Frequency Filtering in Time Dimension?

classic Classic list List threaded Threaded
15 messages Options
Reply | Threaded
Open this post in threaded view
|

Frequency Filtering in Time Dimension?

Jacob Keller-2
Hi All,

is there a way to do high-pass/low-pass filtering in the time dimension? Or
a suggestion for a workaround?

All the best,

Jacob Keller

--
ImageJ mailing list: http://imagej.nih.gov/ij/list.html
Reply | Threaded
Open this post in threaded view
|

Re: Frequency Filtering in Time Dimension?

Herbie
Good day Jacob,

if you mean point-wise filtering in time this can of course be done
either by 1D-convolution or 1D-FFT. The latter is made easy by the
following macro command:

*Array.fourier(array, windowType)*
Calculates and returns the Fourier amplitudes of array. WindowType can
be "none", "Hamming", "Hann", or "flat-top", or may be omitted (meaning
"none"). See the TestArrayFourier macro for an example and more
documentation. Requires 1.49i.

HTH

Herbie

::::::::::::::::::::::::::::::::::::::::::
Am 29.10.18 um 20:51 schrieb Jacob Keller:

> Hi All,
>
> is there a way to do high-pass/low-pass filtering in the time dimension? Or
> a suggestion for a workaround?
>
> All the best,
>
> Jacob Keller
>
> --
> ImageJ mailing list: http://imagej.nih.gov/ij/list.html
>

--
ImageJ mailing list: http://imagej.nih.gov/ij/list.html
Reply | Threaded
Open this post in threaded view
|

Re: Frequency Filtering in Time Dimension?

Michael Schmid-3
In reply to this post by Jacob Keller-2
Hi Jacob,

do you have simply a 1D data set or an image stack, with slices for
different times?

- In the first case, convert your data to an image with n * 1 pixels,
and apply one of the usual image filters.

- In the second case, you can use Process>Filters>Gaussian Blur 3D and
specify the x and y sigma as zero.
For high-pass filtering, you would have to duplicate the stack first and
then subtract the filtered image.

If you rather want moving averages or a median, you can reslice the
stack (Image>Stacks>Reslice) to have the time direction in x or y and
apply the 'Fast Filters' plugin, which can do 1D filtering. Then Reslice
to make time the z axis again.
The 'Fast Filters' plugin also has an option to subtract the filtered
image (i.e., highpass operation)


Michael
________________________________________________________________
On 29.10.18 20:51, Jacob Keller wrote:
> Hi All,
>
> is there a way to do high-pass/low-pass filtering in the time dimension? Or
> a suggestion for a workaround?
>
> All the best,
>
> Jacob Keller

--
ImageJ mailing list: http://imagej.nih.gov/ij/list.html
Reply | Threaded
Open this post in threaded view
|

Re: Frequency Filtering in Time Dimension?

Herbie
In reply to this post by Jacob Keller-2
Jacob,

