I have a bit of a mathematics problem.

I have an image $I_1$ that I have previously blurred (from an image $I_0$) with a Gaussian of standard deviation $\sigma_1$.

I would like to blur image $I_1$ with a Gaussian of standard deviation $\sigma_X$ to produce a new image, image $I_2$. The result, image $I_2$, should be as if $I_0$ has been blurred with a Gaussian of standard deviation $\sigma_2$.

This problem is not difficult, we know that the convolution of two Gaussian filters is the third filter. I can compute the $\sigma_X$ by: $$ \sigma_X = \sqrt{ \sigma_2^2 - \sigma_1^2 } $$ (this works because application of a Gaussian of $\sigma_1$ followed by application of Gaussian of $\sigma_2$, is as if a single Gaussian of $\sigma_3$ had been applied, where $\sigma_3^2 = \sigma_1^2 + \sigma_2^2$).

I can directly apply a Gaussian with $\sigma_X$ to image $I_1$ to achieve my desired image $I_2$.

Furthermore, I would like to decimate the image $I_2$ by a factor of $N$ (e.g., $0.5$ is decimation by half, $0.25$ is by fourth) after the blur. Keep in mind that neither of $1/N$ nor $\log_2(\sigma_2/\sigma_1)$ must be whole numbers, though they will always be positive.

Here is the problem. In the case where $\sigma_X$ is very large, the process becomes computationally expensive. Since I will be decimating the image anyways, I would like to apply the method often used in (dyadic) Gaussian pyramids, where one applies a blur followed by a decimation step, and then repeats this. However, since $\log_2(\sigma_2/\sigma_1)$ and $1/N$ are not guaranteed to be whole numbers, it is not straightforward to break the transform down into smaller blur-decimate steps so that I never have to apply a very large Gaussian convolution.

I have a feeling that this is possible.

Currently, I have:

    double sigmaratio = sigma2 / sigma1;
    double xsigma = sqrt( sigma2*sigma2 - sigma1*sigma1);
    double N = decimation_ratio //(e.g. 0.5 means decimate by half, 0.25 is by fourth, etc.)
    double sig_halves = log2( sigmaratio );
    double N_halves = -1 * log2(N); //will be =1 if we halve it once, =2 if twice, etc.

    //Taking floor to get the whole part...
    int intsighalves = (int) sig_halves;
    int intNhalves = (int) N_halves;
    int num = intsighalves;
    if( intNhalves < num )
       num = intNhalves;
    if( num >= 1 )
       double Ndiff = N_halves - (double)num;
       double sigdiff = sig_halves - (double)num;
       //I need to apply num dyadic steps with some sigmaXprime at each step (getting 2x as large each time) and decimation by half at each step.

My question is: how do I choose the sigmaXprime at each dyadic step for the whole parts, and when/how do I choose the partial sigmaXprimepartial and the partial decimation in order to get an exact implementation of the image $I_1$ convolved with Gaussian $\sigma_X$ followed by decimation by $N$?

Edit: Just a note, this may not be possible due to loss of information at each decimation step, but I think it should be possible since in the case of dyadic Gaussian pyramids, application of a single Gaussian / decimation step with the correct parameters is exactly equal to application of iterated Gaussian / decimation steps (see my answer to https://dsp.stackexchange.com/questions/667/image-pyramid-without-decimation/55654#55654 )

  • 1
    $\begingroup$ It is not clear what you want. Could you please try explaining it clearer? $\endgroup$ – Royi Mar 2 at 18:14
  • $\begingroup$ Thank you for the feedback, I will try. I need to figure out how to write latex here or else it will never be clear I think. $\endgroup$ – rveale Mar 2 at 18:21
  • $\begingroup$ For Math just write as you'd write in LaTeX: For In Line Math: $ <In Line Math Code> $, Display Math: $$ <Display Math Code> $$. $\endgroup$ – Royi Mar 2 at 20:58
  • $\begingroup$ Use FFT for the convolution. Decimation is by an integer factor. Resampling is changing the sample frequency, this is achieved by a convolution with $\frac{\sin(n)}{n}$ similar to that for the upsampling. Since we want non-circular convolution we need zero padding, and $\frac{\sin x}{x}$ is replaced by $\frac{\sin x}{x} w(x)$ for some window $w$. $\endgroup$ – reuns Mar 4 at 20:26

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.