C# Interpolazione lineare

Question

Sto riscontrando problemi nell'interpolazione di un file di dati, che ho convertito da .csv in un array X e un array Y in cui X [0] corrisponde al punto Y [0] ad esempio.

Ho bisogno di interpolare tra i valori per darmi un file liscio alla fine. Sto usando un Picoscope per emettere la funzione, il che significa che ogni riga è equamente distanziata nel tempo, quindi usando solo i valori Y, ed è per questo che sto cercando di farlo in un modo strano quando vedi il mio codice.

Il tipo di valori che deve affrontare sono:

X  Y 
0  0 
2.5 0 
2.5 12000 
7.5 12000 
7.5 3000 
10 3000 
10 6000 
11 6625.254 
12 7095.154 

Così dove 2 valori Y accanto all'altro sono gli stessi sua una linea retta tra loro ma quando variano come da X = 10 sulla guardie diventa un'onda sinusoidale in questo esempio.

Sono stato a guardare le equazioni per l'interpolazione ecc. E altri post qui ma non ho fatto quel tipo di matematica per anni e purtroppo non riesco a capirlo più, perché ci sono 2 incognite e posso Pensa a come programmarlo in C#.

Quello che ho è:

int a = 0, d = 0, q = 0; 
bool up = false; 
double times = 0, change = 0, points = 0, pointchange = 0; 
double[] newy = new double[8192]; 
while (a < sizeoffile-1 && d < 8192) 
{ 
    Console.Write("..."); 
    if (Y[a] == Y[a+1])//if the 2 values are the same add correct number of lines in to make it last the correct time 
    { 
     times = (X[a] - X[a + 1]);//number of repetitions 
     if ((X[a + 1] - X[a]) > times)//if that was a negative number try the other way around 
      times = (X[a + 1] - X[a]);//number of repetitions 
     do 
     { 
      newy[d] = Y[a];//add the values to the newy array which replaces y later on in the program 
      d++;//increment newy position 
      q++;//reduce number of reps in this loop 
     } 
     while (q < times + 1 && d < 8192); 
     q = 0;//reset reps 
    } 
    else if (Y[a] != Y[a + 1])//if the 2 values are not the same interpolate between them 
    { 
     change = (Y[a] - Y[a + 1]);//work out difference between the values 
     up = true;//the waveform is increasing 
     if ((Y[a + 1] - Y[a]) > change)//if that number was a negative try the other way around 
     { 
      change = (Y[a + 1] - Y[a]);//work out difference bwteen values 
      up = false;//the waveform is decreasing 
     } 
     points = (X[a] - X[a + 1]);//work out amount of time between given points 
     if ((X[a + 1] - X[a]) > points)//if that was a negative try other way around 
      points = (X[a + 1] - X[a]);///work out amount of time between given points 
     pointchange = (change/points);//calculate the amount per point in time the value changes by 
     if (points > 1)//any lower and the values cause errors 
     { 
      newy[d] = Y[a];//load first point 
      d++; 
      do 
      { 
       if (up == true)// 
        newy[d] = ((newy[d - 1]) + pointchange); 
       else if (up == false) 
        newy[d] = ((newy[d - 1]) - pointchange); 
       d++; 
       q++; 
      } 
      while (q < points + 1 && d < 8192); 
      q = 0; 
     } 
     else if (points != 0 && points > 0) 
     { 
      newy[d] = Y[a];//load first point 
      d++; 
     } 
    } 
    a++; 
} 

e questo crea una stretta forma d'onda, ma è ancora molto Steppy.

Quindi qualcuno può capire perché questo non è molto preciso? Come migliorare la precisione? Oppure un modo diverso di farlo usando gli array?

Ringraziamenti per lo sguardo :)

Answer 1

provare questo metodo per me:

static public double linear(double x, double x0, double x1, double y0, double y1) 
{ 
    if ((x1 - x0) == 0) 
    { 
     return (y0 + y1)/2; 
    } 
    return y0 + (x - x0) * (y1 - y0)/(x1 - x0); 
} 

Effettivamente si dovrebbe essere in grado di prendere gli array e usarlo in questo modo:

var newY = linear(X[0], X[0], X[1], Y[0], Y[1]); 

ho tirato il codice da here, ma verificato che l'algoritmo corrispondeva alla teoria here, e quindi I t hink è giusto. Tuttavia, probabilmente dovresti considerare di usare l'interpolazione polinomiale se è ancora steppy, per favore nota il link della teoria, mostra che l'interpolazione lineare produce onde steppy.

Così, il primo link che ho dato, dove ho preso questo codice da, ha anche un algoritmo polinomiale:

static public double lagrange(double x, double[] xd, double[] yd) 
{ 
    if (xd.Length != yd.Length) 
    { 
     throw new ArgumentException("Arrays must be of equal length."); //$NON-NLS-1$ 
    } 
    double sum = 0; 
    for (int i = 0, n = xd.Length; i < n; i++) 
    { 
     if (x - xd[i] == 0) 
     { 
      return yd[i]; 
     } 
     double product = yd[i]; 
     for (int j = 0; j < n; j++) 
     { 
      if ((i == j) || (xd[i] - xd[j] == 0)) 
      { 
       continue; 
      } 
      product *= (x - xd[i])/(xd[i] - xd[j]); 
     } 
     sum += product; 
    } 
    return sum; 
} 

Per utilizzare questo quello che si sta andando ad avere per decidere come si vuole intensificare i vostri valori x, quindi diciamo che volevamo farlo trovando il punto medio tra l'iterazione corrente e la successiva:

for (int i = 0; i < X.Length; i++) 
{ 
    var newY = lagrange(new double[] { X[i]d, X[i+1]d }.Average(), X, Y); 
} 