I must admit that using "*Array.fourier" won't help because it doesn't
deal with complex-valued signals. (I wasn't aware of this fact.)

Consequently, you need to follow one of Michael's suggestions.

Good luck

Herbie

::::::::::::::::::::::::::::::::::::::::::
Good day Jacob,

if you mean point-wise filtering in time this can of course be done
either by 1D-convolution or 1D-FFT. The latter is made easy by the
following macro command:

*Array.fourier(array, windowType)*
Calculates and returns the Fourier amplitudes of array. WindowType can
be "none", "Hamming", "Hann", or "flat-top", or may be omitted (meaning
"none"). See the TestArrayFourier macro for an example and more
documentation. Requires 1.49i.

HTH

Herbie

::::::::::::::::::::::::::::::::::::::::::
Am 29.10.18 um 20:51 schrieb Jacob Keller:

> Hi All,
>
> is there a way to do high-pass/low-pass filtering in the time dimension? Or
> a suggestion for a workaround?
>
> All the best,
>
> Jacob Keller
>
> --
> ImageJ mailing list: http://imagej.nih.gov/ij/list.html
>

--
ImageJ mailing list: http://imagej.nih.gov/ij/list.html

--
ImageJ mailing list: http://imagej.nih.gov/ij/list.html
Reply | Threaded
Open this post in threaded view
|

Re: Frequency Filtering in Time Dimension?

Herbie
In reply to this post by Jacob Keller-2
Jacob,

I must admit that using "*Array.fourier" won't help because it doesn't deal
with complex-valued signals. (I wasn't aware of this fact.)

Consequently, you need to follow one of Michael's suggestions.

Good luck

Herbie

::::::::::::::::::::::::::::::::::::::::::
Good day Jacob,

if you mean point-wise filtering in time this can of course be done either
by 1D-convolution or 1D-FFT. The latter is made easy by the following macro
command:

*Array.fourier(array, windowType)*
Calculates and returns the Fourier amplitudes of array. WindowType can be
"none", "Hamming", "Hann", or "flat-top", or may be omitted (meaning
"none"). See the TestArrayFourier macro for an example and more
documentation. Requires 1.49i.

HTH

Herbie

::::::::::::::::::::::::::::::::::::::::::

Jacob Keller-2 wrote

> Hi All,
>
> is there a way to do high-pass/low-pass filtering in the time dimension?
> Or
> a suggestion for a workaround?
>
> All the best,
>
> Jacob Keller
>
> --
> ImageJ mailing list: http://imagej.nih.gov/ij/list.html





--
Sent from: http://imagej.1557.x6.nabble.com/

--
ImageJ mailing list: http://imagej.nih.gov/ij/list.html
Reply | Threaded
Open this post in threaded view
|

Re: Frequency Filtering in Time Dimension?

Jacob Keller-2
In reply to this post by Michael Schmid-3
>
> do you have simply a 1D data set or an image stack, with slices for
> different times?
>

The latter--image stack, one plane over time.

- In the second case, you can use Process>Filters>Gaussian Blur 3D and
> specify the x and y sigma as zero.
> For high-pass filtering, you would have to duplicate the stack first and
> then subtract the filtered image.
>

Ah, this is a great idea--I will try it out. How can I figure out the
relationship between sigma and the desired frequency cutoff? For example,
let's say the stack has 120 frames per period--what sigma value should be
input?

Thanks very much for these suggestions--I think they are going to be
very helpful,

Jacob

If you rather want moving averages or a median, you can reslice the

> stack (Image>Stacks>Reslice) to have the time direction in x or y and
> apply the 'Fast Filters' plugin, which can do 1D filtering. Then Reslice
> to make time the z axis again.
> The 'Fast Filters' plugin also has an option to subtract the filtered
> image (i.e., highpass operation)
>
>
> Michael
> ________________________________________________________________
> On 29.10.18 20:51, Jacob Keller wrote:
> > Hi All,
> >
> > is there a way to do high-pass/low-pass filtering in the time dimension?
> Or
> > a suggestion for a workaround?
> >
> > All the best,
> >
> > Jacob Keller
>
> --
> ImageJ mailing list: http://imagej.nih.gov/ij/list.html
>

--
ImageJ mailing list: http://imagej.nih.gov/ij/list.html
Reply | Threaded
Open this post in threaded view
|

Re: Frequency Filtering in Time Dimension?

Michael Schmid-3
Hi Jacob,

when you are doing a Gaussian Blur, in the Fourier domain it corresponds
to multiplying the amplitudes with a Gaussian as well.

A Gaussian
    exp(-pi * x²/a²)
in real space corresponds to
    exp(-pi * f² * a²)
in Fourier space

If you want to attenuate a frequency component f1 by a factor 1/e =
0.37, you thus need a² = 1/(pi f1²), i.e. exp (pi² * f1² * x²)

In terms of sigma, a Gaussian is given by
   exp(-x²/(2*sigma²))

Thus, attenuation by a factor of 1/e at f1 corresponds to
   sigma = 1/(pi * f1 * sqrt(2))


You can easily check it with the following macro:

period = 16;  // between about 10 and 100
               // best accuracy with powers of 2
f=1/period;   // spatial frequency
newImage("Untitled", "32-bit white", 256, 256, 1);
run("Macro...", "code=v=sin(2*PI*x*"+f+")");  //create sine wave
sigma = 1/(PI*f*sqrt(2));
run("Gaussian Blur...", "sigma="+sigma);
run("Set Measurements...", "min display redirect=None decimal=4");
makeRectangle(25, 0, 200, 256);
run("Measure");


You will see that the maxima and minima of the sine wave (which has had
an original amplitude of 1) gets attenuated to an amplitude of about 0.37.

If you take sigma = 1/(pi * f1) you will get an attenuation factor of
1/e² = 0.13 at the frequency f1.

So far an excursion into the math of Gaussians and their Fourier
transform...


Michael
________________________________________________________________
On 30.10.18 16:23, Jacob Keller wrote:
 >>
 >> do you have simply a 1D data set or an image stack, with slices for
 >> different times?
 >>
 >
 > The latter--image stack, one plane over time.
 >
 > - In the second case, you can use Process>Filters>Gaussian Blur 3D and
 >> specify the x and y sigma as zero.
 >> For high-pass filtering, you would have to duplicate the stack first and
 >> then subtract the filtered image.
 >>
 >
 > Ah, this is a great idea--I will try it out. How can I figure out the
 > relationship between sigma and the desired frequency cutoff? For example,
 > let's say the stack has 120 frames per period--what sigma value should be
 > input?
 >
 > Thanks very much for these suggestions--I think they are going to be
 > very helpful,
 >
 > Jacob
 >
 > If you rather want moving averages or a median, you can reslice the
 >> stack (Image>Stacks>Reslice) to have the time direction in x or y and
 >> apply the 'Fast Filters' plugin, which can do 1D filtering. Then Reslice
 >> to make time the z axis again.
 >> The 'Fast Filters' plugin also has an option to subtract the filtered
 >> image (i.e., highpass operation)
 >>
 >>
 >> Michael
 >> ________________________________________________________________
 >> On 29.10.18 20:51, Jacob Keller wrote:
 >>> Hi All,
 >>>
 >>> is there a way to do high-pass/low-pass filtering in the time
dimension?
 >> Or
 >>> a suggestion for a workaround?
 >>>
 >>> All the best,
 >>>
 >>> Jacob Keller
 >>
 >> --
 >> ImageJ mailing list: http://imagej.nih.gov/ij/list.html
 >>
 >
 > --
 > ImageJ mailing list: http://imagej.nih.gov/ij/list.html
 >

--
ImageJ mailing list: http://imagej.nih.gov/ij/list.html
Reply | Threaded
Open this post in threaded view
|

Re: Frequency Filtering in Time Dimension?

Jacob Keller-2
Thanks Michael--this is really clear and helpful! I checked out some Bode
plots of various filters, and it seems Gaussian is pretty non-brick wall. I
wonder whether a Butterworth or Bessel filter has been implemented?

JPK

On Wed, Oct 31, 2018 at 7:36 AM Michael Schmid <[hidden email]>
wrote:

> Hi Jacob,
>
> when you are doing a Gaussian Blur, in the Fourier domain it corresponds
> to multiplying the amplitudes with a Gaussian as well.
>
> A Gaussian
>     exp(-pi * x²/a²)
> in real space corresponds to
>     exp(-pi * f² * a²)
> in Fourier space
>
> If you want to attenuate a frequency component f1 by a factor 1/e =
> 0.37, you thus need a² = 1/(pi f1²), i.e. exp (pi² * f1² * x²)
>
> In terms of sigma, a Gaussian is given by
>    exp(-x²/(2*sigma²))
>
> Thus, attenuation by a factor of 1/e at f1 corresponds to
>    sigma = 1/(pi * f1 * sqrt(2))
>
>
> You can easily check it with the following macro:
>
> period = 16;  // between about 10 and 100
>                // best accuracy with powers of 2
> f=1/period;   // spatial frequency
> newImage("Untitled", "32-bit white", 256, 256, 1);
> run("Macro...", "code=v=sin(2*PI*x*"+f+")");  //create sine wave
> sigma = 1/(PI*f*sqrt(2));
> run("Gaussian Blur...", "sigma="+sigma);
> run("Set Measurements...", "min display redirect=None decimal=4");
> makeRectangle(25, 0, 200, 256);
> run("Measure");
>
>
> You will see that the maxima and minima of the sine wave (which has had
> an original amplitude of 1) gets attenuated to an amplitude of about 0.37.
>
> If you take sigma = 1/(pi * f1) you will get an attenuation factor of
> 1/e² = 0.13 at the frequency f1.
>
> So far an excursion into the math of Gaussians and their Fourier
> transform...
>
>
> Michael
> ________________________________________________________________
> On 30.10.18 16:23, Jacob Keller wrote:
>  >>
>  >> do you have simply a 1D data set or an image stack, with slices for
>  >> different times?
>  >>
>  >
>  > The latter--image stack, one plane over time.
>  >
>  > - In the second case, you can use Process>Filters>Gaussian Blur 3D and
>  >> specify the x and y sigma as zero.
>  >> For high-pass filtering, you would have to duplicate the stack first
> and
>  >> then subtract the filtered image.
>  >>
>  >
>  > Ah, this is a great idea--I will try it out. How can I figure out the
>  > relationship between sigma and the desired frequency cutoff? For
> example,
>  > let's say the stack has 120 frames per period--what sigma value should
> be
>  > input?
>  >
>  > Thanks very much for these suggestions--I think they are going to be
>  > very helpful,
>  >
>  > Jacob
>  >
>  > If you rather want moving averages or a median, you can reslice the
>  >> stack (Image>Stacks>Reslice) to have the time direction in x or y and
>  >> apply the 'Fast Filters' plugin, which can do 1D filtering. Then
> Reslice
>  >> to make time the z axis again.
>  >> The 'Fast Filters' plugin also has an option to subtract the filtered
>  >> image (i.e., highpass operation)
>  >>
>  >>
>  >> Michael
>  >> ________________________________________________________________
>  >> On 29.10.18 20:51, Jacob Keller wrote:
>  >>> Hi All,
>  >>>
>  >>> is there a way to do high-pass/low-pass filtering in the time
> dimension?
>  >> Or
>  >>> a suggestion for a workaround?
>  >>>
>  >>> All the best,
>  >>>
>  >>> Jacob Keller
>  >>
>  >> --
>  >> ImageJ mailing list: http://imagej.nih.gov/ij/list.html
>  >>
>  >
>  > --
>  > ImageJ mailing list: http://imagej.nih.gov/ij/list.html
>  >
>
> --
> ImageJ mailing list: http://imagej.nih.gov/ij/list.html
>

--
ImageJ mailing list: http://imagej.nih.gov/ij/list.html
Reply | Threaded
Open this post in threaded view
|

Re: Frequency Filtering in Time Dimension?

Herbie
Good day Jacob,

Butterworth or Bessel filters are filters that play a significant role
if the filters have to be implemented electronically as *causal* filters
for real-time processing of temporal signals.

For this purpose, Gaussian or hard-limiting filters are difficult to
realize electronically. Not so for the processing of static signals such
as images that are per se a-causal. Hard-limiting low-pass filters are
easy to implement in the Fourier-domain (Just limit the Fourier-spectral
extent) for a-causal signals and you can do this with ImageJ.

Because here and with ImageJ we deal with images, filter-characteristics
that are typical for the processing of causal time-signals don't play a
role and I don't think that they are or will be implemented for use with
ImageJ.

Regards

Herbie
______

PS:
Th mathematical proof that a Gaussian is a so-called Fourier-pair,
i.e.the fact that a Gaussian corresponds to a Gaussian by
Fourier-transformation, is quite involved.

::::::::::::::::::::::::::::::::::::::::::
Am 31.10.18 um 17:32 schrieb Jacob Keller:

> Thanks Michael--this is really clear and helpful! I checked out some Bode
> plots of various filters, and it seems Gaussian is pretty non-brick wall. I
> wonder whether a Butterworth or Bessel filter has been implemented?
>
> JPK
>
> On Wed, Oct 31, 2018 at 7:36 AM Michael Schmid <[hidden email]>
> wrote:
>
>> Hi Jacob,
>>
>> when you are doing a Gaussian Blur, in the Fourier domain it corresponds
>> to multiplying the amplitudes with a Gaussian as well.
>>
>> A Gaussian
>>      exp(-pi * x²/a²)
>> in real space corresponds to
>>      exp(-pi * f² * a²)
>> in Fourier space
>>
>> If you want to attenuate a frequency component f1 by a factor 1/e =
>> 0.37, you thus need a² = 1/(pi f1²), i.e. exp (pi² * f1² * x²)
>>
>> In terms of sigma, a Gaussian is given by
>>     exp(-x²/(2*sigma²))
>>
>> Thus, attenuation by a factor of 1/e at f1 corresponds to
>>     sigma = 1/(pi * f1 * sqrt(2))
>>
>>
>> You can easily check it with the following macro:
>>
>> period = 16;  // between about 10 and 100
>>                 // best accuracy with powers of 2
>> f=1/period;   // spatial frequency
>> newImage("Untitled", "32-bit white", 256, 256, 1);
>> run("Macro...", "code=v=sin(2*PI*x*"+f+")");  //create sine wave
>> sigma = 1/(PI*f*sqrt(2));
>> run("Gaussian Blur...", "sigma="+sigma);
>> run("Set Measurements...", "min display redirect=None decimal=4");
>> makeRectangle(25, 0, 200, 256);
>> run("Measure");
>>
>>
>> You will see that the maxima and minima of the sine wave (which has had
>> an original amplitude of 1) gets attenuated to an amplitude of about 0.37.
>>
>> If you take sigma = 1/(pi * f1) you will get an attenuation factor of
>> 1/e² = 0.13 at the frequency f1.
>>
>> So far an excursion into the math of Gaussians and their Fourier
>> transform...
>>
>>
>> Michael
>> ________________________________________________________________
>> On 30.10.18 16:23, Jacob Keller wrote:
>>   >>
>>   >> do you have simply a 1D data set or an image stack, with slices for
>>   >> different times?
>>   >>
>>   >
>>   > The latter--image stack, one plane over time.
>>   >
>>   > - In the second case, you can use Process>Filters>Gaussian Blur 3D and
>>   >> specify the x and y sigma as zero.
>>   >> For high-pass filtering, you would have to duplicate the stack first
>> and
>>   >> then subtract the filtered image.
>>   >>
>>   >
>>   > Ah, this is a great idea--I will try it out. How can I figure out the
>>   > relationship between sigma and the desired frequency cutoff? For
>> example,
>>   > let's say the stack has 120 frames per period--what sigma value should
>> be
>>   > input?
>>   >
>>   > Thanks very much for these suggestions--I think they are going to be
>>   > very helpful,
>>   >
>>   > Jacob
>>   >
>>   > If you rather want moving averages or a median, you can reslice the
>>   >> stack (Image>Stacks>Reslice) to have the time direction in x or y and
>>   >> apply the 'Fast Filters' plugin, which can do 1D filtering. Then
>> Reslice
>>   >> to make time the z axis again.
>>   >> The 'Fast Filters' plugin also has an option to subtract the filtered
>>   >> image (i.e., highpass operation)
>>   >>
>>   >>
>>   >> Michael
>>   >> ________________________________________________________________
>>   >> On 29.10.18 20:51, Jacob Keller wrote:
>>   >>> Hi All,
>>   >>>
>>   >>> is there a way to do high-pass/low-pass filtering in the time
>> dimension?
>>   >> Or
>>   >>> a suggestion for a workaround?
>>   >>>
>>   >>> All the best,
>>   >>>
>>   >>> Jacob Keller
>>   >>
>>   >> --
>>   >> ImageJ mailing list: http://imagej.nih.gov/ij/list.html
>>   >>
>>   >
>>   > --
>>   > ImageJ mailing list: http://imagej.nih.gov/ij/list.html
>>   >
>>
>> --
>> ImageJ mailing list: http://imagej.nih.gov/ij/list.html
>>
>
> --
> ImageJ mailing list: http://imagej.nih.gov/ij/list.html
>

--
ImageJ mailing list: http://imagej.nih.gov/ij/list.html
Reply | Threaded
Open this post in threaded view
|

Re: Frequency Filtering in Time Dimension?

Jacob Keller-2
>
> Because here and with ImageJ we deal with images, filter-characteristics
> that are typical for the processing of causal time-signals don't play a
> role and I don't think that they are or will be implemented for use with
> ImageJ.
>

It seems that you are saying that ImageJ is limited in scope to static
images? Huh? I see a lot of folks using it for time series, and there are
quite a few plugins that incorporate time. I have noticed generally a bit
of a split in the community regarding time series, but it seems to me that
time is a critical component of biological processes, and at the level of
microscopy is quite accessible, so why limit imageJ to static images? Have
I misunderstood you?

JPK

--
ImageJ mailing list: http://imagej.nih.gov/ij/list.html
Reply | Threaded
Open this post in threaded view
|

Re: Frequency Filtering in Time Dimension?

Herbie
In reply to this post by Herbie
Jacob,

I didn't want to imply this:

"It seems that you are saying that ImageJ is limited in scope to static
images"

What I had in mind are real-time signals.

"Have I misunderstood you?"

Yes.

ImageJ isn't capable of dealing with real-time temporal signals.

If you have a temporal sequence of images available, the sequence is no
longer temporal, it is a static stack of images and all the problems
occurring with causal temporal processing are irrelevant because the
stack of images is a-causal from an signal processing point of view. Of
course you know that e.g. slice 2 comes before slice 3 but it is no
problem to have e.g. a convolution kernel that extends from slice 1 to
slice 5. In real-time processing this isn't possible as long as slice 5
is available...

I hope this makes things a bit clearer.

Regards

Herbie

::::::::::::::::::::::::::::::::::::::::::
Am 01.11.18 um 16:40 schrieb Jacob Keller:

>     Because here and with ImageJ we deal with images,
>     filter-characteristics
>     that are typical for the processing of causal time-signals don't play a
>     role and I don't think that they are or will be implemented for use
>     with
>     ImageJ.
>
>
> It seems that you are saying that ImageJ is limited in scope to static
> images? Huh? I see a lot of folks using it for time series, and there
> are quite a few plugins that incorporate time. I have noticed generally
> a bit of a split in the community regarding time series, but it seems to
> me that time is a critical component of biological processes, and at the
> level of microscopy is quite accessible, so why limit imageJ to static
> images? Have I misunderstood you?
>
> Jacob
>
>

--
ImageJ mailing list: http://imagej.nih.gov/ij/list.html
Reply | Threaded
Open this post in threaded view
|

Re: Frequency Filtering in Time Dimension?

Michael Schmid-3
In reply to this post by Jacob Keller-2
Hi Jacob,

do I understand it correctly that you want a sharper (more square-like)
frequency response than a Gaussian?
In that case, you will get some ringing in the data, i.e., the step
response will have oscillations around the (smoothed) step, and the
response to a delta function (single peak) will reach below the baseline
for some times.

If you mainly want to suppress *one* frequency, you can use moving
averages, with a kernel length equal to the period you want to subtract
(i.e., a rectangular kernel; this gives you no oscillations). To do it
in time (slice) direction, reslice the stack and use the 'Fast Filters'
plugin. The kernel length is 2*radius+1. Maybe apply some weak Gaussian
filter in addition.

If you want a more square-like filter in the Fourier domain (and you
don't care about ringing), and the desired cutoff frequency is rather
high, you can use Process>Filter>Convolve on the resliced stack. For a
perfectly sharp frequency cutoff you would need a sin(a*x)/(a*x) kernel
(the so-called "sinc" function); but this would give you infinite kernel
size.
   https://en.wikipedia.org/wiki/Sinc_function
Multiply it with a Gaussian to limit the kernel size; this makes the
cutoff a bit smoother (in frequency space, it gets convolved, i.e.,
blurred, by the Fourier transform of the Gaussian that you have used for
multiplying). The narrower the Gaussian you multiply with the sync, the
smaller a kernel you can use (if you neglect very small numbers in the
periphery of the kernel), but the frequency cutoff becomes smoother.
Maybe a Gaussian with a sigma between the first and second zero is a
good starting point (i.e., sigma between a*pi and a*2pi).

If you want to check how it behaves, e.g. create a 32-bit image with a
single vertical line, apply the kernel (in x direction, i.e., have a
single line with all the numbers in the text box of the ImageJ
convolver), and look at the power spectrum that you get with the ImageJ
FFT (you need to enable the power spectrum in the FFT options). Maybe
you find examples on some web site (I did not find any in a very quick
search).

In contrast to Bessel or Butterworth filters, the response to a delta
function of a sync*Gaussian filter is symmetric in time. Thus, the time
of any symmetric peak or (anti-symmetric) step function, etc. does not
get not shifted by the filter.

Bessel or Butterworth filters can be implemented in real time (e.g. with
analog electronics), so due to causality, their response is always
*after* the signal. This means that peaks, steps, etc. get shifted to
later times.
This is just an other way to express what Herbie has talked in his
recent mails as "filtering in (real) time" vs. "filtering an image
stack".


Michael
________________________________________________________________


On 2018-10-31 17:32, Jacob Keller wrote:

> Thanks Michael--this is really clear and helpful! I checked out some
> Bode
> plots of various filters, and it seems Gaussian is pretty non-brick
> wall. I
> wonder whether a Butterworth or Bessel filter has been implemented?
>
> JPK
>
> On Wed, Oct 31, 2018 at 7:36 AM Michael Schmid
> <[hidden email]>
> wrote:
>
>> Hi Jacob,
>>
>> when you are doing a Gaussian Blur, in the Fourier domain it
>> corresponds
>> to multiplying the amplitudes with a Gaussian as well.
>>
>> A Gaussian
>>     exp(-pi * x²/a²)
>> in real space corresponds to
>>     exp(-pi * f² * a²)
>> in Fourier space
>>
>> If you want to attenuate a frequency component f1 by a factor 1/e =
>> 0.37, you thus need a² = 1/(pi f1²), i.e. exp (pi² * f1² * x²)
>>
>> In terms of sigma, a Gaussian is given by
>>    exp(-x²/(2*sigma²))
>>
>> Thus, attenuation by a factor of 1/e at f1 corresponds to
>>    sigma = 1/(pi * f1 * sqrt(2))
>>
>>
>> You can easily check it with the following macro:
>>
>> period = 16;  // between about 10 and 100
>>                // best accuracy with powers of 2
>> f=1/period;   // spatial frequency
>> newImage("Untitled", "32-bit white", 256, 256, 1);
>> run("Macro...", "code=v=sin(2*PI*x*"+f+")");  //create sine wave
>> sigma = 1/(PI*f*sqrt(2));
>> run("Gaussian Blur...", "sigma="+sigma);
>> run("Set Measurements...", "min display redirect=None decimal=4");
>> makeRectangle(25, 0, 200, 256);
>> run("Measure");
>>
>>
>> You will see that the maxima and minima of the sine wave (which has
>> had
>> an original amplitude of 1) gets attenuated to an amplitude of about
>> 0.37.
>>
>> If you take sigma = 1/(pi * f1) you will get an attenuation factor of
>> 1/e² = 0.13 at the frequency f1.
>>
>> So far an excursion into the math of Gaussians and their Fourier
>> transform...
>>
>>
>> Michael
>> ________________________________________________________________
>> On 30.10.18 16:23, Jacob Keller wrote:
>>  >>
>>  >> do you have simply a 1D data set or an image stack, with slices
>> for
>>  >> different times?
>>  >>
>>  >
>>  > The latter--image stack, one plane over time.
>>  >
>>  > - In the second case, you can use Process>Filters>Gaussian Blur 3D
>> and
>>  >> specify the x and y sigma as zero.
>>  >> For high-pass filtering, you would have to duplicate the stack
>> first
>> and
>>  >> then subtract the filtered image.
>>  >>
>>  >
>>  > Ah, this is a great idea--I will try it out. How can I figure out
>> the
>>  > relationship between sigma and the desired frequency cutoff? For
>> example,
>>  > let's say the stack has 120 frames per period--what sigma value
>> should
>> be
>>  > input?
>>  >
>>  > Thanks very much for these suggestions--I think they are going to
>> be
>>  > very helpful,
>>  >
>>  > Jacob
>>  >
>>  > If you rather want moving averages or a median, you can reslice the
>>  >> stack (Image>Stacks>Reslice) to have the time direction in x or y
>> and
>>  >> apply the 'Fast Filters' plugin, which can do 1D filtering. Then
>> Reslice
>>  >> to make time the z axis again.
>>  >> The 'Fast Filters' plugin also has an option to subtract the
>> filtered
>>  >> image (i.e., highpass operation)
>>  >>
>>  >>
>>  >> Michael
>>  >> ________________________________________________________________
>>  >> On 29.10.18 20:51, Jacob Keller wrote:
>>  >>> Hi All,
>>  >>>
>>  >>> is there a way to do high-pass/low-pass filtering in the time
>> dimension?
>>  >> Or
>>  >>> a suggestion for a workaround?
>>  >>>
>>  >>> All the best,
>>  >>>
>>  >>> Jacob Keller
>>  >>
>>  >> --
>>  >> ImageJ mailing list: http://imagej.nih.gov/ij/list.html
>>  >>
>>  >
>>  > --
>>  > ImageJ mailing list: http://imagej.nih.gov/ij/list.html
>>  >
>>
>> --
>> ImageJ mailing list: http://imagej.nih.gov/ij/list.html
>>
>
> --
> ImageJ mailing list: http://imagej.nih.gov/ij/list.html

--
ImageJ mailing list: http://imagej.nih.gov/ij/list.html
Reply | Threaded
Open this post in threaded view
|

Re: Frequency Filtering in Time Dimension?

Fred Damen
In reply to this post by Herbie
Herbie...

There are more things in heaven and earth, Horatio,
Than are dreamt of in your philosophy.
- Hamlet (1.5.167-8), Hamlet to Horatio

Depending on the time frame that you are concerned with will dictate if you
look at your data as temporal.  From the image acquisition perspective your
view has merit.  Although, from the study of what is contained in these
sequential volumes of images, e.g., data points, they may very well be
temporal.

In MRI:

fMRI and perfusion are real time data that monitors physiological processes as
they occur, and, are analyzed as temporal signals.

MRE(lastography) are real time data in which the subject of the images is
perturbed in real time and the response, i.e., the images, are analyzed as
temporal responses.

B0 map, T1 map, T2 map, the independent variable is time, but not real time,
albeit, temporal analysis tools are applicable.

Diffusion data is varied by magnetic gradient field strength and timing
parameters, temporal analysis tools are applicable.

I could go on, but suffice it to say that ImageJ is very useful to examine
this data in the frame dimension, and there is a lot on untapped functionality
that can be provided here...

Thanks for listening,

Fred


On Thu, November 1, 2018 11:50 am, Herbie wrote:

> Jacob,
>
> I didn't want to imply this:
>
> "It seems that you are saying that ImageJ is limited in scope to static
> images"
>
> What I had in mind are real-time signals.
>
> "Have I misunderstood you?"
>
> Yes.
>
> ImageJ isn't capable of dealing with real-time temporal signals.
>
> If you have a temporal sequence of images available, the sequence is no
> longer temporal, it is a static stack of images and all the problems
> occurring with causal temporal processing are irrelevant because the
> stack of images is a-causal from an signal processing point of view. Of
> course you know that e.g. slice 2 comes before slice 3 but it is no
> problem to have e.g. a convolution kernel that extends from slice 1 to
> slice 5. In real-time processing this isn't possible as long as slice 5
> is available...
>
> I hope this makes things a bit clearer.
>
> Regards
>
> Herbie
>
> ::::::::::::::::::::::::::::::::::::::::::
> Am 01.11.18 um 16:40 schrieb Jacob Keller:
>>     Because here and with ImageJ we deal with images,
>>     filter-characteristics
>>     that are typical for the processing of causal time-signals don't play a
>>     role and I don't think that they are or will be implemented for use
>>     with
>>     ImageJ.
>>
>>
>> It seems that you are saying that ImageJ is limited in scope to static
>> images? Huh? I see a lot of folks using it for time series, and there
>> are quite a few plugins that incorporate time. I have noticed generally
>> a bit of a split in the community regarding time series, but it seems to
>> me that time is a critical component of biological processes, and at the
>> level of microscopy is quite accessible, so why limit imageJ to static
>> images? Have I misunderstood you?
>>
>> Jacob
>>
>>
>
> --
> ImageJ mailing list: http://imagej.nih.gov/ij/list.html
>

--
ImageJ mailing list: http://imagej.nih.gov/ij/list.html
Reply | Threaded
Open this post in threaded view
|

Re: Frequency Filtering in Time Dimension?

Herbie
Fred,

here and with ImageJ we deal with signal//l/image processing not in
real-time. That's the only fact we need to understand what I've written
before.

An image stack is a static entity and as such there is no causality,
except in its semantic interpretation. That's what you like to do but it
doesn't help much with mathematics.

"From the image acquisition perspective your view has merit."

I don't agree, because it is just the acquisition process that is
temporal. However, if it is finished and you've stored the data, then it
is static.

Regards

Herbie

::::::::::::::::::::::::::::::::::::::::::::
Am 01.11.18 um 20:30 schrieb [hidden email]:

> Herbie...
>
> There are more things in heaven and earth, Horatio,
> Than are dreamt of in your philosophy.
> - Hamlet (1.5.167-8), Hamlet to Horatio
>
> Depending on the time frame that you are concerned with will dictate if you
> look at your data as temporal.  From the image acquisition perspective your
> view has merit.  Although, from the study of what is contained in these
> sequential volumes of images, e.g., data points, they may very well be
> temporal.
>
> In MRI:
>
> fMRI and perfusion are real time data that monitors physiological processes as
> they occur, and, are analyzed as temporal signals.
>
> MRE(lastography) are real time data in which the subject of the images is
> perturbed in real time and the response, i.e., the images, are analyzed as
> temporal responses.
>
> B0 map, T1 map, T2 map, the independent variable is time, but not real time,
> albeit, temporal analysis tools are applicable.
>
> Diffusion data is varied by magnetic gradient field strength and timing
> parameters, temporal analysis tools are applicable.
>
> I could go on, but suffice it to say that ImageJ is very useful to examine
> this data in the frame dimension, and there is a lot on untapped functionality
> that can be provided here...
>
> Thanks for listening,
>
> Fred
>
>
> On Thu, November 1, 2018 11:50 am, Herbie wrote:
>> Jacob,
>>
>> I didn't want to imply this:
>>
>> "It seems that you are saying that ImageJ is limited in scope to static
>> images"
>>
>> What I had in mind are real-time signals.
>>
>> "Have I misunderstood you?"
>>
>> Yes.
>>
>> ImageJ isn't capable of dealing with real-time temporal signals.
>>
>> If you have a temporal sequence of images available, the sequence is no
>> longer temporal, it is a static stack of images and all the problems
>> occurring with causal temporal processing are irrelevant because the
>> stack of images is a-causal from an signal processing point of view. Of
>> course you know that e.g. slice 2 comes before slice 3 but it is no
>> problem to have e.g. a convolution kernel that extends from slice 1 to
>> slice 5. In real-time processing this isn't possible as long as slice 5
>> is available...
>>
>> I hope this makes things a bit clearer.
>>
>> Regards
>>
>> Herbie
>>
>> ::::::::::::::::::::::::::::::::::::::::::
>> Am 01.11.18 um 16:40 schrieb Jacob Keller:
>>>      Because here and with ImageJ we deal with images,
>>>      filter-characteristics
>>>      that are typical for the processing of causal time-signals don't play a
>>>      role and I don't think that they are or will be implemented for use
>>>      with
>>>      ImageJ.
>>>
>>>
>>> It seems that you are saying that ImageJ is limited in scope to static
>>> images? Huh? I see a lot of folks using it for time series, and there
>>> are quite a few plugins that incorporate time. I have noticed generally
>>> a bit of a split in the community regarding time series, but it seems to
>>> me that time is a critical component of biological processes, and at the
>>> level of microscopy is quite accessible, so why limit imageJ to static
>>> images? Have I misunderstood you?
>>>
>>> Jacob
>>>
>>>
>>
>> --
>> ImageJ mailing list: http://imagej.nih.gov/ij/list.html
>>
>
>

--
ImageJ mailing list: http://imagej.nih.gov/ij/list.html
Reply | Threaded
Open this post in threaded view
|

Re: Frequency Filtering in Time Dimension?

Jacob Keller-2
In reply to this post by Michael Schmid-3
If you mainly want to suppress *one* frequency, you can use moving
> averages, with a kernel length equal to the period you want to subtract
> (i.e., a rectangular kernel; this gives you no oscillations). To do it
> in time (slice) direction, reslice the stack and use the 'Fast Filters'
> plugin. The kernel length is 2*radius+1. Maybe apply some weak Gaussian
> filter in addition.
>

My data are responses to a given stimulus frequency, so I want to extract
responses at just that frequency. I have actually been doing just what you
suggest with a plugin from Jay Unruh which subtracts a moving window
average--I set the window to the period, output both the moving window and
the subtracted stacks. It's not perfect, however, since the moving window
introduces its own artifacts, as one can see from Bode plots of moving
windows. I guess the ideal might be simply to fit the data with the right
equations/parameters, then plot those parameters over time. But that's not
as agnostic as a filter.

If you want a more square-like filter in the Fourier domain (and you

> don't care about ringing), and the desired cutoff frequency is rather
> high, you can use Process>Filter>Convolve on the resliced stack. For a
> perfectly sharp frequency cutoff you would need a sin(a*x)/(a*x) kernel
> (the so-called "sinc" function); but this would give you infinite kernel
> size.
>    https://en.wikipedia.org/wiki/Sinc_function
> Multiply it with a Gaussian to limit the kernel size; this makes the
> cutoff a bit smoother (in frequency space, it gets convolved, i.e.,
> blurred, by the Fourier transform of the Gaussian that you have used for
> multiplying). The narrower the Gaussian you multiply with the sync, the
> smaller a kernel you can use (if you neglect very small numbers in the
> periphery of the kernel), but the frequency cutoff becomes smoother.
> Maybe a Gaussian with a sigma between the first and second zero is a
> good starting point (i.e., sigma between a*pi and a*2pi).
>

I was thinking of trying to use Lanczos filtering, which is two sinc's
superimposed. It's supposed to be a good compromise.

This is just an other way to express what Herbie has talked in his
> recent mails as "filtering in (real) time" vs. "filtering an image
> stack".
>

Okay, I get what Herbie was saying now. It would seem that the static case
would be much easier, then. I am thinking that fitting with some restraints
might be the best way to do this.

All the best,

Jacob



>
> Michael
> ________________________________________________________________
>
>
> On 2018-10-31 17:32, Jacob Keller wrote:
> > Thanks Michael--this is really clear and helpful! I checked out some
> > Bode
> > plots of various filters, and it seems Gaussian is pretty non-brick
> > wall. I
> > wonder whether a Butterworth or Bessel filter has been implemented?
> >
> > JPK
> >
> > On Wed, Oct 31, 2018 at 7:36 AM Michael Schmid
> > <[hidden email]>
> > wrote:
> >
> >> Hi Jacob,
> >>
> >> when you are doing a Gaussian Blur, in the Fourier domain it
> >> corresponds
> >> to multiplying the amplitudes with a Gaussian as well.
> >>
> >> A Gaussian
> >>     exp(-pi * x²/a²)
> >> in real space corresponds to
> >>     exp(-pi * f² * a²)
> >> in Fourier space
> >>
> >> If you want to attenuate a frequency component f1 by a factor 1/e =
> >> 0.37, you thus need a² = 1/(pi f1²), i.e. exp (pi² * f1² * x²)
> >>
> >> In terms of sigma, a Gaussian is given by
> >>    exp(-x²/(2*sigma²))
> >>
> >> Thus, attenuation by a factor of 1/e at f1 corresponds to
> >>    sigma = 1/(pi * f1 * sqrt(2))
> >>
> >>
> >> You can easily check it with the following macro:
> >>
> >> period = 16;  // between about 10 and 100
> >>                // best accuracy with powers of 2
> >> f=1/period;   // spatial frequency
> >> newImage("Untitled", "32-bit white", 256, 256, 1);
> >> run("Macro...", "code=v=sin(2*PI*x*"+f+")");  //create sine wave
> >> sigma = 1/(PI*f*sqrt(2));
> >> run("Gaussian Blur...", "sigma="+sigma);
> >> run("Set Measurements...", "min display redirect=None decimal=4");
> >> makeRectangle(25, 0, 200, 256);
> >> run("Measure");
> >>
> >>
> >> You will see that the maxima and minima of the sine wave (which has
> >> had
> >> an original amplitude of 1) gets attenuated to an amplitude of about
> >> 0.37.
> >>
> >> If you take sigma = 1/(pi * f1) you will get an attenuation factor of
> >> 1/e² = 0.13 at the frequency f1.
> >>
> >> So far an excursion into the math of Gaussians and their Fourier
> >> transform...
> >>
> >>
> >> Michael
> >> ________________________________________________________________
> >> On 30.10.18 16:23, Jacob Keller wrote:
> >>  >>
> >>  >> do you have simply a 1D data set or an image stack, with slices
> >> for
> >>  >> different times?
> >>  >>
> >>  >
> >>  > The latter--image stack, one plane over time.
> >>  >
> >>  > - In the second case, you can use Process>Filters>Gaussian Blur 3D
> >> and
> >>  >> specify the x and y sigma as zero.
> >>  >> For high-pass filtering, you would have to duplicate the stack
> >> first
> >> and
> >>  >> then subtract the filtered image.
> >>  >>
> >>  >
> >>  > Ah, this is a great idea--I will try it out. How can I figure out
> >> the
> >>  > relationship between sigma and the desired frequency cutoff? For
> >> example,
> >>  > let's say the stack has 120 frames per period--what sigma value
> >> should
> >> be
> >>  > input?
> >>  >
> >>  > Thanks very much for these suggestions--I think they are going to
> >> be
> >>  > very helpful,
> >>  >
> >>  > Jacob
> >>  >
> >>  > If you rather want moving averages or a median, you can reslice the
> >>  >> stack (Image>Stacks>Reslice) to have the time direction in x or y
> >> and
> >>  >> apply the 'Fast Filters' plugin, which can do 1D filtering. Then
> >> Reslice
> >>  >> to make time the z axis again.
> >>  >> The 'Fast Filters' plugin also has an option to subtract the
> >> filtered
> >>  >> image (i.e., highpass operation)
> >>  >>
> >>  >>
> >>  >> Michael
> >>  >> ________________________________________________________________
> >>  >> On 29.10.18 20:51, Jacob Keller wrote:
> >>  >>> Hi All,
> >>  >>>
> >>  >>> is there a way to do high-pass/low-pass filtering in the time
> >> dimension?
> >>  >> Or
> >>  >>> a suggestion for a workaround?
> >>  >>>
> >>  >>> All the best,
> >>  >>>
> >>  >>> Jacob Keller
> >>  >>
> >>  >> --
> >>  >> ImageJ mailing list: http://imagej.nih.gov/ij/list.html
> >>  >>
> >>  >
> >>  > --
> >>  > ImageJ mailing list: http://imagej.nih.gov/ij/list.html
> >>  >
> >>
> >> --
> >> ImageJ mailing list: http://imagej.nih.gov/ij/list.html
> >>
> >
> > --
> > ImageJ mailing list: http://imagej.nih.gov/ij/list.html
>
> --
> ImageJ mailing list: http://imagej.nih.gov/ij/list.html
>

--
ImageJ mailing list: http://imagej.nih.gov/ij/list.html