Decimation vs. Collapsing and their spectral effects

Thanks for the expansion Martin. Collapsing in 1 dimension is equivalent to
first convolving your sequence or scan with a two tap 'boxcar' filter then
throwing away every second sample. You get some reduction in the high
frequency content of the scan because of the convolution so you can hope
that , if there was no very large high frequency component, then the low
pass effect will be enough to stop the aliased components from messing up
your new, low resolution image. When you do it time and time again the
aliased components must start to build up. Do you have a set of criteria for
deciding what is significant in your image and what constitutes a
noticeable/agravating distortion?

Best of luck - Mike

"Martin J. Stumpf" <mjs@neuron.rad.jhmi.edu> wrote in message
news:slrncl0p59.1hs.mjs@neuron.rad.jhmi.edu...
> Thanks for your response Mike,
>
> I am collapseing the samples on purpose to double the pixel size of
> the image.  This is necessary because I am looking at the effects of
> different projection bin sizes on iterative image reconstruction from
> projections. I know what you are saying in your first response and
> appreciate you saying "What are you worried about?". I am worried that
> when we collapse an image by two say to change 1mm pixels into 2mm
> pixels that we are also changing the frequencies in that image.
> ie. The noise will be attenuated just by virtue of the collapsing. My
> goal is to be able to describe or predict the outcome of changing the
> pixel sizes using DSP theory. Another more direct way of asking the
> question is this: Is collapsing in this manner considered resampling?
> And if it is can I predict the effect on the resulting power spectrum?
> I guess what I'm 'worried about' is whether or not I am changing the
> spatial resolution in the image by collapsing. I guess in the time
> domain the analogy would be that by adding together every two samples
> one would be doubling the time per sample and therefore affecting the
> time resolution.
>
> Thank you very much for your time
>
> -Martin
>
> In article <cipmaf$c8c$1@hercules.btinternet.com>, Mike Yarwood wrote:
> > Another way of looking at it is as a 1/2 sample shift due to using a
> > trapezoidal interpolator then decimation by two.  Why don't you just
throw
> > every other sample away instead?
> >
> > Best of Luck - Mike
> >
> > "Mike Yarwood" <mpyarwood@btopenworld.com> wrote in message
> > news:cipll6$ru2$1@titan.btinternet.com...
> >> Hello Martin :  it looks as though you are redefining the spatial
sampling
> >> interval and doing some sort of averaging.   If the original space
between
> >> your samples was dphi your new samples are spaced 2*dphi and offset
from
> > the
> >> original samples by dphi/2  (either +/- depending on which edge you
start
> > to
> >> collapse from).  By adding them together it's just like taking the
> >> averaqge*2.   You must have already had these thoughts though so what
are
> >> you really worried about?
> >>
> >> Mike
> >>
> >> "christophe grimault" <christophe.grimault@novagrid.com> wrote in
message
> >> news:4150497D.5040402@novagrid.com...
> >> >
> >> >
> >> > Martin J. Stumpf wrote:
> >> >
> >> > >Hello all,
> >> > >
> >> > >I am struggling with what happens to the time axis in a sampled
> >> > >sequence when it is collapsed by adding every two samples together,
> >> > >with or wihout averaging.
> >> > >I am working in the spatial domain but believe the principles should
> >> > >be the same as with time varying signals.
> >> > >
> >> > >I understand decimation and how it effectively resamples at a lower
> >> > >sampling frequency and has the potential for inducing aliasing if
the
> >> > >sequence is not first band limited. I understand this because it is
> >> > >obvious that in decimation you are just using a different sampling
> >> > >period and resampling. But I am not sure what actually takes place
in
> >> > >collapsing.
> >> > >
> >> > >What is tripping me up is that it seems like when adding two samples
> >> > >together the sample is not instantaneous but rather the integral
over
> >> > >the duration of the two samples used. I have been working with
matlab
> >> > >to demonstrate to myself what happens to a sampled sequence with
> >> > >repeated interpolations and collapsings. When I collapse a sequence
> >> > >say of N=64 to N=32 I have combined each pair of the original
samples
> >> > >and thus must have doubled the time/space of each sample. The
collapse
> >> > >acts like a lowpass filter and attenuates the high frequencies but
it
> >> > >is not directly analagous I don't think because the bin width has
been
> >> > >doubled. With a FIR filter the bin width doesn't change.
> >> > >
> >> > >It is difficult for me to frame the actual question and I apologize
> >> > >for that. I want to be able to understand the theoretical
> >> > >underpinnings of just what happens during a collapse. We do this all
> >> > >the time in medical imaging and I want to be able to explain it with
> >> > >DSP.
> >> > >
> >> > >Any insights greatly appreciated.
> >> > >
> >> > >-Martin
> >> > >
> >> > >
> >> >
> >>
> >>
> >
> >

