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
Revision 0:9087ea0ad689, committed 2015-03-10
- Comitter:
- GerritPathuis
- Date:
- Tue Mar 10 12:14:31 2015 +0000
- Commit message:
- Program works
Changed in this revision
diff -r 000000000000 -r 9087ea0ad689 dsp.lib --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/dsp.lib Tue Mar 10 12:14:31 2015 +0000 @@ -0,0 +1,1 @@ +http://mbed.org/users/simon/code/dsp/#208cacc9d789
diff -r 000000000000 -r 9087ea0ad689 main.cpp --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/main.cpp Tue Mar 10 12:14:31 2015 +0000 @@ -0,0 +1,88 @@ +/////////////////////////////////////// +// 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 ); + } +}
diff -r 000000000000 -r 9087ea0ad689 mbed.bld --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/mbed.bld Tue Mar 10 12:14:31 2015 +0000 @@ -0,0 +1,1 @@ +http://mbed.org/users/mbed_official/code/mbed/builds/7e07b6fb45cf \ No newline at end of file