<p>I tried to follow Nathan Reed post <a href="http://www.reedbeta.com/blog/antialiasing-to-splat-or-not/" rel="nofollow noreferrer">Antialiasing: To Splat Or Not</a> to implement splat method. The image I am getting thorugh the Lanczos filter is very different than Nathan. I do not understand where the problem is.</p>

<p>According to <a href="https://en.wikipedia.org/wiki/Lanczos_resampling" rel="nofollow noreferrer">Wikipedia page</a> Lanczos filter is defined as either of below formulas:</p>

<p><a href="https://i.sstatic.net/UEPNb.png" rel="nofollow noreferrer"><img src="https://i.sstatic.net/UEPNb.png" alt="enter image description here"></a></p>

<p><a href="https://i.sstatic.net/oqq67.png" rel="nofollow noreferrer"><img src="https://i.sstatic.net/oqq67.png" alt="enter image description here"></a></p>

<p><a href="https://gist.github.com/Reedbeta/893b63390160e33ddb3c" rel="nofollow noreferrer">Nathan implementation</a> seems to have missed a division by pi**2 and a multiplication by filter radius.</p>

<pre><code>elif options.kernel == 'lanczos':
        weight = 
          (np.sin(np.pi * sampleOffsetX) * np.sin(np.pi * sampleOffsetY) *
           np.sin(np.pi * sampleOffsetX / options.radius) * 
           np.sin(np.pi * sampleOffsetY / options.radius) /                           
          (sampleOffsetX**2 * sampleOffsetY**2))
</code></pre>

<p>Here is my code which borrowed most of the filter implementation from pbrt code:</p>

<pre><code>private void btn_Render_Click(object sender, EventArgs e)
{
    Random rand = new Random();
    int samples = 16;
    ImageFilm film = new ImageFilmSplat(pictureBox1.Width, pictureBox1.Height);
    film.Filter = new SincFilter(tau: 4, width: 4);

    for (int i = 0; i &lt; film.Width; i++)
    {
        for (int j = 0; j &lt; film.Height; j++)
        {
            EvaluatePixel(i, j, samples, rand, film);
        }
    }

    for (int i = 0; i &lt; film.Width; i++)
    {
        for (int j = 0; j &lt; film.Height; j++)
        {
            Pixel px = film.GetPixel(i, j);

            if (px.filterWeightSum != 0)
            {
                px.radiance = px.radiance / px.filterWeightSum;

                film.SetPixel(i, j, px);
            }
        }
    }

    pictureBox1.Image = film.GetImage();
}

void EvaluatePixel(int i, int j, int nSamples, Random rand, ImageFilm film)
{
    Vector halfpixel = new Vector(.5 / film.Width, .5 / film.Height, 0);

    //all samples between [0,1]
    Sample[] samples = sampler.Random(rand, nSamples);

    for (int n = 0; n &lt; nSamples; n++)
    {
        //set sample coord at centre of the pixel
        double cx = ((double)i / film.Width) + halfpixel.X;
        double cy = ((double)j / film.Height) + halfpixel.Y;

        //jitter the centre
        cx += ((samples[n].X * 2) - 1) * halfpixel.X
        cy += ((samples[n].Y * 2) - 1) * halfpixel.Y

        //sample the function
        double f = func(cx, cy);

        film.AddToPixel(i, j, new Vector(f));
    }
}

//ImageFilm class
public void AddToPixel(int i, int j, Vector radiance)
{ 
    //splat the value to neighbor pixels
    //find neighbour pixels based on filter width
    int x0 = (int)Math.Max(0, i - Filter.Width);
    int x1 = (int)Math.Min(Width - 1, i + Filter.Width);
    int y0 = (int)Math.Max(0, j - Filter.Width);
    int y1 = (int)Math.Min(Height - 1, j + Filter.Width);

    for (int x = x0; x &lt;= x1; x++)
    {
        for (int y = y0; y &lt;= y1; y++)
        {
            double weight = this.Filter.eval(x - i, y - j);

            Pixel pxl = GetPixel(x, y);

            pxl.radiance += radiance * weight;
            pxl.filterWeightSum += weight;

            SetPixel(x, y, pxl);
        }
    }
}

public class SincFilter
{
    private double tau;

    public SincFilter(double tau1, double width)
    {
        Width = width;
        tau = tau1;
    }
    double Sinc(double x)
    {
        x = Math.Abs(x);
        if (x &lt; 1e-5) return 1;
        return Math.Sin(Math.PI * x) / (Math.PI * x);
    }
    double WindowedSinc(double x, double radius)
    {
        x = Math.Abs(x);
        if (x &gt; radius) return 0;
        return Sinc(x) * Sinc(x / tau);
    }
    public override double eval(double x, double y)
    {
        return WindowedSinc(x, Width) * WindowedSinc(y, Width);
    }
}