Jerry Avins wrote:
> 
> jim wrote:
> 
>    ...
> 
> > Hi Martin
> >       The frequency response of your filter [1,1] is just cos function.
> 
>    ...
> 
> That seems very strange to me; I would have thought it is a sinc. Can
> you explain?

e^(i(x+a)) + e^(i(x-a)) = 2*cos(a)*e^ix

The frequency response is a sinc (sin(x)/x) only in the limit if you
have an infinite number of ones in your filter. For a finite length
filter its what's often called an aliased or periodic sinc. A periodic
sinc can also be expressed as a summation of geometric progression of
cosines (you need only the above equation to derive that fact). In this
particular case the summation only has one term.

-jim


-----= Posted via Newsfeeds.Com, Uncensored Usenet News =-----
http://www.newsfeeds.com - The #1 Newsgroup Service in the World!
-----==  Over 100,000 Newsgroups - 19 Different Servers! =-----

Bernhard Holzmayer wrote:

> > So it's a very simple, more-or-less brute force 2:1 decimating
> > filter. It lacks some elegance and isn't very high performance,
> > but it's awfully simple to compute.
> 
> Can you imagine a better version of anti-aliasing filter when
> decimation ratio is variable as N(t):1 ?
> 

Your filter is essentially sampling a rectangular window. To get a
better filter simply sample a better window. For instance sampling a
gaussian would be a simple way to make a better anti-aliasing filter.
You could refer to pascal's triangle to get your coefficients.

	But I caution you again to find out what you might be filtering. From
your previous statements there is no reason, that I can see, not to
believe that it is just as likely the signal might vary with speed and
that the noise might have a linear relationship to the speed of the
process.

-jim


