StochasticAnalytics / emClarity

GNU Lesser General Public License v3.0
41 stars 6 forks source link

Masks, filters and optional inputs. #111

Closed thomasfrosio closed 4 years ago

thomasfrosio commented 4 years ago

1) OPTIONAL inputs

Managing inputs and especially inputs with default values is not very straightforward in MATLAB. I chose to add a parameter (OPTIONAL ; cell|struct) at the end of functions that have optional inputs. I deliberately didn't use varargin because it embeds the inputs into a cell, which makes everything more complicated when passing this OPTIONAL structure from one function to another. As a result, if one wants to use the default values, an empty cell or empty structure should be specified as the OPTIONAL parameter (ex: EMC_resize(image, limits, 'cpu', {} )). It's kinda ugly but at least it is clear. Then, if one wants to change the origin and deactivate the taper for example: EMC_resize(image, limits, 'cpu', {'origin',1 ; 'taper', false}). To be sure that all parameters are known and used by the function, EMC_extract_optional, called by the functions during checkIN, is going to make sure all the fields are known and will raise an error if a field is not expected by the function. I also added a small function, named EMC_is3d: [flg3d, ndim] = EMC_is3d(SIZE). As commented in my last pull request by Ben, z=1 is considered as 2d in every EMC functions, from now.

2) BH_padZeros3d

BH_padZeros3d is changed to EMC_resize:

taper_one_side

pad_fft_output

--The taper is not correctly applied with fourierOverSample: it is applied to the edges and not the center, which doesn't make sense to me because the center is padded in this case.

taper_fftoutput

3) BH_bandpass3d

Equivalent to EMC_bandpass. One general question is: Should we start to roll off at the specified cutoff or start the roll off slightly before to have the cutoff value at a given threshold (like 0.8)? At 80%, the frequencies are still playing their part is the transform and it would allow to make 'more accurate' bandpass. By that I mean that it would let less unwanted frequencies go through while keeping most of the wanted frequencies. Just an idea... Also, when the cutoff is at Nyquist, should we start to roll off at Nyquist (like we currently do), start the roll off before to end up at Nyquist or start the roll off before to have Nyquist at a given threshold (like 50%), considering that frequencies below this threshold won't be 'visible' in the reconstruction? Of course, this should be tested, it's worth mentioning though.

bandpass_comparaison

extended_rolloff

4) BH_mask3d

BH_mask3d.m is spitted into two different functions: EMC_shapeMask, that manages 'shape' masks, so 'cylinder', 'sphere' and 'rectangle'. For binary masks, EMC_binaryMask will be the function to call. EMC_binaryMask will be the subject of another pull request.

Future pull requests in December 2019.

bHimes commented 4 years ago

This all looks nice! I'll review more closely on Monday. Regarding the bandpass for non square images, maybe double check the definition of the discrete fft...I'm not convinced there is bug.

Ben

On Sun, Dec 15, 2019, 8:36 AM ffyr2w notifications@github.com wrote:

1) OPTIONAL inputs

Managing inputs and especially inputs with default values is not very straightforward in MATLAB. I chose to add a parameter (OPTIONAL ; cell|struct) at the end of functions that have optional inputs. I deliberately didn't use varargin because it embeds the inputs into a cell, which makes everything more complicated when passing this OPTIONAL structure from one function to another. As a result, if one wants to use the default values, an empty cell or empty structure should be specified as the OPTIONAL parameter (ex: EMC_resize(image, limits, 'cpu', {} )). It's kinda ugly but at least it is clear. Then, if one wants to change the origin and deactivate the taper for example: EMC_resize(image, limits, 'cpu', {'origin',1 ; 'taper', false}). To be sure that all parameters are known and used by the function, EMC_extract_optional, called by the functions during checkIN, is going to make sure all the fields are known and will raise an error if a field is not expected by the function. I also added a small function, named EMC_is3d: [flg3d, ndim] = EMC_is3d(SIZE). As commented in my last pull request by Ben, z=1 is considered as 2d in every EMC functions, from now. 1) BH_padZeros3d

BH_padZeros3d is changed to EMC_resize:

  • Taper type: the user can specify the type (cosine or linear) and length of the taper (ex: 'taper', {'cosine', 7}). In the next pull request I will add Gaussian tapers. One can also send their own taper, but this is probably temporary and might be removed once Gaussian taper will be added.
  • The taper is applied only to the edges that are padded. It makes more sense to me, but it might not be a good idea.

