curve fitting problem with arduino

Hi,

I try to read the data from infrared sensors and curve fitting the data by using matrix. I don’t know why the serial monitor cannot show me the value.
I want put time to array_x and put the infrared sensor value to array_y.
x-axis: time(second)
y-axis: the value of infrared sensor(cm)

Here is my code:

const int ii=10;
const int jj=4;
double a0, a1, a2, a3;
double array_a[4][1]={a0,a1,a2,a3};
double array_y[ii][1];              //the infrared sensor's value
double array_x[ii][jj];
double array_xt[jj][ii];              //transpose x
double array_yc[ii][1];              //change the y value
double A[jj][jj];                         //x* xt
double Y[jj][jj];                         //Inverse(A)
double array_Y_xt[jj][ii];          //Inverse(A) * xt
double array_Y_xt_xt[jj][ii];
long indxc[4], indxr[4], ipiv[4];    
long i, icol, irow, j, ir, ic;
float big, dum, pivinv, temp, bb;



void setup() {   
        Serial.begin(9600);   
        Invert();
        delay(100);   
}   
   
void loop() 
{
       for (int p=0; p<4; p++ ) 
       {
          for (int k=0; k<10;k++ )                 
           {
            array_x[k][p]=(k+1)^p;                     //3-order
            array_y[k][1]=read_value(0);             //put the infrared sensor's value to the 'y'matrix
           }
       }
       
       for (int p=0; p<4; p++ ) 
        {
          for (int k=0; k<10;k++ )                 
          {
           array_xt[p][k]=array_x[k][p];          //transpose of array_x
          }
        }           
       
       //array_xt * array_x
        for (int i=0;i<4;i++)
        {
           for (int j=0;j<4;j++)
           {
             for(int k=0, s=0; k<10; k++)
             A[i][j]=s+(array_xt[i][k])*(array_x[k][j]);
            
           }
        }
          
        //Invert(array_xt * array_x)=Y
        Invert();
        
        // Y * array_xt
        for (int i=0;i<4;i++)
        {
           for (int j=0;j<10;j++)
           {
             for(int k=0, s=0; k<4; k++)
             array_Y_xt[i][j]=s+(Y[i][k])*(array_xt[k][j]);
           }
        }
        
        //array_Y_xt * array_y
         for (int i=0;i<4;i++)
        {
           for (int j=0;j<1;j++)
           {
             for(int k=0, s=0; k<10; k++)
             array_a[i][j]=s+(array_Y_xt[i][k])*(array_y[k][j]);
           Serial.print("a");
           Serial.print(i);
           Serial.print(" ");
           Serial.print(array_a[i][j]);
           } 
         }
                 
   
       // array_x * array_a 
          for (int i=0;i<10;i++)
        {
           for (int j=0;j<1;j++)
           {
             for(int k=0, s=0; k<4; k++)
             { 
               array_yc[i][j]=s+(array_x[i][k])*(array_a[k][j]);
               Serial.print("  Sensor_ ");
               Serial.print(array_y[k][j]);
               Serial.print("   Change ");
               Serial.println(array_yc[k][j]);
              }
            }
          }
        
     
       
       delay(1000);       
 
 }
 
double read_value(byte pin) { 
        double tmp; 
        double value;
        tmp = analogRead(pin); 
        value = (6787.0 /(tmp - 3.0)) - 4.0;                                                   //DP2D12 According to http://www.acroname.com/robotics/info/articles/irlinear/irlinear.html
        if (tmp < 4.0 ) return 0;                                                              // invalid Value
        else if (value < 10.0 || value > 80.0) return 0.0;                                     // out of the range of infrared
        return value;                                                              
} 


 //calculate Matrix Inverse of A  
 void Invert(){
    ipiv[0] = -1;
    ipiv[1] = -1;
    ipiv[2] = -1;
    ipiv[3] = -1;
    Y[0][0] = A[0][0];
    Y[1][0] = A[1][0];
    Y[2][0] = A[2][0];
    Y[3][0] = A[3][0];
    Y[0][1] = A[0][1];
    Y[1][1] = A[1][1];
    Y[2][1] = A[2][1];
    Y[3][1] = A[3][1];
    Y[0][2] = A[0][2];
    Y[1][2] = A[1][2];
    Y[2][2] = A[2][2];
    Y[3][2] = A[3][2];
    Y[0][3] = A[0][3];
    Y[1][3] = A[1][3];
    Y[2][3] = A[2][3];
    Y[3][3] = A[3][3];
    for (i = 0; i < 4; i++) {
        big = 0.0f;
        for (j = 0; j < 4; j++) {
            if (ipiv[j] != 0) {
                if (ipiv[0] == -1) {
                    if ((bb = (float) fabs(Y[j][0])) > big) {
                        big = bb;
                        irow = j;
                        icol = 0;
                    }
                } else if (ipiv[0] > 0) {
                    return;
                }
                if (ipiv[1] == -1) {
                    if ((bb = (float) fabs((float) Y[j][1])) > big) {
                        big = bb;
                        irow = j;
                        icol = 1;
                    }
                } else if (ipiv[1] > 0) {
                    return;
                }
                if (ipiv[2] == -1) {
                    if ((bb = (float) fabs((float) Y[j][2])) > big) {
                        big = bb;
                        irow = j;
                        icol = 2;
                    }
                } else if (ipiv[2] > 0) {
                    return;
                }
                if (ipiv[3] == -1) {
                    if ((bb = (float) fabs((float) Y[j][3])) > big) {
                        big = bb;
                        irow = j;
                        icol = 3;
                    }
                } else if (ipiv[3] > 0) {
                    return;
                }
            }
        }
        ++(ipiv[icol]);
        if (irow != icol) {
            temp = Y[irow][0];
            Y[irow][0] = Y[icol][0];
            Y[icol][0] = temp;
            temp = Y[irow][1];
            Y[irow][1] = Y[icol][1];
            Y[icol][1] = temp;
            temp = Y[irow][2];
            Y[irow][2] = Y[icol][2];
            Y[icol][2] = temp;
            temp = Y[irow][3];
            Y[irow][3] = Y[icol][3];
            Y[icol][3] = temp;
        }
        indxr[i] = irow;
        indxc[i] = icol;
        if (Y[icol][icol] == 0.0) {
            return;
        }
        pivinv = 1.0f / Y[icol][icol];
        Y[icol][icol] = 1.0f;
        Y[icol][0] *= pivinv;
        Y[icol][1] *= pivinv;
        Y[icol][2] *= pivinv;
        Y[icol][3] *= pivinv;
        if (icol != 0) {
            dum = Y[0][icol];
            Y[0][icol] = 0.0f;
            Y[0][0] -= Y[icol][0] * dum;
            Y[0][1] -= Y[icol][1] * dum;
            Y[0][2] -= Y[icol][2] * dum;
            Y[0][3] -= Y[icol][3] * dum;
        }
        if (icol != 1) {
            dum = Y[1][icol];
            Y[1][icol] = 0.0f;
            Y[1][0] -= Y[icol][0] * dum;
            Y[1][1] -= Y[icol][1] * dum;
            Y[1][2] -= Y[icol][2] * dum;
            Y[1][3] -= Y[icol][3] * dum;
        }
        if (icol != 2) {
            dum = Y[2][icol];
            Y[2][icol] = 0.0f;
            Y[2][0] -= Y[icol][0] * dum;
            Y[2][1] -= Y[icol][1] * dum;
            Y[2][2] -= Y[icol][2] * dum;
            Y[2][3] -= Y[icol][3] * dum;
        }
        if (icol != 3) {
            dum = Y[3][icol];
            Y[3][icol] = 0.0f;
            Y[3][0] -= Y[icol][0] * dum;
            Y[3][1] -= Y[icol][1] * dum;
            Y[3][2] -= Y[icol][2] * dum;
            Y[3][3] -= Y[icol][3] * dum;
        }
    }
    if (indxr[3] != indxc[3]) {
        ir = indxr[3];
        ic = indxc[3];
        temp = Y[0][ir];
        Y[0][ir] = Y[0][ic];
        Y[0][ic] = temp;
        temp = Y[1][ir];
        Y[1][ir] = Y[1][ic];
        Y[1][ic] = temp;
        temp = Y[2][ir];
        Y[2][ir] = Y[2][ic];
        Y[2][ic] = temp;
        temp = Y[3][ir];
        Y[3][ir] = Y[3][ic];
        Y[3][ic] = temp;
    }
    if (indxr[2] != indxc[2]) {
        ir = indxr[2];
        ic = indxc[2];
        temp = Y[0][ir];
        Y[0][ir] = Y[0][ic];
        Y[0][ic] = temp;
        temp = Y[1][ir];
        Y[1][ir] = Y[1][ic];
        Y[1][ic] = temp;
        temp = Y[2][ir];
        Y[2][ir] = Y[2][ic];
        Y[2][ic] = temp;
        temp = Y[3][ir];
        Y[3][ir] = Y[3][ic];
        Y[3][ic] = temp;
    }
    if (indxr[1] != indxc[1]) {
        ir = indxr[1];
        ic = indxc[1];
        temp = Y[0][ir];
        Y[0][ir] = Y[0][ic];
        Y[0][ic] = temp;
        temp = Y[1][ir];
        Y[1][ir] = Y[1][ic];
        Y[1][ic] = temp;
        temp = Y[2][ir];
        Y[2][ir] = Y[2][ic];
        Y[2][ic] = temp;
        temp = Y[3][ir];
        Y[3][ir] = Y[3][ic];
        Y[3][ic] = temp;
    }
    if (indxr[0] != indxc[0]) {
        ir = indxr[0];
        ic = indxc[0];
        temp = Y[0][ir];
        Y[0][ir] = Y[0][ic];
        Y[0][ic] = temp;
        temp = Y[1][ir];
        Y[1][ir] = Y[1][ic];
        Y[1][ic] = temp;
        temp = Y[2][ir];
        Y[2][ir] = Y[2][ic];
        Y[2][ic] = temp;
        temp = Y[3][ir];
        Y[3][ir] = Y[3][ic];
        Y[3][ic] = temp;
    }
    // end of calculate matrix inverse*/
}