double func(double x, double y)
{
    double minPeriod = 2e-5;
    double maxPeriod = 0.2;

    double period = minPeriod + (maxPeriod - minPeriod) * (y * y);
    double phase = x / period;
    phase -= Math.Floor(phase);

    return Math.Round(phase);
}
</code></pre>

<p>Here is Nathan image : </p>

<p><a href="https://i.sstatic.net/R5VlB.png" rel="nofollow noreferrer"><img src="https://i.sstatic.net/R5VlB.png" alt="enter image description here"></a></p>

<p>And here is mine for lanczos filter width 4, tau 4</p>

<p><a href="https://i.sstatic.net/e3thP.png" rel="nofollow noreferrer"><img src="https://i.sstatic.net/e3thP.png" alt="enter image description here"></a></p>

<p>My implementation is quite noisy for some reason. It doesn't look any better than box filter.</p>

CC BY-SA 4.0
2019-03-19T06:08:51.940Z
I notice in Nathan's version there is a linear->sRGB conversion which may account for his version looking brighter & flatter. I would suggest looking at stratifying the samples too  to lower graininess, unless your sampler.Random() already does that ?
PaulHK
CC BY-SA 4.0
2019-03-19T10:10:28.180Z
Yes I have tried stratified samples too which improves it but nothing like Nathan's. Also I wonder the second formula in Wiki page should be modified to "a**2" in the numerator instead of "a"?
ali
CC BY-SA 4.0
2019-03-19T10:19:02.987Z
Ideally a minimum of 256 samples should be used so that it can converge on the full 0-255 RGB output range accurately. As most of your samples are going to be either black or white in your input image, 16 average samples are going to give you a none-smooth / quantised range. Although your filter is going to negate most of that, I didn't study too closely how your filter is working. I still think the linear->sRGB is responsible for most of Nathans 'flatness'
PaulHK
CC BY-SA 4.0
2019-03-19T13:59:14.527Z
Thank Paul. Could explain a bit more how  linear->sRGB works? Do you refer to tone-mapping?
ali
CC BY-SA 4.0
2019-03-19T16:31:24.770Z
BTW, the factor $a/\pi^2$ in the weights is not needed because the final image is normalized by the sum of weights anyway, so any constant factor will cancel out.
Nathan Reed
CC BY-SA 4.0
2019-03-19T16:36:59.303Z
Re: sRGB, no it's not tone-mapping, it's a conventional nonlinear encoding of pixel values in images to obtain better precision in the dark values where our eyes are more sensitive. You can read here for more: https://en.wikipedia.org/wiki/SRGB#The_sRGB_transfer_function_(%22gamma%22)
Nathan Reed
CC BY-SA 4.0
2019-03-19T22:30:44.710Z
Paul, your guess was correct. After applying the conversion to sRGB it soomths out the image. Thanks.
ali
<image-processing><interpolation><filtering>
Lanczos filter implementation issue
ali
# Lanczos filter implementation issue

(Asked by: username_0 on 2019-03-18)

