function to calculate RMS value

Hello Experts,
I am trying to write a quick function to calculate the RMS values, without using any blocking calls(mainly loops)
This is what I have so far, it is not optimized for speed for which I need your help:
Here is the function:

Code:

---

/*
 * How to calculate the RMS
 * 0: True RMS - Square all the inputs within the time duration and add them, then take mean and finally square root the value.
 * 1: Approx RMS - Find the peak value within the time duration and divide by sqrt(2)
 */
#define MOTOR_RMS_TECHNEQUE 0    //This can be 0 or 1


#define NONBLOCKING_EXECUTION(a) (uint32_t)(__divf32(a,0.450))        //Timing for all functions running in 450ms Task


typedef struct _PhaseRMS_Vars_t_
{
    /// Temporary Phase U offset
    float32_t phaseU_offset;
    /// Temporary Phase V offset
    float32_t phaseV_offset;
    /// Temporary Phase W offset
    float32_t phaseW_offset;
    /// Phase U RMS
    float32_t phaseU_RMS;
    /// Phase V RMS
    float32_t phaseV_RMS;
    /// Phase W RMS
    float32_t phaseW_RMS;
    /// Phase U Peak
    float32_t phaseU_Peak;
    /// Phase V Peak
    float32_t phaseV_Peak;
    /// Phase W Peak
    float32_t phaseW_Peak;
    /// RMS Counter
    uint32_t rmsCounter;
    /// Sample aquition time
    float32_t time;
    /// phaseU squared and summed
    float32_t phaseUSumSquared;
    /// phaseV squared and summed
    float32_t phaseVSumSquared;
    /// phaseW squared and summed
    float32_t phaseWSumSquared;
} PhaseRMS_Vars_t;
extern PhaseRMS_Vars_t MotorCurrent[2];


