Example of a least square caculation, Matrix calculation based on CMSIS-DSP library, x(t)=a+b*t+c*t^2.

Dependencies:   dsp mbed

Comment out line 226 of "arm_math.h" to prevent a lot of compiler warnings

Revision:
0:9087ea0ad689
--- /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 );
+    }
+}