[image: taper_one_side] https://user-images.githubusercontent.com/51373921/70863107-8c443300-1f3c-11ea-8419-7f9445c988e3.png

  • The algorithm for cropping is changed. Now padding and cropping can be performed at the same time (3d or 2d). In term of performance, EMC_resize is as fast or faster than BH_padZeros3d, depending if a taper should be applied, or if only padding or cropping is performed.
  • LIMITS: BH_padZeros3d uses padlow and padtop to define the edges that should be cropped or padded. EMC_resize uses LIMITS (see EMC_multi_limits, analogue to BH_multi_padVal). LIMITS is as follow: [x_low, x_high, y_low, y_high], or for 3d: [x_low, x_top, y_low, y_top, z_low, z_top]. So LIMITS = [1,2,3,4] is equivalent to PADLOW=[1,3]; PADTOP=[2,4].
  • But fixes (?): --BH_padZeros3d has the flag 'fourierOverSample' which triggers padding on fft outputs (zero frequency first). EMC functions uses the 'origin' field-name to specify the origin (origin=-1 specifies zero frequency first for example). Anyway, I am not sure about the output of fourierOverSample (it looks like it is shifted by one pixel and is not correctly re-centered by fftshift). This shift is also there in BH_multi_gridCoordinates with flgShiftOrigin=0.

[image: pad_fft_output] https://user-images.githubusercontent.com/51373921/70863117-b269d300-1f3c-11ea-8ed9-f4b399e40f4d.png

--The taper is not correctly applied with fourierOverSample: it is applied to the edges and not the center, which doesn't make sense to me because the center is padded in this case.

[image: taper_fftoutput] https://user-images.githubusercontent.com/51373921/70863121-b990e100-1f3c-11ea-84db-40c5f0c11540.png 2) BH_bandpass3d

Equivalent to EMC_bandpass. One general question is: Should we start to roll off at the specified cutoff or start the roll off slightly before to have the cutoff value at a given threshold (like 0.8)? At 80%, the frequencies are still playing their part is the transform and it would allow to make 'more accurate' bandpass. By that I mean that it would let less unwanted frequencies go through while keeping most of the wanted frequencies. Just an idea... Also, when the cutoff is at Nyquist, should we start to roll off at Nyquist (like we currently do), start the roll off before to end up at Nyquist or start the roll off before to have Nyquist at a given threshold (like 50%), considering that frequencies below this threshold won't be 'visible' in the reconstruction? Of course, this should be tested, it's worth mentioning though.

  • The biggest change in EMC_bandpass compared to BH_bandpass3d is the gaussian roll off used. It uses a roll off that is invariant to the dimensions of the image/volume. BH_bandpass3d will always generate a taper of ~7pixels, which is almost useless (to prevent the ringing effect) for dimensions larger than ~100/200 pixels.

[image: bandpass_comparaison] https://user-images.githubusercontent.com/51373921/70863124-c57ca300-1f3c-11ea-8b3a-14c3ecf1a04f.png

  • BH_bandpass3d allows 'extended' roll offs: EMC_bandpass do the same to the exception that the roll off stops at Nyquist for the high frequency extended roll off. Although, if the cutoff is very close to Nyquist, to the point where the extended roll off is smaller that the default roll off, it switches back to the default Gaussian taper to prevent ringing effects.

[image: extended_rolloff] https://user-images.githubusercontent.com/51373921/70863132-d6c5af80-1f3c-11ea-9ab4-0c22819701e1.png

-

One can of course change the strength of the Gaussian roll and switch to an extended roll. lowroll and highroll are independent from each other and accessible as OPTIONAL inputs.

I changed the naming convention: LOWCUT is the cutoff at lower frequency/resolution and HIGHCUT at higher frequencies (higher than LOWCUT). They can be 'nan', which will compute lowpass/highpass filters.

Bug fix (?): BH_bandpass3d doesn't correctly work with uneven dimensions. EMC_bandpass fixes this by using 'isotropic' radial grids (see EMC_multi_vectorCoordinates for more details).

3) BH_mask3d

BH_mask3d.m is spitted into two different functions: EMC_shapeMask, that manages 'shape' masks, so 'cylinder', 'sphere' and 'rectangle'. For binary masks, EMC_binaryMask will be the function to call. EMC_binaryMask will be the subject of another pull request.

  • EMC_shapeMask uses the same gaussian rolls that EMC_bandpass, which is also much more efficient that what is currently implemented in BH_mask3d. The algorithm to generate the fixed ~7pixels taper in BH_mask3d is very inefficient and become almost unusable if the taper is ~100 pixels. As with EMC_bandpass, the gaussian roll off in EMC_shapeMask is invariant to the length of the dimension and can be easily adjusted with OPTIONAL.
  • Rectangular masks now have the exact same gaussian taper as spherical and cylindrical masks.
  • Currently, 3d cylinders are not supported: I am finishing a new algorithm that should efficiently deal with cylinders and will also remove the need to call EMC_resize for rectangular masks.

