Quansight-Labs / numpy.net

A port of NumPy to .Net
BSD 3-Clause "New" or "Revised" License
131 stars 14 forks source link

np.interp not implemented #15

Closed Karthick47v2 closed 2 years ago

Karthick47v2 commented 2 years ago

Hi, i saw that, its mentioned that np.interp not working in one of commits.. Did you fix it?. I also tried to implement it myself by viewing numpy's implementation. But they mentioned like its a wrapper for interpand complex_interp which are PyCFunction and i couldn't find more details about that. Could you share any resources so that i can catch up?

@KevinBaselinesw @teoliphant

KevinBaselinesw commented 2 years ago

You are correct. interp is not implemented in numpydotnet. It was a huge task to do that.

Instead, I recommend that you use this package instead which produces the same results:

https://numerics.mathdotnet.com/

That is what I did in the following unit tests.

` [TestMethod] public void test_interp1d() { var x = np.arange(2, 12, 2); var y = np.arange(1, 6, 1);

        double[] a = new double[] { 2, 4, 6, 8, 10 };
        double[] b = new double[] { 0.51341712, 0.26359714, 0.13533528, 0.06948345, 0.03567399 };

        var ynew = MathNet.Numerics.Interpolate.Linear(a, b);

        print(x);
        print(y);
        print(ynew);

       // AssertArray(ynew, new double[] { 0.51341712, 0.26359714, 0.13533528, 0.06948345, 0.03567399 });

    }

    [TestMethod]
    public void test_interp1d_2()
    {
        ndarray xarray = np.arange(0, 10, dtype: np.Float64);
        ndarray yarray = np.array(new double[0]);
        ndarray xnewarray = np.arange(0, 9, 0.1, dtype: np.Float64);

        double[] x = xarray.rawdata().datap as double[];
        double[] y = new double[] { 1.0, 0.71653131, 0.51341712, 0.36787944, 0.26359714, 0.1888756,
                                         0.13533528, 0.09697197, 0.06948345, 0.04978707 };

        var f = MathNet.Numerics.Interpolate.Linear(x, y);

        double[] xnew = xnewarray.rawdata().datap as double[];
        double[] ynew = new double[xnew.Length];

        int index = 0;
        foreach (var xn in xnew)
        {
            ynew[index] = Math.Round(f.Interpolate(xn), 8, MidpointRounding.AwayFromZero);
            index++;
        }

        ndarray ynewarray = np.array(ynew);

        print(xarray);
        print(yarray);
        print(xnewarray);
        print(ynewarray);
    }`
Karthick47v2 commented 2 years ago

@KevinBaselinesw Thank you

rainyl commented 1 year ago

Although this issue has been closed and Kevin provided mathnet to solve it, I want to talk more about it. TLDR: the interpolation method of mathnet and numpy seems to be different, so I think it is, some how, necessary to implement np.interp() in NumpyDotNet.

For numpy, points besides xp and fp will use the minimum or maximum value directly (or left and right parameter if provided). But MathNet will extrapolate the y values using spline, in some cases, it will be useful, but not intuitive and friendly for a numpy user.

Here is the test example:

public void test_interp_1()
        {
            var xp = new float[] { 1, 2, 3 };
            var fp = new float[] { 3, 2, 0 };

            var a = np.interp(2.5, xp, fp);
            print("a: ", a);  // 1.0

            // 3.0, 3.0, 2.5, 0.56, 0.0
            var xb = new double[] { 0, 1, 1.5, 2.72, 3.14 };
            var b_real = np.array(new double[] { 3.0, 3.0, 2.5, 0.56, 0.0 });
            var b = np.interp(xb, xp, fp);
            // mathnet
            var intp = MathNet.Numerics.Interpolate.Linear(
                Array.ConvertAll<float, double>(xp, (x)=>(double)x),
                Array.ConvertAll<float, double>(fp, (x)=>(double)x));
            var b2 = new double[xb.Length];
            for (int i = 0; i < xb.Length; i++)
                b2[i] = intp.Interpolate(xb[i]);
            //print
            print("b (np.interp): ", b);
            print("b (mathnet): ", b2);
            b = np.round(b, 2);
            Assert.IsTrue(b.Equals(b_real).All().AsBoolArray().First());

            float UNDEF = -99.0f;
            var c = np.interp(new float[] { 3.14f, -1f }, xp, fp, left: UNDEF, right: UNDEF);
            print("c (left, right): ", c);  // -99.0, -99.0

            // with period
            var xd = new float[] { -180, -170, -185, 185, -10, -5, 0, 365 };
            var xpp = new float[] { 190, -190, 350, -350 };
            var fpp = new float[] { 5, 10, 3, 4 };
            // [7.5, 5., 8.75, 6.25, 3., 3.25, 3.5, 3.75]
            var d_real = np.array(new double[] { 7.5, 5.0, 8.75, 6.25, 3.0, 3.25, 3.5, 3.75 });
            var d = np.interp(xd, xpp, fpp, period: 360);
            Assert.IsTrue(d.Equals(d_real).All().AsBoolArray().First());
            print("d (period): ", d);
        }

and the result: image

So, I implemented the np.interp() according to the source code of Numpy (https://github.com/numpy/numpy/blob/main/numpy/core/src/multiarray/compiled_base.c from line 497, https://github.com/numpy/numpy/blob/main/numpy/lib/function_base.py Line 1462), but only for real numbers, I have almost returned all my knowledge about complex numbers to my teacher, so, someone more familiar with it may do some help.

I have opened a pull request, #50

Finally, the implementation of mine may be inefficient, if any one can improve it, just do it.

KevinBaselinesw commented 1 year ago

thank you. When I get a little bit more time, I will try to release this in an official build.

KevinBaselinesw commented 1 year ago

I just uploaded 0.9.83.8 with your implementation. I modified the API to use doubles instead of floats. Most data types will be automatically upscaled to double by .NET.

Thank you very much for your contribution in getting np.interp implemented!

rainyl commented 1 year ago

Great! Hoping this project will get more attention 😄