prega di notare che ci sarà più di questo ciclo, come garantire che v'è un i+1 e tale , ma volevo vedere se c Potrei darti un inizio.

Answer 2

Theoretical base at Wolfram

La soluzione qui di seguito calcola le medie dei valori Y di punti dati con stesso X, così come la funzione polyfit Matlab fa

Linq e NET versione quadro> 3.5 sono obbligatorie per questa rapida API. Commenti all'interno del codice.

using System; 
using System.Collections.Generic; 
using System.Linq; 


/// <summary> 
/// Linear Interpolation using the least squares method 
/// <remarks>http://mathworld.wolfram.com/LeastSquaresFitting.html</remarks> 
/// </summary> 
public class LinearLeastSquaresInterpolation 
{ 
    /// <summary> 
    /// point list constructor 
    /// </summary> 
    /// <param name="points">points list</param> 
    public LinearLeastSquaresInterpolation(IEnumerable<Point> points) 
    { 
     Points = points; 
    } 
    /// <summary> 
    /// abscissae/ordinates constructor 
    /// </summary> 
    /// <param name="x">abscissae</param> 
    /// <param name="y">ordinates</param> 
    public LinearLeastSquaresInterpolation(IEnumerable<float> x, IEnumerable<float> y) 
    { 
     if (x.Empty() || y.Empty()) 
      throw new ArgumentNullException("null-x"); 
     if (y.Empty()) 
      throw new ArgumentNullException("null-y"); 
     if (x.Count() != y.Count()) 
      throw new ArgumentException("diff-count"); 

     Points = x.Zip(y, (unx, uny) => new Point { x = unx, y = uny }); 
    } 

    private IEnumerable<Point> Points; 
    /// <summary> 
    /// original points count 
    /// </summary> 
    public int Count { get { return Points.Count(); } } 

    /// <summary> 
    /// group points with equal x value, average group y value 
    /// </summary> 
                public IEnumerable<Point> UniquePoints 
{ 
    get 
    { 
     var grp = Points.GroupBy((p) => { return p.x; }); 
     foreach (IGrouping<float,Point> g in grp) 
     { 
      float currentX = g.Key; 
      float averageYforX = g.Select(p => p.y).Average(); 
      yield return new Point() { x = currentX, y = averageYforX }; 
     } 
    } 
} 
    /// <summary> 
    /// count of point set used for interpolation 
    /// </summary> 
    public int CountUnique { get { return UniquePoints.Count(); } } 
    /// <summary> 
    /// abscissae 
    /// </summary> 
    public IEnumerable<float> X { get { return UniquePoints.Select(p => p.x); } } 
    /// <summary> 
    /// ordinates 
    /// </summary> 
    public IEnumerable<float> Y { get { return UniquePoints.Select(p => p.y); } } 
    /// <summary> 
    /// x mean 
    /// </summary> 
    public float AverageX { get { return X.Average(); } } 
    /// <summary> 
    /// y mean 
    /// </summary> 
    public float AverageY { get { return Y.Average(); } } 

    /// <summary> 
    /// the computed slope, aka regression coefficient 
    /// </summary> 
    public float Slope { get { return ssxy/ssxx; } } 

    // dotvector(x,y)-n*avgx*avgy 
    float ssxy { get { return X.DotProduct(Y) - CountUnique * AverageX * AverageY; } } 
    //sum squares x - n * square avgx 
    float ssxx { get { return X.DotProduct(X) - CountUnique * AverageX * AverageX; } } 

    /// <summary> 
    /// computed intercept 
    /// </summary> 
    public float Intercept { get { return AverageY - Slope * AverageX; } } 

    public override string ToString() 
    { 
     return string.Format("slope:{0:F02} intercept:{1:F02}", Slope, Intercept); 
    } 
} 

/// <summary> 
/// any given point 
/// </summary> 
public class Point 
{ 
    public float x { get; set; } 
    public float y { get; set; } 
} 

/// <summary> 
/// Linq extensions 
/// </summary> 
public static class Extensions 
{ 
    /// <summary> 
    /// dot vector product 
    /// </summary> 
    /// <param name="a">input</param> 
    /// <param name="b">input</param> 
    /// <returns>dot product of 2 inputs</returns> 
    public static float DotProduct(this IEnumerable<float> a, IEnumerable<float> b) 
    { 
     return a.Zip(b, (d1, d2) => d1 * d2).Sum(); 
    } 
    /// <summary> 
    /// is empty enumerable 
    /// </summary> 
    /// <typeparam name="T"></typeparam> 
    /// <param name="a"></param> 
    /// <returns></returns> 
    public static bool Empty<T>(this IEnumerable<T> a) 
    { 
     return a == null || a.Count() == 0; 
    } 
} 

utilizzarlo come:

var llsi = new LinearLeastSquaresInterpolation(new Point[] 
    { 
     new Point {x=1, y=1}, new Point {x=1, y=1.1f}, new Point {x=1, y=0.9f}, 
     new Point {x=2, y=2}, new Point {x=2, y=2.1f}, new Point {x=2, y=1.9f}, 
     new Point {x=3, y=3}, new Point {x=3, y=3.1f}, new Point {x=3, y=2.9f}, 
     new Point {x=10, y=10}, new Point {x=10, y=10.1f}, new Point {x=10, y=9.9f}, 
     new Point {x=100, y=100}, new Point{x=100, y=100.1f}, new Point {x=100, y=99.9f} 
    }); 

Oppure:

var llsi = new LinearLeastSquaresInterpolation(
    new float[] { 1, 2, 3, 4, 5, 6, 7, 8, 9, 1, 2, 3, 4, 5, 6, 7, 8, 9 }, 
    new float[] { 1, 2, 3, 4, 5, 6, 7, 8, 9, 1, 2, 3, 4, 5, 6, 7, 8, 9 });