Condition Monitoring with frequency analysis
This tutorial configures a complete monitoring application, based on the TwinCAT 3 Condition Monitoring API. It facilitates creation of a workflow for Condition Monitoring applications, including data collection and adding high-performance analysis algorithms. The block diagram below illustrates the arrangement of the application. For a better understanding of the programming tasks, the document is subdivided into design steps.
The sample is available for download from here:
FrequencyAnalysis_Tutorial.zip
It can be modified and extended as required.
Block diagram
Step 1: Application specification
The first step for the design of a condition monitoring application is to determine the main aims of the application, e.g. automatic warning in the event of excessive vibrations or in the event of a malfunction in the bearing, based on frequency analysis. It is also important to consider other technical factors such as measuring sensors, the sampling rate of the controller and the expected accuracy. The aim of this tutorial is to detect small and large errors in the input signal with the aid of the magnitude spectrum and its quantile distribution. In addition, a classifier is used for predicting the general state as "normal state", "warning state" or "alarm state". The table below shows a list of the function blocks used in this tutorial.
Function block |
---|
FB_CMA_Souce |
FB_CMA_Sink |
FB_CMA_MagnitudeSpectrum |
FB_CMA_Quantiles |
FB_CMA_DiscreteClassification |
For a more detailed description of the algorithm selection for specific issues such as bearing analysis, gear unit analysis or frequency analysis, we refer to the solutions described elsewhere. Since the aim of the tutorial is to detect general changes in the input signal, a magnitude spectrum with a resolution of 4096 lines is sufficient. The 50 % and 90 % quantile of the spectral values are calculated, and the result is classified as "normal state", “warning state” or "alarm state“.
Step 2: Configuration of the PLC tasks
Since condition monitoring and analysis is comprised of a data acquisition stage, a calculation stage and an analysis stage, the task has to be structured according to the calculation requirements for each step. You can find additional information on this topic at the task setting. The aim of this tutorial is to calculate the magnitude spectrum of 4096 lines, for which approx. 3200 data samples are required. The means that, during the data collection stage, a source multi-array has to handle 1600 data samples, considering overlapping. With 10x oversampling, 160 cycles are required, or 160 ms with 1 ms trigger, to fill a single multi-array. The following setting is therefore recommended for the calculation task:
Calculation cycle time < (data collection cycle time * buffer size / oversampling factor)*0.5
For the tutorial the calculation cycle time is set to 10 ms. For the actual application it is important to consider the computation load, which is affected by other tasks that run simultaneously on the same controller, such as visualization or network communication. Further information on task settings can be found here in the Task cycle time section. Make sure that adequate router memory capacity is allocated before starting to build a Condition Monitoring application. This tutorial was set for working with a router memory capacity adjusted to 32 MB.
Step 3: Configuration of the function blocks
In this step the function blocks listed above are configured according to the application requirements. As already mentioned, the source multi-array collects 1600 data samples for calculating a spectrum. The aDimSizes parameter is therefore set to 1600. Since the tutorial only considers one channel, nDims is set to 1.
PROGRAM CM_Worker
VAR CONSTANT
cInitSourceSpectrum : ST_MA_MultiArray_InitPars := ( eTypeCode := eMA_TypeCode_LREAL, nDims := 1, aDimSizes := [1600]);
END_VAR
In the calculation task the magnitude spectrum for calculating a spectrum of 4096 lines is configured, indicated by cFFTLength. A so-called window function is used, since the spectrum calculation is associated with periodic processing of discrete segments of a continuous signal. A correctly selected window function improves the signal transformation efficiency, reduces fluctuations thanks to the overlap-add method and improves the spectral resolution. In practical applications the window function also reduces the leakage effect near critical frequencies. In the tutorial Hann window was selected. The magnitude spectrum function block offers a wide range of scaling options, out of which the RMS value was selected. The reason is that for time-varying physical signals, an RMS value is a preferred indicator of the mean signal power, in contrast to the peak value, for example. In the vibration acceleration spectrum, individual lines indicate the effective values of the vibrations at the corresponding frequency and can be expressed directly in the corresponding units such as mm/s² or g.
PROGRAM MAIN_CM
VAR CONSTANT
cInitSpectrum : ST_CM_MagnitudeSpectrum_InitPars := ( nFFT_Length := 4096,
nWindowLength := 2*1600,
bTransformToDecibel:= FALSE,
eWindowType := eCM_HannWindow,
eScalingType := eCM_RMS);
END_VAR
The result of the magnitude spectrum is copied to an array via a sink function block, with specified array length of nFFT_Length/2+1. In the next step of the analysis chain, a quantile function block for calculating the 50 % and 90 % quantiles of the spectral values is configured. In many cases the spectral values fluctuate strongly, so that an evaluation is difficult if the values are too low or too high. Using the quantiles it is possible to determine the maximum, minimum or indeed the average value over a specified time interval. This type of range-based evaluation is often more reliable and easier to handle. A 50% quantile value Q0.5 is a value for which almost 50% of the values of a distribution are less than Q0.5. It is also referred to as median value. Similarly, a 90% quantile Q0.9 indicates a value for which 90% of the values of a distribution are less than Q0.9.
VAR CONSTANT
cInitQuantiles : ST_CM_Quantiles_InitPars := ( nChannels := (4096/2+1),
fMinBinned := -10,
fMaxBinned := 10,
nBins := 100,
nMaxQuantiles := 2);
END_VAR
In the program the quantiles are configured as follows:
(*--------- Configure quantile args ---------*)
IF bConfigureQuantile THEN
FOR nChannel := 1 TO (cFFTLength/2+1) DO
aQuantilesArg[nChannel,1]:= 0.50; // 50% quantile
aQuantilesArg[nChannel,2]:= 0.90; // 90% quantile
END_FOR
fbQuantiles.Configure( pArg := ADR(aQuantilesArg), nArgSize := SIZEOF(aQuantilesArg));
bConfigureQuantile := FALSE;
END_IF
Here you can find a more detailed description of the function block FB_CMA_Quantiles. Note that the parameters fMinBinned and fMaxBinned define the expected input signal range and nBins indicates the number of Bins into which the signal range is divided. These parameters depend on the respective input signal. The signal state is classified based on the quantiles information. The discrete function block can process several channels simultaneously, therefore the quantile output is sent directly to the block. The classifier is set to distinguish between three states and to display the corresponding state via the nMaxClasses parameter.
VAR CONSTANT
cInitClassification : ST_CM_DiscreteClassification_InitPars := (nChannels:= (4096/2+1),
nMaxClasses := 3);
END_VAR
The output of the quantiles function block is a 2D array, which in this case is the number of spectral lines over the number of quantiles. But the discrete classifier only allows a one-dimensional array, which contains the number of input channels. In order to avoid a dimension conflict, the buffer converter of FB_CMA_BufferConverting should be used. This function block converts a two-dimensional multi-array to a one-dimensional array without any data loss. The code snippet describes the corresponding application.
VAR CONSTANT
cInitBuffer : ST_MA_MultiArray_InitPars := ( eTypeCode := eMA_TypeCode_LREAL,
nDims := 1,
aDimSizes := [(4096/2+1)]);
END_VAR
VAR
fbBufferConverter : FB_CMA_BufferConverting := (stInitPars := cInitBuffer, nOwnID := eID_BufferConverter, aDestIDs := [eID_Classify]);
END_VAR
The buffer converter calls a method:
fbBufferConverter.Copy1D(nWorkDimIn := 0,
nWorkDimOut := 0,
nElements := 0,
pStartIndexIn := 0,
pStartIndexOut := 0,
nOptionPars := 0);
Further information on this function block can be found under FB_CMA_BufferConverting. To complete the function block configuration, each sink function block must be linked to PLC arrays with correct dimensions.
Step 4: Fine-tuning of the application parameters
Before starting the analysis, it is important to configure the discrete classifier with regard to its limit values. A classification limit or threshold value enables the discrete classifier to monitor incoming channels continuously and determine whether one of the input channels exceeds this threshold value. The threshold values depend on the respective application, the accuracy requirements, the permitted detection tolerances, etc. The aim of this tutorial is to detect small errors, which are comparable to random noise, and also large-sized errors, which occur at a specific frequency (e.g. 200 Hz). The threshold values fWarning and fAlarm are determined. If the amplitude of the input channels exceeds fWarning, the general state switches to warning state. If fAlarm is exceeded, an alarm message is issued. If the threshold values are not exceeded, the channel state is in the normal range.
(*--------- Configure classifier args ---------*)
IF bConfigureClassifier THEN
fWarning := (fMonitoringLevel/100)*1.5;
fAlarm := (fMonitoringLevel/100)*2.5;
fbTeachTimer(IN := TRUE, PT := T#15S);
IF fbTeachTimer.Q THEN
fbTeachTimer(IN := FALSE);
FOR nChannel := 1 TO (cFFTLength/2+1) DO
aClassArgs[nChannel, 1] := (fMonitoringLevel/100)*aQuantilesCopy[nChannel,1];
aClassArgs[nChannel, 2] := fWarning*aQuantilesCopy[nChannel,1];
aClassArgs[nChannel, 3] := fAlarm*aQuantilesCopy[nChannel,1];
END_FOR
fbClassification.Configure(pArg := ADR(aClassArgs), nArgSize := SIZEOF(aClassArgs));
bConfigureClassifier := FALSE;
END_IF
END_IF
The code snippet above describes the configuration of this discrete classifier, so that a timer block allows a normal operating window to pass through a so-called teaching phase, during which the discrete classifier is configured. It is assumed that the input signal behaves normally during this time, i.e. within the permissible range. The warning threshold is 150 % of the “normal” 50 % quantile, the alarm threshold is 250 % of the normal 50 % quantile. Since the 50% quantile describes the average behavior, this threshold value configuration is suitable for applications whose inputs only have few outliers. The 90 % quantile can also be determined as threshold value, if it is assumed that the input signal is likely to fluctuate strongly. It is also possible to configure another variable, fMonitoringLevel, which can be used to apply a certain tolerance range around the permissible value, in order to control the number of false alarms. This parameter can be used to fine-tune the threshold values. Note that the threshold values for the discrete classifier can be specified individually for all input channels.
Step 5: Starting the application
Compile the code, download it to the target system and start the PLC, in order to execute the tutorial. A small prepared visualization, referred to as Dashboard, can be found in the Solution Explorer under the VISU node, which can be used for a quick test. For the simulation the input signal is linked to a function generator, which was configured for generating a sinusoidal 50 Hz signal with an amplitude of 5. Other available signals such as pulse, triangle or saw tooth, or indeed a hardware module such as EL3632, can be applied to the input. Once the application has been started, the display fields show the maximum amplitude, the RMS amplitude and the frequency at the maximum amplitude of the PLC in real-time.
The diagram illustrates that the state of the application is shown in the corresponding display field. A small-sized error can be simulated by pressing the Add Fault button. You can see how the RMS value of the input signal slowly increases beyond the threshold value and how the state changes accordingly. To simulate a large-sized error, press the Small/Large button. Similar to the previous error the RMS value will increase, but this time the “Fault Frequency” field shows the frequency of the fault signal, in this case 200 Hz.
Step 6: Monitoring
Once the PLC has started, the display fields show the current values. Initially the Reconfigure button is shown as pressed, and the signal in the right-hand corner is disabled. The means that the limit values are in the process of being specified for the discrete classifier. Once the configuration is complete, the Reconfigure button resets itself, and the machine status is shown as “normal state”. The signal switches to green, which has the same meaning.
To simulate an error, leave the option field at “Small Fault” and press the Activate Disturbance button. The machine will switch between “Normal” and “Warning” state, and the signal switches between green and orange. If a large error is simulated by switching the option field, the machine state switches to “Alarm” state, and the signal switches to red. To prevent the fault, release the Activate Disturbance button. The signal state returns to green. Note that a change in the signal amplitude also results in an error state. If this is undesirable, press the Reconfigure button again, in order to adjust the discrete classifier to this new signal state.
Requirements
Development environment |
Target system type |
PLC libraries to include |
---|---|---|
TwinCAT v3.1.4013 |
PC or CX (x86, x64) |
Tc3_CM, Tc3_CM_Base, Tc3_MultiArray |