Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Masks, filters and optional inputs. #111

Merged
merged 13 commits into from
Dec 19, 2019

Conversation

thomasfrosio
Copy link

@thomasfrosio thomasfrosio commented Dec 15, 2019

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 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.

taper_one_side

  • 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.

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.

  • 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.

bandpass_comparaison

  • 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.

extended_rolloff

  • 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).

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.

  • 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

…ag (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.
…pic vectors are now directly made by EMC_multi_vectorCoordinates. TYPE supports low or upper cases.
…opping and padding with origin=-1|1|2 and allows to specify your own taper.
…aussian roll off invariant to the filter size.
…here' and 'rectangle'. Gaussian taper is invariant to the image/volume size. Asymmetric restriction not supported yet, as well as 3d cylinders.
@bHimes
Copy link
Collaborator

bHimes commented Dec 15, 2019 via email

@thomasfrosio
Copy link
Author

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
Copy link
Author

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

@thomasfrosio
Copy link
Author

thomasfrosio commented Dec 16, 2019

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.

… not padded or edges that are cropped with the new OPTIONAL input: 'force_taper' (bool).

Fix a stupid bug with origin=-1 with cropping: I forgot that MATLAB had linear indexing...
Now, even if no padding, no cropping and no taper, the function will change the precision if asked with the OPTIONAL input 'precision'.
@thomasfrosio
Copy link
Author

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
Copy link
Collaborator

bHimes commented Dec 19, 2019

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)

@bHimes bHimes merged commit e55f101 into StochasticAnalytics:development Dec 19, 2019
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants