Example of a least square caculation, Matrix calculation based on CMSIS-DSP library, x(t)=a+b*t+c*t^2.
Comment out line 226 of "arm_math.h" to prevent a lot of compiler warnings
main.cpp
- Committer:
- GerritPathuis
- Date:
- 2015-03-10
- Revision:
- 0:9087ea0ad689
File content as of revision 0:9087ea0ad689:
/////////////////////////////////////// // FRDM-K22F // // Least square calculation (Matrix) // // Based on the DSP libtary // // Digital Signal Processing // // // // Results can be read via Tera Term // /////////////////////////////////////// #include "mbed.h" #include "arm_math.h" #define NUMSAMPLES 51 #define NUMUNKNOWS 3 Serial pc(USBTX, USBRX); // tx, rx // time is evenly spaced bit this in not required float32_t tData[NUMSAMPLES] = { 0.0f, 0.1f, 0.2f, 0.3f, 0.4f, 0.5f, 0.6f, 0.7f, 0.8f, 0.9f, 1.0f, 1.1f, 1.2f, 1.3f, 1.4f, 1.5f, 1.6f, 1.7f, 1.8f, 1.9f, 2.0f, 2.1f, 2.2f, 2.3f, 2.4f, 2.5f, 2.6f, 2.7f, 2.8f, 2.9f, 3.0f, 3.1f, 3.2f, 3.3f, 3.4f, 3.5f, 3.6f, 3.7f, 3.8f, 3.9f, 4.0f, 4.1f, 4.2f, 4.3f, 4.4f, 4.5f, 4.6f, 4.7f, 4.8f, 4.9f, 5.0f }; float32_t xData[NUMSAMPLES] = { 07.4213f, 21.7231f, -7.2828f, 21.2254f, 20.2221f, 10.3585f, 20.2033f, 29.2690f, 57.7152f, 53.6075f, 22.8209f, 59.8714f, 43.1712f, 38.4436f, 46.0499f, 39.8803f, 41.5188f, 55.2256f, 55.1803f, 55.6495f, 49.8920f, 34.8721f, 50.0859f, 57.0099f, 47.3032f, 50.8975f, 47.4671f, 38.0605f, 41.4790f, 31.2737f, 42.9272f, 24.6954f, 23.1770f, 22.9120f, 3.2977f, 35.6270f, 23.7935f, 12.0286f, 25.7104f, -2.4601f, 06.7021f, 01.6804f, 02.0671f, -2.2891f, -16.2070f, -14.2204f, -20.1780f, -18.9303f, -20.4859f, -25.8338f, -47.2892f }; float32_t AData[NUMSAMPLES * NUMUNKNOWS]; float32_t ATData[NUMSAMPLES * NUMUNKNOWS]; float32_t ATAData[NUMUNKNOWS * NUMUNKNOWS]; float32_t invATAData[NUMUNKNOWS * NUMUNKNOWS]; float32_t BData[NUMSAMPLES * NUMSAMPLES]; float32_t cData[NUMUNKNOWS]; arm_matrix_instance_f32 t= {NUMSAMPLES, 1, tData}; arm_matrix_instance_f32 x= {NUMSAMPLES, 1, xData}; arm_matrix_instance_f32 A= {NUMSAMPLES, NUMUNKNOWS, AData}; arm_matrix_instance_f32 AT= {NUMUNKNOWS, NUMSAMPLES, ATData}; arm_matrix_instance_f32 ATA= {NUMUNKNOWS, NUMUNKNOWS, ATAData}; arm_matrix_instance_f32 invATA= {NUMUNKNOWS, NUMUNKNOWS, invATAData}; arm_matrix_instance_f32 B= {NUMUNKNOWS, NUMSAMPLES, BData}; arm_matrix_instance_f32 c= {NUMUNKNOWS, 1, cData}; main(void) { int i, time; float y; pc.printf("/n/nMatrix calculation (Least squares) with 3 Unknowns \n\r"); pc.printf("We calculates a, b and c from; x(t)= a + b*t + c*t*t \n\r"); pc.printf("Results approx a= 8.701225, b= 38.880058, c= -9.792965\n\r"); y= sqrtf(xData[0]); cData[0]= y; for(i=0; i<NUMSAMPLES; i++) { AData[i*NUMUNKNOWS +0] = 1.0f; AData[i*NUMUNKNOWS +1] = tData[i]; AData[i*NUMUNKNOWS +2] = tData[i] * tData[i]; } /* calculation of A transpose */ arm_mat_trans_f32(&A, &AT); arm_mat_mult_f32(&AT, &A, &ATA); arm_mat_inverse_f32(&ATA, &invATA); arm_mat_mult_f32(&invATA, &AT, &B); arm_mat_mult_f32(&B, &x, &c); pc.printf("Results a= %f, b= %f, c= %f \n\r", cData[0], cData[1], cData[2]); for(time=0; time<6; time++) { y= cData[0] + cData[1]* time + cData[2] * time * time; pc.printf("time= %d, x= %f \n\r", time, y ); } }