> >
> > I'd opine that this isn't a very good thing to do from a DSP
> > perspective,
> So do I, and indeed, I have a very bad feeling, especially,
> because I'm doing it in a measurement application where aliasing
> might be a critical effect.
> 
> > but since what you're really doing is image
> > processing it might be suitable for whatever it is that you're
> > trying to do (which I don't completely understand).
> >
> 
> > Hope that helps a little, anyway.
> 
> Helps me to see, that others deal with similar issues.
> Gives hope, that there's a better approach already...
> 
> Bernhard


-----= Posted via Newsfeeds.Com, Uncensored Usenet News =-----
http://www.newsfeeds.com - The #1 Newsgroup Service in the World!
-----==  Over 100,000 Newsgroups - 19 Different Servers! =-----

Bernhard Holzmayer wrote:

   ...

> Can you imagine a better version of anti-aliasing filter when 
> decimation ratio is variable as N(t):1 ?

There was a time when I thought I understood your problem. It's clear
now that I don't. Why do you need to decimate?

Jerry
-- 
Engineering is the art of making what you want from things you can get.
�����������������������������������������������������������������������

jim wrote:

   ...

> In this particular case the summation only has one term.
> 
> -jim

Neat!. Thanks.

Jerry
-- 
Engineering is the art of making what you want from things you can get.
�����������������������������������������������������������������������

On Wed, 22 Sep 2004 08:23:30 -0500, jim <"N0sp"@m.sjedging@mwt.net>
wrote:

>
>
>Jerry Avins wrote:
>> 
>> jim wrote:
>> 
>>    ...
>> 
>> > Hi Martin
>> >       The frequency response of your filter [1,1] is just cos function.
>> 
>>    ...
>> 
>> That seems very strange to me; I would have thought it is a sinc. Can
>> you explain?
>
>e^(i(x+a)) + e^(i(x-a)) = 2*cos(a)*e^ix
>
>The frequency response is a sinc (sin(x)/x) only in the limit if you
>have an infinite number of ones in your filter. For a finite length
>filter its what's often called an aliased or periodic sinc. A periodic
>sinc can also be expressed as a summation of geometric progression of
>cosines (you need only the above equation to derive that fact). In this
>particular case the summation only has one term.
>
>-jim

Hi Jim,

   I never realized that the absolute value of the 
DFT of [1,1] was a cosine function.  I think you're 
right, although I don't know from where your 
equation came.  (If "x" represents frequency, 
what is "a"?).

You made me realize that a system with an 
impulse response of [1,1] is a special kind of 
delay-comb filter and I happen to know that 
those filters have a freq magnitude response 
that's sinusoidal in shape with a peak value 
of 2.

So I believe the magnitude of the freq response 
of a two-tap (unity-valued coefficients) is:

  |X(w)| = 2cos(w/2)

where w is in the range 0 -to- 2pi, 
and represents frequency.  

The discrete time Fourier transform of the 
sequence [1,1] is:

  X(w) = 1 + cos(w) -jsin(w).

What I'm unable to do is prove to myself 
(algebraically) that the magnitude of 

     1 + cos(w) -jsin(w)

equals 2cos(w/2).

By the way, as far as I've seen, the DFT of 
a finite number of ones is called a 
"Dirichlet kernel".

See Ya,
[-Rick-]

Jerry Avins wrote:

> Bernhard Holzmayer wrote:
> 
>    ...
> 
>> Can you imagine a better version of anti-aliasing filter when
>> decimation ratio is variable as N(t):1 ?
> 
> There was a time when I thought I understood your problem. It's
> clear now that I don't. Why do you need to decimate?
> 
> Jerry

My samples have been taken at a rate of 48kS/s (A), at an earlier 
stage.
Then I resample the signal at a lower rate which is variable,
in fact between 0...48kS/s, let's assume 32S/s...8kS/s for now 
(to avoid other issues), rate B.
Although my variable rate B is special, the concept seems to be  
that of a decimation process.

First approach:  
at every clock of rate B I pick the current sample out of rate A.
Works well, but lots of aliasing...

Second approach:
at every clock of rate B I pick the average of the samples at rate A 
since the last take at rate B. That's the filter with h=1,1,... 
coeffs. Works better, since less aliasing effects.

Guess:
the filter cannot be optimal, because the averaging isn't done 
continuously, but in "chunks". The steps on the borders are 
suboptimal. 
Therefore my query for a better "decimation filter".

Maybe, my current solution is already good enough, and cannot be 
improved much without big effort...

Bernhard
 

Rick Lyons wrote:

   ...

> What I'm unable to do is prove to myself 
> (algebraically) that the magnitude of 
> 
>      1 + cos(w) -jsin(w)
> 
> equals 2cos(w/2).

Let M = magnitude
M� = [1 + cos(w) - jsin(w)]� � [1 + cos(w) + jsin(w)]�
= [1 + cos(w)]� + sin�(w)
= 1 + 2cos(w) + cos�(w) + sin�(w)
= 2 + 2cos(w) = 2[1+ cos(w)]
Now we need the square root. Let's cheat, and square 2cos(w/2).
Do you see it now?

Jerry
-- 
... they proceeded on the sound principle that the magnitude of a lie
always contains a certain factor of credibility, ... and that therefor
... they more easily fall victim to a big lie than to a little one ...
                                                                   A. H.

�����������������������������������������������������������������������

On Thu, 23 Sep 2004 10:23:43 -0400, Jerry Avins <jya@ieee.org> wrote:

>Rick Lyons wrote:
>
>   ...
>
>> What I'm unable to do is prove to myself 
>> (algebraically) that the magnitude of 
>> 
>>      1 + cos(w) -jsin(w)
>> 
>> equals 2cos(w/2).

Hi Jerry,

>Let M = magnitude
>M&#4294967295; = [1 + cos(w) - jsin(w)]&#4294967295; &#4294967295; [1 + cos(w) + jsin(w)]&#4294967295;

Humm, I don't know about the above, 
but I like the following.

>= [1 + cos(w)]&#4294967295; + sin&#4294967295;(w)
>= 1 + 2cos(w) + cos&#4294967295;(w) + sin&#4294967295;(w)
>= 2 + 2cos(w) = 2[1+ cos(w)]

Yes, yes.

Next we say: M&#4294967295; = 4[1+cos(w)]/2.

Then we use the "half-angle" trig identity 
of: 

    sqrt([1+cos(w)]/2) = cos(w/2), 

to state

   sqrt(4[1+cos(w)]/2) = sqrt(4)cos(w/2) = 2cos(w/2).

Good Jerry!
---
Aarrggh!  I see another way:

The complex freq response of

   1 + cos(w) -jsin(w) 

can be written

   1 + exp(-jw)
 
     = exp(-jw/2)[exp(jw/2) + exp(-jw/2)].

Those terms inside the brackets are equal to

     2cos(w/2).  

So the magnitude of the complex freq response is

    |exp(-jw/2)[2cos(w/2)]|  =  2cos(w/2).

Sheece, how simple!
  
(The above is surely what Jim was saying when he 
wrote: "e^(i(x+a)) + e^(i(x-a)) = 2*cos(a)*e^ix".)

*** NOW ***

Someone might read our thread and say: "Why the heck 
are Jerry and Rick screwing around with this 
silly simple two-coefficient filter's magnitude 
response?"

The reason is: the above sinusoidal shaped magnitude 
response *also* applies to a transversal FIR filter 
whose coefficients are

    h1(k) = 1,0,0,0,0,0,0,1.

This h1(k) is an example of a comb filter used in 
some "frequency sampling filters". (Section 7.1 of 
my book).

In addition, that sinusoidal shaped magnitude response
*also* applies to a transversal FIR filter whose 
coefficients are

    h1(k) = 1,0,0,0,0,0,0,-1

which is an example of the comb portion of a cascaded 
integrator-comb (CIC) filter used in many digital 
communications systems.

Neat huh?

[-Rick-]

Hi Rick
	My newsreader seems to have missed some posts. This is the first one I
see since my last post. (see more below)

Rick Lyons wrote:
> 
> On Thu, 23 Sep 2004 10:23:43 -0400, Jerry Avins <jya@ieee.org> wrote:
> 
> >Rick Lyons wrote:

> 
> The complex freq response of
> 
>    1 + cos(w) -jsin(w)
> 
> can be written
> 
>    1 + exp(-jw)
> 
>      = exp(-jw/2)[exp(jw/2) + exp(-jw/2)].
> 
> Those terms inside the brackets are equal to
> 
>      2cos(w/2).
> 
> So the magnitude of the complex freq response is
> 
>     |exp(-jw/2)[2cos(w/2)]|  =  2cos(w/2).
> 
> Sheece, how simple!
> 
> (The above is surely what Jim was saying when he
> wrote: "e^(i(x+a)) + e^(i(x-a)) = 2*cos(a)*e^ix".)

Yes, glad you waded through the math (and didn't make me do it).

	Remember the OP's problem was with image processing and pixels, so
there's no need to be concerned with causality. Therefore one is free to
express the frequency response with reference to the center of the
filter rather than the end as you do (that's where your exp(-jw/2) comes
from). So it could be said like this: averaging 2 pixels produces the
cosine magnitude response at the point halfway in between.

> 
> *** NOW ***
> 
> Someone might read our thread and say: "Why the heck
> are Jerry and Rick screwing around with this
> silly simple two-coefficient filter's magnitude
> response?"
> 
> The reason is: the above sinusoidal shaped magnitude
> response *also* applies to a transversal FIR filter
> whose coefficients are
> 
>     h1(k) = 1,0,0,0,0,0,0,1.
> 
> This h1(k) is an example of a comb filter used in
> some "frequency sampling filters". (Section 7.1 of
> my book).
> 
> In addition, that sinusoidal shaped magnitude response
> *also* applies to a transversal FIR filter whose
> coefficients are
> 
>     h1(k) = 1,0,0,0,0,0,0,-1

Right, but in that case the response will be a sine function.

> 
> which is an example of the comb portion of a cascaded
> integrator-comb (CIC) filter used in many digital
> communications systems.
> 
> Neat huh?

Yes, and if you constrain your filters to purely symmetric or purely
anti-symmetric FIR (as people inevitably do in image processing), then
any such filter can be expressed as a simple sum of sines or cosines.
For example, an even length box car filter can be considered to be the
sum of the filters:

 [1,1] + [1,0,0,1] + [1,0,0,0,0,1] .....

And therefore the magnitude response is just the sum of the respective
cosine functions. 

-jim


-----= Posted via Newsfeeds.Com, Uncensored Usenet News =-----
http://www.newsfeeds.com - The #1 Newsgroup Service in the World!
-----==  Over 100,000 Newsgroups - 19 Different Servers! =-----