Perlin Noise

void Main()
{
    // Uses Oxyplot.WindowsForms.    
    var noise = new PerlinNoise();

    var results = new List<double>();
    for(double i = 0.0d; i < 100; i += 0.01d)
    {
        results.Add(noise.Noise(i));
    }
    
    var scatterSeries = new ScatterSeries();
    scatterSeries.Points.AddRange(results.Select((r, i) => new ScatterPoint(i, r, 1d)));
    
    var chart = new PlotModel();
    chart.Series.Add(scatterSeries);
    
    var view = new PlotView()
    {
	    Model = chart,
    };
	view.Dump();
    
    CreateRandomTexture(1366, 300).Dump();
}

public static class OxyPlotExtensions
{
   public static void AddScatterSeries(this PlotModel model, IEnumerable<double> xSeries, IEnumerable<double> ySeries)
   {
       model.AddScatterSeries(xSeries, ySeries, OxyColors.Automatic);
   }
 
   public static void AddScatterSeries(this PlotModel model, IEnumerable<double> xSeries, IEnumerable<double> ySeries, OxyColor color)
   {
       var scatterSeries = new ScatterSeries()
       {
           MarkerFill = color,
           MarkerSize = 1,
       };
 
       foreach (var item in xSeries.Zip(ySeries, (x, y) => new { x, y }))
       {
           scatterSeries.Points.Add(new ScatterPoint(item.x, item.y));
       }
 
       model.Series.Add(scatterSeries);
   }
}

public Bitmap CreateRandomTexture(int width, int height)
{
    var noise = new PerlinNoise();

    var bmp = new Bitmap(width, height);
    
    var xOffset = 30d;
    
    for(int x = 0; x < width; x ++)
    {
        var yOffset = 500d;
    
        for(int y = 0; y < height; y ++)
        {
            var noiseValue = noise.Noise(xOffset, yOffset);
        
            var rgb = (int)Math.Abs(Map(noiseValue, -1d, 1d, 100, 255));
        
            bmp.SetPixel(x, y, Color.FromArgb(rgb, rgb, rgb));
            
            yOffset += 0.01d;
        }
        
        xOffset += 0.01d;
    }
    
    return bmp;
}

public double Map(double value, double sourceRangeMinimum, double sourceRangeMaximum, double targetRangeMinimum, double targetRangeMaximum)
{
    if ((sourceRangeMaximum - sourceRangeMinimum) == 0)
    {
        return (targetRangeMinimum + targetRangeMaximum) / 2;
    }
    return targetRangeMinimum + (value - sourceRangeMinimum) * (targetRangeMaximum - targetRangeMinimum) / (sourceRangeMaximum - sourceRangeMinimum);
}

public class PerlinNoise
{
    /* coherent noise function over 1, 2 or 3 dimensions */
    /* (copyright Ken Perlin) */
    const int B = 0x100; // 256
    const int BM = 0xff;
    
    const int N = 0x1000;
    const int NP = 12;   /* 2^N */
    const int NM = 0xfff;
    
    int[] p = new int[B + B + 2];
    double[][] g3 = new double[B + B + 2][]; // each child needs 3 
    double[][] g2 = new double[B + B + 2][]; // each child needs 2
    double[] g1 = new double[B + B + 2];
    
    private Random rnd = new Random();
    
    public PerlinNoise()
    {
        int i, j, k;
    
        SizeArray(g2, 2);
        SizeArray(g3, 3);
    
        for (i = 0 ; i < B ; i++) {
            p[i] = i;
    
            g1[i] = (double)((rnd.Next() % (B + B)) - B) / B;

            // Size the array.
            for (j = 0 ; j < 2 ; j++)
                g2[i][j] = (float)((rnd.Next() % (B + B)) - B) / B;
            Normalize(g2[i]);
    
            g3[i] = new double[3];
            for (j = 0 ; j < 3 ; j++)
                g3[i][j] = (float)((rnd.Next() % (B + B)) - B) / B;
            Normalize(g3[i]);
        }
    
        while (--i > 0) {
            k = p[i];
            p[i] = p[j = rnd.Next() % B];
            p[j] = k;
        }
    
        for (i = 0 ; i < B + 2 ; i++) {
            p[B + i] = p[i];
            g1[B + i] = g1[i];
            for (j = 0 ; j < 2 ; j++)
                g2[B + i][j] = g2[i][j];
            for (j = 0 ; j < 3 ; j++)
                g3[B + i][j] = g3[i][j];
        }
    }
    
    private static double s_curve(double t) 
    {
        return t * t * (3.0d - 2.0d * t);
    }
    
    private static double lerp(double t, double a, double b)
    {
        return a + t * (b - a);
    }