void rmsCalc(PhaseRMS_Vars_t *iMotor, MOTOR_Vars_t *pMotor)
{
    //Mid point of bipolar range(0 to 4096)
    const uint16_t offset = 2048;


    //Get the data from the ADC and subtract the offset.
    /*
    * TODO Use the ADC PPB values in future.
    * The PPB values already have the offset subtracted from them.
    * So the offset value subtraction can be eliminated.
    */
    if(pMotor->MotorNum == 0)    //Motor A
    {
        iMotor->phaseU_offset = fabsf(M1_IFB_U - offset);
        iMotor->phaseV_offset = fabsf(M1_IFB_V - offset);
        iMotor->phaseW_offset = fabsf(M1_IFB_W - offset);
    }
    else if(pMotor->MotorNum == 1)    //Motor B
    {
        iMotor->phaseU_offset = fabsf(M2_IFB_U - offset);
        iMotor->phaseV_offset = fabsf(M2_IFB_V - offset);
        iMotor->phaseW_offset = fabsf(M2_IFB_W - offset);
    }


    //If speed is zero, div/0 error will occur, so make sure the RMS will only be calculated once the motor is running.
    if(pMotor->speed.SpeedRpm > 50)
    {
        //N = 120*F/P
        //From frequency(F) calculate the time period.
        iMotor->time = __divf32(1.0,__divf32(pMotor->speed.SpeedRpm * 4,120.0)) * 1000;  //time in ms
    }
    else
    {
        iMotor->time = 0;
    }
    //Convert time is milliseconds to counter ticks
    const uint16_t maxCount = NONBLOCKING_EXECUTION(iMotor->time);


    //Increment the counter
    iMotor->rmsCounter++;
#if(MOTOR_RMS_TECHNEQUE == 0)
    //Find the peak Value
    if(iMotor->rmsCounter < maxCount)
    {
        if(iMotor->phaseU_offset > iMotor->phaseU_Peak)
        {
            iMotor->phaseU_Peak = iMotor->phaseU_offset;
        }
        if(iMotor->phaseV_offset > iMotor->phaseV_Peak)
        {
            iMotor->phaseV_Peak = iMotor->phaseV_offset;
        }
        if(iMotor->phaseW_offset > iMotor->phaseW_Peak)
        {
            iMotor->phaseW_Peak = iMotor->phaseW_offset;
        }
    }
    else
    {
        /* Calculate RMS
        * This is still the RMS in terms of ADC count.
        * Still need to convert count to voltage, and finally voltage to actual current(To be done later)
        * */
        iMotor->phaseU_RMS = iMotor->phaseU_Peak * 0.7072f;
        iMotor->phaseV_RMS = iMotor->phaseV_Peak * 0.7072f;
        iMotor->phaseW_RMS = iMotor->phaseW_Peak * 0.7072f;
        iMotor->rmsCounter = 0;
        iMotor->phaseU_Peak = 0;
        iMotor->phaseV_Peak = 0;
        iMotor->phaseW_Peak = 0;
    }
#elif(MOTOR_RMS_TECHNEQUE == 1)
    if(iMotor->rmsCounter < maxCount)
    {
        //pow(iMotor->phaseU_offset,2) uses a lot of stack
        iMotor->phaseUSumSquared += (iMotor->phaseU_offset * iMotor->phaseU_offset);
        iMotor->phaseVSumSquared += (iMotor->phaseV_offset * iMotor->phaseV_offset);
        iMotor->phaseWSumSquared += (iMotor->phaseW_offset * iMotor->phaseW_offset);
    }
    else
    {
        /* Calculate RMS
        * This is still the RMS in terms of ADC count.
        * Still need to convert count to voltage, and finally voltage to actual current(To be done later)
        * */


        //The mean of the squared values is calculated and then the square root is found.
        iMotor->phaseU_RMS = __sqrt(__divf32(iMotor->phaseUSumSquared, maxCount));
        iMotor->phaseV_RMS = __sqrt(__divf32(iMotor->phaseVSumSquared, maxCount));
        iMotor->phaseW_RMS = __sqrt(__divf32(iMotor->phaseWSumSquared, maxCount));
        iMotor->phaseUSumSquared = 0;
        iMotor->phaseVSumSquared = 0;
        iMotor->phaseWSumSquared = 0;
        iMotor->rmsCounter = 0;
    }
#endif
}

---

That's not the definition of Root Mean Square (RMS). This is:
Attachment 16641
Something like this:

Code:

---

float rms( float *p, unsigned int samples )
{
  double sum;
  float *q;

  if ( ! samples )
    return 0.0;

  q = p + samples;

  sum = 0.0;
  while ( p < q )
  {
    sum += (double)*p * *p;
    p++;
  }

  return sqrt( sum / samples );
}

---

If you don't want to calculate the RMS of all values ever seen (e.g. only the last 10 samples), and you want to do it efficently, then this is how I would do it:

Code:

---

#include <stdio.h>
#include <math.h>
#include <stdlib.h>


struct rms_state {
  int count;
  int max_values;
  int index;
  int total;
  int *data_squared;  // History of values
};


struct rms_state *rms_init(int max_values) {
  struct rms_state *s = malloc(sizeof(struct rms_state));
  if(s == NULL) {
      return NULL;
  }
  s->count = 0;
  s->max_values = max_values;
  s->total = 0.0;
  s->index = 0;
  s->data_squared = malloc(sizeof(double)*max_values);
  if(s->data_squared == NULL) {
      free(s);
      return NULL;
  }
  return s;
}


double rms_update(struct rms_state *s, int value) {
  value *= value;

  // See if we need to remove the oldest value from the running total
  if(s->count != s->max_values)
      s->count++;
  else {
      s->total -= s->data_squared[s->index];
  }

  // Add the squared new value into the running total
  s->data_squared[s->index] = value;
  s->total += value;


  s->index++;
  if(s->index == s->max_values)
      s->index = 0;


  return sqrt(s->total / s->count);
}

