Gaussian approximation of a rectangular pulse

Posted on by

Some image processing methods inherently smooth the images with an ideal filter, that is, a rectangular pulse. Comparing those methods with others that use a more typical Gaussian kernel can thus be challenging, as even if the full-width at half-maximum (FWHM) of the Gaussian filter is made equal to the width of the rectangular pulse, the resulting smoothnesses are not similar.

What would be the FWHM of a Gaussian filter that produces smoothness as similar as possible as that of a rectangular pulse of width ww? This must depend only on the choice of the standard deviation of the Gaussian filter, which is the only parameter we can choose, versus the width of the ideal filter, which is also the only parameter that can be controlled (when it even can). If you ever wondered about this important question, here you will find the answer.

Note: This is the first time I use AI for an article (Claude 4.6, specifically). It certainly saved some time. I believe I am in good company in that.

Definitions: Consider the following two real-valued functions each normalized to have unit area:

1) A rectangular pulse of width w>0w>0:

    r(x)={1w,|x|≤w2,0,|x|>w2r(x) = \begin{cases}\displaystyle\frac{1}{w},\qquad |x|\le \dfrac{w}{2},\\ 0, \qquad |x|\gt \dfrac{w}{2}\end{cases}

    2. A Gaussian centered at the origin with standard deviation σ>0\sigma > 0:

    g(x)=1σ2πexp⁡(−x22σ2)g(x) = \frac{1}{\sigma\sqrt{2\pi}}\exp\left(-\frac{x^{2}}{2\sigma^{2}}\right)

    We want to find the value of σ\sigma as a function of ww that minimizes the integral of squared difference between the two functions:

    D(σ)=∫−∞∞(g(x)−r(x))2dxD(\sigma) = \int_{-\infty}^{\infty}\left(g(x)-r(x)\right)^{2}dx

    Solution: We start by expanding the integrand into 3 separate terms:

    D(σ)=∫−∞∞g2(x)dx⏟I1−2∫−∞∞g(x)r(x)dx⏟I2+∫−∞∞r2(x)dx⏟I3.D(\sigma) = \underbrace{\int_{-\infty}^{\infty} g^{2}(x)dx}_{I_{1}} -2\underbrace{\int_{-\infty}^{\infty} g(x)r(x)dx}_{I_{2}} +\underbrace{\int_{-\infty}^{\infty} r^{2}(x)dx}_{I_{3}}.

    For I1I_1, we substitute u=x/σu=x/\sigma so dx=σdu:dx=\sigma du:

    I1=∫−∞∞12πσ2e−x2/σ2dx=12πσ2⋅σ∫−∞∞e−u2du=12πσ⋅π=12σπI_{1} = \int_{-\infty}^{\infty}\frac{1}{2\pi\sigma^{2}} e^{-x^{2}/\sigma^{2}} dx = \frac{1}{2\pi\sigma^{2}}\cdot\sigma \int_{-\infty}^{\infty}e^{-u^{2}} du = \frac{1}{2\pi\sigma}\cdot\sqrt{\pi} = \boxed{\frac{1}{2\sigma\sqrt{\pi}}}

    For I3I_3, since r(x)=1/wr(x)=1/w on an interval of length ww and zero elsewhere, we have:

    I3=∫−w/2w/21w2dx=ww2=1wI_{3} = \int_{-w/2}^{w/2}\frac{1}{w^{2}} dx = \frac{w}{w^{2}} = \boxed{\frac{1}{w}}

    which does not depend on σ\sigma.

    For I2I_2, given that r(x)=1/wr(x)=1/w only on the interval [−w/2,w/2][-w/2, w/2], we have:

    I2=1w∫−w/2w/21σ2πe−x2/(2σ2)dxI_{2} = \frac{1}{w}\int_{-w/2}^{w/2} \frac{1}{\sigma\sqrt{2\pi}} e^{-x^{2}/(2\sigma^{2})} dx

    Substituting u=x/(2σ)u = x/(\sqrt{2} \sigma) so that dx=2σdudx=\sqrt{2} \sigma du, and defining β≡w22σ\beta \equiv \frac{w}{2\sqrt{2} \sigma}, thus mapping the integration limits to ±β\pm\beta:

    I2=1wσ2π⋅2σ∫−ββe−u2du=1wπ∫−ββe−u2duI_{2} = \frac{1}{w \sigma\sqrt{2\pi}}\cdot\sqrt{2} \sigma \int_{-\beta}^{\beta}e^{-u^{2}} du = \frac{1}{w\sqrt{\pi}} \int_{-\beta}^{\beta}e^{-u^{2}} du

    Noting that ∫−ββe−u2du=πerf⁡(β)\displaystyle\int_{-\beta}^{\beta}e^{-u^{2}} du = \sqrt{\pi} \operatorname{erf}(\beta), we obtain:

    I2=erf⁡(β)wI_{2} = \frac{\operatorname{erf}(\beta)}{w}

    Hence:

    D=12σπ−2erf⁡(β)w+1wD = \frac{1}{2\sigma\sqrt{\pi}} – \frac{2 \operatorname{erf}(\beta)}{w} + \frac{1}{w}

    We need now to differentiate DD with respect to σ\sigma and set it to zero. The derivative of the first term is:

    ddσ(12σπ)=−12σ2π\frac{d}{d\sigma}\left(\frac{1}{2\sigma\sqrt{\pi}}\right) = – \frac{1}{2\sigma^{2}\sqrt{\pi}}

    For the second term, using the chain rule with ddβerf⁡(β)=2πe−β2\displaystyle\frac{d}{d\beta}\operatorname{erf}(\beta)=\frac{2}{\sqrt{\pi}} e^{-\beta^{2}} and dβdσ=−w22σ2\displaystyle\frac{d\beta}{d\sigma}=-\frac{w}{2\sqrt{2} \sigma^{2}}, we have:

    ddσ[−2erf⁡(β)w]=−2w⋅2πe−β2⋅(−w22σ2)=2e−β22πσ2=2e−β2σ2π\frac{d}{d\sigma}\left[\frac{-2 \operatorname{erf}(\beta)}{w}\right] = \frac{-2}{w}\cdot\frac{2}{\sqrt{\pi}} e^{-\beta^{2}} \cdot\left(-\frac{w}{2\sqrt{2} \sigma^{2}}\right) = \frac{2 e^{-\beta^{2}}}{\sqrt{2\pi} \sigma^{2}} = \frac{\sqrt{2} e^{-\beta^{2}}}{\sigma^{2}\sqrt{\pi}}

    Setting dD/dσ=0dD/d\sigma=0:

    −12σ2π+2e−β2σ2π=0– \frac{1}{2\sigma^{2}\sqrt{\pi}} +\frac{\sqrt{2} e^{-\beta^{2}}}{\sigma^{2}\sqrt{\pi}} = 0

    and solving for β\beta:

    e−β2=122e^{-\beta^{2}} = \frac{1}{2\sqrt{2}}
    β2=ln⁡(22)=ln⁡2+12ln⁡2=32ln⁡2\beta^{2} = \ln\bigl(2\sqrt{2}\bigr) = \ln 2 + \tfrac{1}{2}\ln 2 = \frac{3}{2}\ln 2
    β=32ln⁡2\beta = \sqrt{\tfrac{3}{2}\ln 2}

    Recalling we defined β=w/(22σ)\beta = w/(2\sqrt{2} \sigma) and solving for σ\sigma:

    σ=w12ln⁡2\sigma = \frac{w}{\sqrt{12\ln 2}}

    From a previous article, we know that the FWHM of a Gaussian filter can be expressed as σ8ln⁡2\sigma \sqrt{8 \ln 2}.

    Final result: Hence, we have:

    FWHM=23w\text{FWHM} = \sqrt{\frac{2}{3}}w

    Which means that the FWHM of the Gaussian filter is approximately 0.8164965809277 times the width of the rectangular filter.

    Example: Suppose you did some image processing that involved averaging over a sphere with radius 30 mm. The width ww of the rectangular function is 60 mm. The FWHM of the closest Gaussian filter that yields similar degree of smoothing is approximately 48.99 mm.

    Leave a comment