Future pull requests in December 2019.

  • Binary masks
  • Update EMC_shapeMask to deal with cylinders and add asymmetric restriction.
  • unit tests

You can view, comment on, or merge this pull request online at:

https://github.com/bHimes/emClarity/pull/111 Commit Summary

  • Allow z=1 to be processed as 2d; Fix: TRANS can be a structure and still checked-in
  • EMC_multi_vectorCoordinates.m: Use OPTIONAL and added the isotrope flag (used by bandpass). Bug fixed when computing half vectors for origin=0 and even nb of pixel. Shifts are now not allowed for origin=0 and half=true. All the vectors are now single floats.
  • coordinates/EMC_multi_gridMasks.m: Now use the OPTIONAL input. Isotropic vectors are now directly made by EMC_multi_vectorCoordinates. TYPE supports low or upper cases.
  • masking/EMC_resize.m: EMC version of BH_padZeros3d. Fully supports cropping and padding with origin=-1|1|2 and allows to specify your own taper.
  • masking/EMC_bandpass.m: EMC version of BH_bandpass3d. Isotropic and gaussian roll off invariant to the filter size.
  • masking/EMC_shapeMask.m: EMC version of BH_mask3d for 'cylinder', 'sphere' and 'rectangle'. Gaussian taper is invariant to the image/volume size. Asymmetric restriction not supported yet, as well as 3d cylinders.
  • coordinates/EMC_multi_gridCoordinates.m: flags are now in TRANS
  • EMC_multi_taper.m: create cosine and linear taper.
  • Utilitaries functions for OPTIONAL inputs and checkIN.
  • Merge branch 'mask3d' into development
  • I messed up with the merge on this one. It's fixed now.
  • coordinates/EMC_multi_gridCoordinates.m: fix error message (direction in lower case...)

File Changes

Patch Links:

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/bHimes/emClarity/pull/111?email_source=notifications&email_token=AAXPLGDUBRH34CKZFH5TRR3QYYXD7A5CNFSM4J3AEZD2YY3PNVWWK3TUL52HS4DFUVEXG43VMWVGG33NNVSW45C7NFSM4IASA26Q, or unsubscribe https://github.com/notifications/unsubscribe-auth/AAXPLGBKEL235W3EOGEZZDTQYYXD7ANCNFSM4J3AEZDQ .

thomasfrosio commented 4 years ago

I am almost done with cylinder masks, so I'll update EMC_shapeMask.m tomorrow.

For the fft, I am following this (see p.10): https://www.iap.uni-jena.de/iapmedia/de/Lecture/Computational+Photonics1472680800/CPho16_Seminar7_FFT_Primer.pdf It matches with the output I have on matlab: with real inputs, the hermitian symmetry makes it easy to spot the zero frequency. But it doesn't mean that what you do is not correct so I'll test it with one image to see if there is a difference between what you do and what I do.

Thomas

thomasfrosio commented 4 years ago

Hum, I might have misunderstood. For non square images, if you have doubt I will definitely check...

thomasfrosio commented 4 years ago

I just realized that the 'rectangle' mask is not correct if the shifts push the rectangle to the edge of the image or if, radius+(2*length(taper))+1, is equal/larger than the mask size. The reason is EMC_resize is not going to apply the taper to the edges that are not padded. I am currently adding an option to EMC_resize ({'force_taper', true}) to apply the taper even when padding is not required. To follow BH_padZeros3d, this option will tape AFTER cropping.

thomasfrosio commented 4 years ago

After looking at the discrete fft equation and comparing results, your version of BH_bandpass3d.m is correct for non-squared image. I have adjusted EMC_bandpass and some other fixes in EMC_shapeMask and EMC_resize. EMC_shapeMask can deal with cylinders now. I'll update everything during the evening.

bHimes commented 4 years ago

This all sounds good after a second read through. I'm glad you worked out the DFT issue. I guess it was probably that the frequency step of a Fourier voxel is not isotropic for a non-square image, since it is related to the size of the image - i.e. for M x N each fourier pixel is 1/(Mapix) and 1/(Napix)