void rms_free(struct rms_state *s) {
  free(s->data_squared);
  free(s);
}

int main(void) {
  // Allocate state info
  struct rms_state *s = rms_init(10);
  if(s == NULL) {
      printf("Out of mem?\n");
      return 0;
  }

  // Add 100 numbers to RMS value, one at a time.
  for(int i = 0; i < 100; i++) {
      int new_value = rand()%100;
      double rms_value = rms_update(s, new_value);
      printf("%8d %8f\n", new_value, rms_value);
  }
  // Flush out the RMS buffer - should drop back to zero.
  for(int i = 0; i < 15; i++) {
      int new_value = 0;
      double rms_value = rms_update(s, new_value);
      printf("%8d %8f\n", new_value, rms_value);
  }

  // Free resources
  rms_free(s);
  return 0;
}

---

You have to be careful that the the maximum number of values, at their maximum ranges to not accumulate into a value larger than the total can hold. (.e.g this is good for calculating the RMS of 127 12-bit values as 4095*4095*127 is just under 2^31. You could make it better by using an unsigned total, and should really use fixed-size integers (eg. int32_t).

If you use a floating point type for the total, then you are in for a world of hurt with accumulated errors.

isnt rms 0.707ie sine of 45 degrees

Quote:

---

Originally Posted by hamster_nz

If you don't want to calculate the RMS of all values ever seen (e.g. only the last 10 samples), and you want to do it efficently, then this is how I would do it:

Code:

---

#include <stdio.h>
#include <math.h>
#include <stdlib.h>


struct rms_state {
  int count;
  int max_values;
  int index;
  int total;
  int *data_squared;  // History of values
};


struct rms_state *rms_init(int max_values) {
  struct rms_state *s = malloc(sizeof(struct rms_state));
  if(s == NULL) {
      return NULL;
  }
  s->count = 0;
  s->max_values = max_values;
  s->total = 0.0;
  s->index = 0;
  s->data_squared = malloc(sizeof(double)*max_values);
  if(s->data_squared == NULL) {
      free(s);
      return NULL;
  }
  return s;
}


double rms_update(struct rms_state *s, int value) {
  value *= value;

  // See if we need to remove the oldest value from the running total
  if(s->count != s->max_values)
      s->count++;
  else {
      s->total -= s->data_squared[s->index];
  }

  // Add the squared new value into the running total
  s->data_squared[s->index] = value;
  s->total += value;


  s->index++;
  if(s->index == s->max_values)
      s->index = 0;


  return sqrt(s->total / s->count);
}

void rms_free(struct rms_state *s) {
  free(s->data_squared);
  free(s);
}

int main(void) {
  // Allocate state info
  struct rms_state *s = rms_init(10);
  if(s == NULL) {
      printf("Out of mem?\n");
      return 0;
  }

  // Add 100 numbers to RMS value, one at a time.
  for(int i = 0; i < 100; i++) {
      int new_value = rand()%100;
      double rms_value = rms_update(s, new_value);
      printf("%8d %8f\n", new_value, rms_value);
  }
  // Flush out the RMS buffer - should drop back to zero.
  for(int i = 0; i < 15; i++) {
      int new_value = 0;
      double rms_value = rms_update(s, new_value);
      printf("%8d %8f\n", new_value, rms_value);
  }

  // Free resources
  rms_free(s);
  return 0;
}

---

You have to be careful that the the maximum number of values, at their maximum ranges to not accumulate into a value larger than the total can hold. (.e.g this is good for calculating the RMS of 127 12-bit values as 4095*4095*127 is just under 2^31. You could make it better by using an unsigned total, and should really use fixed-size integers (eg. int32_t).

If you use a floating point type for the total, then you are in for a world of hurt with accumulated errors.

---

Thanks for the code , the controller which we am using has a hardware floating point unit , so we am thinking the accumulated errors will not be much..