    public double Noise(double x)
    {     
    	double t = x + N;
        int bx0 = ((int)t) & BM;
        int bx1 = (bx0+1) & BM;
        double rx0 = t - (int)t;
        double rx1 = rx0 - 1.0d;
        
        double sx = s_curve(rx0);
    
        double u = rx0 * g1[ p[ bx0 ] ];
        double v = rx1 * g1[ p[ bx1 ] ];
    
        return lerp(sx, u, v);
    }
    
    private static double at2(double[] q, double rx, double ry)
    {
        return rx * q[0] + ry * q[1];
    }
    
    public double Noise(double x, double y)
    {
        // setup(0, bx0,bx1, rx0,rx1);
    	double t = x + N;
        int bx0 = ((int)t) & BM;
        int bx1 = (bx0+1) & BM;
        double rx0 = t - (int)t;
        double rx1 = rx0 - 1.0d;
        
        //setup(1, by0,by1, ry0,ry1);
    	t = y + N;
        int by0 = ((int)t) & BM;
        int by1 = (by0+1) & BM;
        double ry0 = t - (int)t;
        double ry1 = ry0 - 1.0d;
    
        int i = p[ bx0 ];
        int j = p[ bx1 ];
    
        int b00 = p[ i + by0 ];
        int b10 = p[ j + by0 ];
        int b01 = p[ i + by1 ];
        int b11 = p[ j + by1 ];
    
        double sx = s_curve(rx0);
        double sy = s_curve(ry0);
    
        double[] q;
        q = g2[b00]; 
        double u = at2(q, rx0,ry0);
        q = g2[ b10 ] ; 
        double v = at2(q, rx1,ry0);
        double a = lerp(sx, u, v);
    
        q = g2[ b01 ] ; 
        u = at2(q, rx0,ry1);
        q = g2[ b11 ] ; 
        v = at2(q, rx1,ry1);
        double b = lerp(sx, u, v);
    
        return lerp(sy, a, b);
    }
    
    private static double at3(double[] q, double rx, double ry, double rz) 
    {
        return rx * q[0] + ry * q[1] + rz * q[2];
    }
    
    public double Noise(double x, double y, double z)
    {
        int bx0, bx1, by0, by1, bz0, bz1, b00, b10, b01, b11;
        double rx0, rx1, ry0, ry1, rz0, rz1, sy, sz, a, b, c, d, t, u, v;
        int i, j;
        
        // setup(0, bx0,bx1, rx0,rx1);
        t = x + N;
        bx0 = ((int)t) & BM;
        bx1 = (bx0+1) & BM;
        rx0 = t - (int)t;
        rx1 = rx0 - 1.0d;
        
        //setup(1, by0,by1, ry0,ry1);
    	t = y + N;
        by0 = ((int)t) & BM;
        by1 = (by0+1) & BM;
        ry0 = t - (int)t;
        ry1 = ry0 - 1.0d;
        
        // setup(2, bz0,bz1, rz0,rz1);
        t = z + N;
        bz0 = ((int)t) & BM;
        bz1 = (bz0+1) & BM;
        rz0 = t - (int)t;
        rz1 = rz0 - 1.0d;
    
        i = p[ bx0 ];
        j = p[ bx1 ];
    
        b00 = p[ i + by0 ];
        b10 = p[ j + by0 ];
        b01 = p[ i + by1 ];
        b11 = p[ j + by1 ];
    
        t  = s_curve(rx0);
        sy = s_curve(ry0);
        sz = s_curve(rz0);
        
        double[] q;
        q = g3[ b00 + bz0 ] ; u = at3(q, rx0,ry0,rz0);
        q = g3[ b10 + bz0 ] ; v = at3(q, rx1,ry0,rz0);
        a = lerp(t, u, v);
    
        q = g3[ b01 + bz0 ] ; u = at3(q, rx0,ry1,rz0);
        q = g3[ b11 + bz0 ] ; v = at3(q, rx1,ry1,rz0);
        b = lerp(t, u, v);
    
        c = lerp(sy, a, b);
    
        q = g3[ b00 + bz1 ] ; u = at3(q, rx0,ry0,rz1);
        q = g3[ b10 + bz1 ] ; v = at3(q, rx1,ry0,rz1);
        a = lerp(t, u, v);
    
        q = g3[ b01 + bz1 ] ; u = at3(q, rx0,ry1,rz1);
        q = g3[ b11 + bz1 ] ; v = at3(q, rx1,ry1,rz1);
        b = lerp(t, u, v);
    
        d = lerp(sy, a, b);
    
        return lerp(sz, c, d);
    }
    
    private static void Normalize(double[] vec)
    {
        var sumOfSquares = vec.Select(v => v * v).Sum();
        var squareRootOfSumOfSquares = Math.Sqrt(sumOfSquares);
        
        for(int i = 0; i < vec.Length; i ++)
        {
            vec[i] = vec[i] / squareRootOfSumOfSquares;
        }
    }
    
    private static void SizeArray(double[][] array, int childDimensions)
    {
        for(int i = 0; i < array.Length; i++)
        {
            array[i] = new double[childDimensions];
        }
    }
}