thanks a lot

You've probably run out of RAM:

double array_x[ii][jj];
double array_xt[jj][ii];              //transpose x
double array_yc[ii][1];              //change the y value
double A[jj][jj];                         //x* xt
double Y[jj][jj];                         //Inverse(A)
double array_Y_xt[jj][ii];          //Inverse(A) * xt
double array_Y_xt_xt[jj][ii];

ii * jj * 4 = 160 bytes, and you've got 6 of these. Try a 328 or Mega.

double a0, a1, a2, a3;
double array_a[4][1]={a0,a1,a2,a3};

Not sure why you're doing this.

I want put time to array_x and put the infrared sensor value to array_y.

What you are storing in array_x in the first for loop in the loop function is not the time.

I don't know why the serial monitor cannot show me the value.

Do you see anything being written to the serial monitor?

What you are storing in array_x in the first for loop in the loop function is not the time.

I want to put time to array_x and put the data to array_y. after 10 data, I can get the new data to fit curve.
could you please teach me how to modify the code to what I want?

double a0, a1, a2, a3;
double array_a[4][1]={a0,a1,a2,a3};

It’s not necessary.
maybedouble array_a[4][1];is enough

I want to put time to array_x

The current time, as in it's 4:09 AM PST, or the elapsed time since some event? If it's elapsed time, relative to what event?

Sorry, I didn't describe it very clear. The "time" is the delay time I set. for example, x y 1 18.93 2 20.00 3 17.98

You've declared over 1100 bytes of storage on a processor with only 1024 bytes of RAM, assuming you're using a 168 Arduino.

When the x value is the array_y index + 1, why is it necessary to store it? Why is it necessary to store any of this data in a two dimensional array? Why is it necessary to store integers in a double array?

What are you doing with the "curve" after you've computed it?

Another thing:

array_x[k][p]=(k+1)^p;

^ is the exclusive OR operator - not what you want here, I think.

When the x value is the array_y index + 1, why is it necessary to store it? Why is it necessary to store any of this data in a two dimensional array? Why is it necessary to store integers in a double array? What are you doing with the "curve" after you've computed it?

(I use cubic equation to do the curve fitting) If I don't store it, I couldn't do the curve fitting. I only know the way to do curve fitting is using matrix. I want to use curve fitting do self-controlling. If I get the "curve", I can also get inflection point. I want to start the motor after the inflection point.

^ is the exclusive OR operator

Would you please teach me how to describe "oder" in arduino?

I'd probably use two multiply operations.