I tried to follow username_2 username_2 post [Antialiasing: To Splat Or Not](http://www.reedbeta.com/blog/antialiasing-to-splat-or-not/) to implement splat method. The image I am getting thorugh the Lanczos filter is very different than username_2. I do not understand where the problem is.

According to [Wikipedia page](https://en.wikipedia.org/wiki/Lanczos_resampling) Lanczos filter is defined as either of below formulas:

[![enter image description here](https://i.sstatic.net/UEPNb.png)](https://i.sstatic.net/UEPNb.png)

[![enter image description here](https://i.sstatic.net/oqq67.png)](https://i.sstatic.net/oqq67.png)

[username_2 implementation](https://gist.github.com/Reedbeta/893b63390160e33ddb3c) seems to have missed a division by pi\*\*2 and a multiplication by filter radius.

```
elif options.kernel == 'lanczos':
        weight = 
          (np.sin(np.pi * sampleOffsetX) * np.sin(np.pi * sampleOffsetY) *
           np.sin(np.pi * sampleOffsetX / options.radius) * 
           np.sin(np.pi * sampleOffsetY / options.radius) /                           
          (sampleOffsetX**2 * sampleOffsetY**2))

```

Here is my code which borrowed most of the filter implementation from pbrt code:

```
private void btn_Render_Click(object sender, EventArgs e)
{
    Random rand = new Random();
    int samples = 16;
    ImageFilm film = new ImageFilmSplat(pictureBox1.Width, pictureBox1.Height);
    film.Filter = new SincFilter(tau: 4, width: 4);

    for (int i = 0; i < film.Width; i++)
    {
        for (int j = 0; j < film.Height; j++)
        {
            EvaluatePixel(i, j, samples, rand, film);
        }
    }

    for (int i = 0; i < film.Width; i++)
    {
        for (int j = 0; j < film.Height; j++)
        {
            Pixel px = film.GetPixel(i, j);

            if (px.filterWeightSum != 0)
            {
                px.radiance = px.radiance / px.filterWeightSum;

                film.SetPixel(i, j, px);
            }
        }
    }

    pictureBox1.Image = film.GetImage();
}

void EvaluatePixel(int i, int j, int nSamples, Random rand, ImageFilm film)
{
    Vector halfpixel = new Vector(.5 / film.Width, .5 / film.Height, 0);

    //all samples between [0,1]
    Sample[] samples = sampler.Random(rand, nSamples);

    for (int n = 0; n < nSamples; n++)
    {
        //set sample coord at centre of the pixel
        double cx = ((double)i / film.Width) + halfpixel.X;
        double cy = ((double)j / film.Height) + halfpixel.Y;

        //jitter the centre
        cx += ((samples[n].X * 2) - 1) * halfpixel.X
        cy += ((samples[n].Y * 2) - 1) * halfpixel.Y

        //sample the function
        double f = func(cx, cy);

        film.AddToPixel(i, j, new Vector(f));
    }
}

//ImageFilm class
public void AddToPixel(int i, int j, Vector radiance)
{ 
    //splat the value to neighbor pixels
    //find neighbour pixels based on filter width
    int x0 = (int)Math.Max(0, i - Filter.Width);
    int x1 = (int)Math.Min(Width - 1, i + Filter.Width);
    int y0 = (int)Math.Max(0, j - Filter.Width);
    int y1 = (int)Math.Min(Height - 1, j + Filter.Width);

    for (int x = x0; x <= x1; x++)
    {
        for (int y = y0; y <= y1; y++)
        {
            double weight = this.Filter.eval(x - i, y - j);

            Pixel pxl = GetPixel(x, y);

            pxl.radiance += radiance * weight;
            pxl.filterWeightSum += weight;

            SetPixel(x, y, pxl);
        }
    }
}

public class SincFilter
{
    private double tau;

    public SincFilter(double tau1, double width)
    {
        Width = width;
        tau = tau1;
    }
    double Sinc(double x)
    {
        x = Math.Abs(x);
        if (x < 1e-5) return 1;
        return Math.Sin(Math.PI * x) / (Math.PI * x);
    }
    double WindowedSinc(double x, double radius)
    {
        x = Math.Abs(x);
        if (x > radius) return 0;
        return Sinc(x) * Sinc(x / tau);
    }
    public override double eval(double x, double y)
    {
        return WindowedSinc(x, Width) * WindowedSinc(y, Width);
    }
}

double func(double x, double y)
{
    double minPeriod = 2e-5;
    double maxPeriod = 0.2;

    double period = minPeriod + (maxPeriod - minPeriod) * (y * y);
    double phase = x / period;
    phase -= Math.Floor(phase);

    return Math.Round(phase);
}

```

Here is username_2 image :

[![enter image description here](https://i.sstatic.net/R5VlB.png)](https://i.sstatic.net/R5VlB.png)

And here is mine for lanczos filter width 4, tau 4

[![enter image description here](https://i.sstatic.net/e3thP.png)](https://i.sstatic.net/e3thP.png)

My implementation is quite noisy for some reason. It doesn't look any better than box filter.

**Tags:** image-processing, interpolation, filtering

> **Comments:**
> **username_1**: I notice in username_2's version there is a linear->sRGB conversion which may account for his version looking brighter & flatter. I would suggest looking at stratifying the samples too  to lower graininess, unless your sampler.Random() already does that ?  
> **username_0**: Yes I have tried stratified samples too which improves it but nothing like username_2's. Also I wonder the second formula in Wiki page should be modified to "a**2" in the numerator instead of "a"?  
> **username_1**: Ideally a minimum of 256 samples should be used so that it can converge on the full 0-255 RGB output range accurately. As most of your samples are going to be either black or white in your input image, 16 average samples are going to give you a none-smooth / quantised range. Although your filter is going to negate most of that, I didn't study too closely how your filter is working. I still think the linear->sRGB is responsible for most of Nathans 'flatness'  
> **username_0**: Thank Paul. Could explain a bit more how  linear->sRGB works? Do you refer to tone-mapping?  
> **username_2**: BTW, the factor $a/\pi^2$ in the weights is not needed because the final image is normalized by the sum of weights anyway, so any constant factor will cancel out.  
> **username_2**: Re: sRGB, no it's not tone-mapping, it's a conventional nonlinear encoding of pixel values in images to obtain better precision in the dark values where our eyes are more sensitive. You can read here for more: https://en.wikipedia.org/wiki/SRGB#The_sRGB_transfer_function_(%22gamma%22)  
> **username_0**: Paul, your guess was correct. After applying the conversion to sRGB it soomths out the image. Thanks.  

---

## 0 answers

