1 ////////////////////////////////////////////////////////////////////// 2 // 连分式等距插值 3 ////////////////////////////////////////////////////////////////////// 4 static float GetValuePqs(const void* valuesPtr, int stride, int n, float t) 5 {  6     int i,j,k,m, l; 7     float z,hh,xi,xj; 8     float b[8]; 9 10     // 初值11     z = 0.0;12 13     // 特例处理14     if (n < 1) 15     {16         return(z);17     }18     if (n == 1) 19     { 20         z = YfGetFloatValue(valuesPtr, stride, 0);  21         return(z);22     }23     if (n == 2)24     { 25         float y0 = YfGetFloatValue(valuesPtr, stride, 0); 26         float y1 = YfGetFloatValue(valuesPtr, stride, 1); 27         z = y0 + (y1 - y0)*t;28         return(z);29     }30 31     float xStep = 1.0f/(n - 1);32 33     // 连分式插值34     if (n <= 8) 35     { 36         k = 0; 37         m = n;38     }39     else if (t < (4.0f*xStep)) 40     { 41         k=0; 42         m=8;43     }44     else if (t > ((n-5)*xStep)) 45     { 46         k = n-8; 47         m = 8;48     }49     else 50     { 51         k = (int)(t/xStep)-3; 52         m = 8;53     }54 55     b[0] = YfGetFloatValue(valuesPtr, stride, k); 56     for (i = 2; i <= m; i++)57     { 58         hh = YfGetFloatValue(valuesPtr, stride, i+k-1); 59         l = 0; 60         j = 1;61 62         while ((l == 0) && (j <= i-1))63         { 64             if (fabs(hh-b[j-1])+1.0f == 1.0f) 65             {66                 l = 1;67             }68             else69             { 70                 xi = (i+k-1)*xStep;71                 xj = (j+k-1)*xStep;72                 hh = (xi-xj)/(hh-b[j-1]);73             }74 75             j = j+1;76         }77 78         b[i-1]=hh;79         if (l != 0) 80         {81             b[i-1] = 1.0e+35F;82         }83     }84 85     z = b[m-1];86     for (i = m-1; i >= 1; i--)87     {88         z = b[i-1]+(t-(i+k-1)*xStep)/z;89     }90 91     return(z);92 }