Preprocess Induction Motor Data and Autogenerate Current Controller Calibration Tables

This example shows how to preprocess raw induction motor data generated by Finite Element Analysis and fit point-by-point models to generate optimal calibrations for an induction motor. First, rearrange the data into the Ids-Iqs coordinate system, define the speed operating points, and create a torque contour. Next, store the data in a table object. Then, use Model-Based Calibration Toolbox™ (MBC) functions to fit the point-by-point models for the stator q-axis current, Iq, and the modulation voltage, Vs, using the processed data. Finally, generate current controller calibrations for the induction motor.

You can also import the table object with the preprocessed data to the MBC Model Fitting app to fit models, and generate an optimal calibration using the MBC Optimization app.

Import Raw Data

Import the raw induction motor data. This example uses a data set based on a joint project between MathWorks® and Motor-CAD™ References. For more information, see Import Induction Machine (Squirrel Cage) Flux Linkage Data from Motor-CAD (Simscape Electrical).

dataFile = 'IMFluxMotorCADData.mat';
load(dataFile);

Visualize Raw Data

Visualize the flux linkages with respect to the d-axis and q-axis currents.

stem3(Isd_Lambda_mat,Isq_Lambda_mat,LambdaSD_Lambda_mat, ...
    "LineStyle","none");
title("Flux d [Wb]");
xlabel("Ids [A]");
ylabel("Iqs [A]");

stem3(Isd_Lambda_mat,Isq_Lambda_mat,LambdaSQ_Lambda_mat, ...
    "LineStyle","none");
title("Flux q [Wb]");
xlabel("Ids [A]");
ylabel("Iqs [A]");

Visualize the torque data with respect to the d-axis and q-axis currents.

stem3(Isd_Lambda_mat,Isq_Lambda_mat,Tdq_mat, ...
    "LineStyle","none");
title("Measured Torque [Nm]");
xlabel("Ids [A]");
ylabel("Iqs [A]");

Define Induction Motor Specification

Specify these parameters to align with your induction motor specification:

• DC voltage, Vdc

• Number of pole pairs, pp

• Stator leakage inductance, Lls

• Rotor leakage inductance, Llr

• Stator phase resistance, Rs

• Rotor phase resistance, Rr

• Phase peak current, IsMax

• Maximum modulation voltage, VsMax

Vdc = 400;
pp = 3;    
Lls = 0.0231;
Llr = 0.0403;
Rs = 0.0094;
Rr = mean(RrdVec);
IsMax = 500*sqrt(2);
VsMax = Vdc*0.95/sqrt(3);

Calculate Gridded Data for Id and Iq Axis Currents

Calculate the gridded data, [IDgrid,IQgrid], using the ndgrid function to create a grid for idAxis and iqAxis, bound by IsMax. The software uses gridded data to calculate torque in a later step and create the torque contour plot.

idAxis = linspace(0,IsMax,51); 
iqAxis = linspace(0,IsMax,52); 
[IDgrid,IQgrid] = ndgrid(idAxis,iqAxis);

Define the breakpoints in a cell array.

currentBreakpoints = {idAxis,iqAxis};

Transform Raw Data into Columns

Transform the raw data to column-based data so you can then use in the fitlookupn function.

Note that Fr_col is the rotor electrical frequency in Hz.

id_col = reshape(Isd_Lambda_mat,[numel(Isd_Lambda_mat),1]);
iq_col = reshape(Isq_Lambda_mat,[numel(Isq_Lambda_mat),1]);
flux_d_col = reshape(LambdaSD_Lambda_mat,[numel(LambdaSD_Lambda_mat),1]);
flux_q_col = reshape(LambdaSQ_Lambda_mat,[numel(LambdaSQ_Lambda_mat),1]);
trq_col = reshape(Tdq_mat,[numel(Tdq_mat),1]);
Lphi_col = reshape(Lphi_mat,[numel(Lphi_mat),1]);
Lt_col = reshape(Lt_mat,[numel(Lt_mat),1]);
Fr_col = repmat(fr_vec_fine',[numel(Is_vec_fine),1]); 

Transform Column Data to Cartesian Coordinate System

Use the fitlookupn function to transform the column data from polar coordinates to the Cartesian coordinate system for these parameters:

• Slip frequency, SlipFreq

• d-axis flux linkage, Lambdad

• q-axis flux linkage, Lambdaq

• Measured Torque, Trq

• Magnetizing Inductance, Lphi

• Transient Inductance, Lt

For each transformation, use the current breakpoints defined in the previous section, the current column data for id and iq, and the corresponding column data for the parameter being transformed.

Transform Slip Frequency Data

Transform the slip frequency data and convert the data from hertz to radians per second. Plot the result.

SlipFreq = fitlookupn( ...
    currentBreakpoints,[id_col,iq_col],Fr_col);
ws = 2*pi*SlipFreq;
surf(IDgrid,IQgrid,ws);
title("Slip Speed [rad/s]");
xlabel("Ids [A]")
ylabel("Iqs [A]");

Transform Flux Linkages Data

Transform the flux linkages data and plot the results.

Lambdad = fitlookupn( ...
    currentBreakpoints,[id_col,iq_col],flux_d_col,ShowPlot=true);
title("Flux d [Wb]");
xlabel("Ids [A]")
ylabel("Iqs [A]");

Lambdaq = fitlookupn( ...
    currentBreakpoints,[id_col,iq_col],flux_q_col,ShowPlot=true);
title("Flux Q [Wb]");
xlabel("Ids [A]")
ylabel("Iqs [A]");

Transform Torque Data

Transform the torque data and plot the result.

Trq = fitlookupn( ...
    currentBreakpoints,[id_col,iq_col],trq_col,ShowPlot=true);
title("Torque [Nm]");
xlabel("Ids [A]")
ylabel("Iqs [A]");

Transform Inductances Data

Transform the inductances data and plot the results.

Lphi = fitlookupn( ...
    currentBreakpoints,[id_col,iq_col],Lphi_col,ShowPlot=true);
title("Lphi [H]");
xlabel("Ids [A]")
ylabel("Iqs [A]");

Lt = fitlookupn( ...
    currentBreakpoints,[id_col,iq_col],Lt_col,ShowPlot=true);
title("Lt [H]");
xlabel("Ids [A]")
ylabel("Iqs [A]");

Calculate Torque Error

Calculate the torque error from the calculated torque and plot the results. A small torque error indicates the data was transformed correctly using the fitlookupn function.

Trq_calc = 1.5*pp*(Lambdad.*IQgrid-Lambdaq.*IDgrid);
Trq_err = Trq - Trq_calc;
surf(IDgrid,IQgrid,Trq_calc);
title("Calculated Torque [Nm]");
xlabel("Ids [A]")
ylabel("Iqs [A]");

surf(IDgrid,IQgrid,Trq_err);
title("Torque Error [Nm]");
xlabel("Ids [A]")
ylabel("Iqs [A]");

Define Operating Points

Specify the number of torque operating points, numTorqueLevels, and the low torque levels, minTorqueLevels. The remaining torque levels are evenly spaced up to the maximum measured torque. Torque is specified in Nm.

numTorqueLevels = 50;
minTorqueLevels = [3 7 10 15];

Downsample the data by a factor of 5 to reduce the computational load in MBC.

downsampleFactor = 5;

Specify minimum speed speedMin, maximum speed speedMax, and speed step size speedStep, to define the speed operating points. Speed is not required to be evenly spaced, but some steps in model-based calibration will be easier with evenly spaced speed operating points. Speed is specified in rpm.

If speed is a lookup table input, then match operating points to breakpoints.

speedMin = 500;
speedMax = 12000;
speedStep = 500;

speedBreakpoints = [speedMin:speedStep:4000 ...
                    4200:200:6000 ...
                    6500:speedStep:speedMax];

Rearrange Data in Torque-Speed Operating Points

Introduce torque as an operating point in addition to speed by rearranging the data to have both torque and speed information.

Calculate Torque Contour Levels

Initialize the torque contours at the levels in minTorqueLevels. Then, evenly space the levels up to the maximum measured torque.

TrqMax = max(abs(Trq_calc),[],"all");
numPoints = numTorqueLevels-length(minTorqueLevels)+1;
torqueLevels = ...
    [minTorqueLevels(1:end-1) linspace(minTorqueLevels(end), ...
    TrqMax,numPoints)];

Create and Plot Torque Contour

Create and plot the torque contour using the contour function.

trqContour is a 2-by-N array containing sets of contour data for each contour level.

$\left[\begin{array}{cc}{Z}_i& X_i\\ N_i& Y_i\end{array}\right]$

trqContour = contour(IDgrid,IQgrid,Trq_calc,torqueLevels,"ShowText","on");
title("Torque Contour  [Nm]");
xlabel("Ids [A]");
ylabel("Iqs [A]");

Rearrange Data by Torque Levels

For each torque level, extract the Trq, Ids, and Iqs values. Next, downsample the data for each torque level. If you have the Signal Processing Toolbox™, use the downsample (Signal Processing Toolbox) function to calculate dataDownsample. Then, expand each torque level to include all speeds, attaching every speed operating point to each (ids, iqs) on the torque contour.

start = 1;
accumulatedData = []; 
while start < size(trqContour,2)
    % process each torque levels (the contour levels). 
    % There may be some contour lines missing.

    currentTorque = trqContour(1,start);

    % The number of points in the contour for this torque.
    ni = trqContour(2,start);
    stop = start + ni;
    
    torqueData = [repelem(currentTorque,ni) ; 
    trqContour(:,start+1:stop)]'; % [Id, Iq] points for current torque
    
    if false % Check to use downsample function.
        dataDownsample = downsample(torqueData,downsampleFactor);
    else
        dataDownsample = torqueData(1:downsampleFactor:end,:);
    end
    
    if ( dataDownsample(end,2)~= torqueData(end,2) ) 
        % Make sure we have the last point in contour
        torqueData = [dataDownsample;torqueData(end,:)];
    else
        torqueData = dataDownsample;
    end
    ni = size(torqueData,1);
    
    % Repeat torque contours for each speed and accumulate data.
    % Make sure speed is a column vector [Trq, Ids, Iqs, n]
    accumulatedData=[accumulatedData;
        repmat( ...
        torqueData,[length(speedBreakpoints),1])...
        repelem(speedBreakpoints(:),ni)];

    % Start of next torque contour.
    start = stop+1;
end

Create Table Object with Preprocessed Data

Store accumulatedData with the variables Trq, Id, Iq, and n to the MBCdata table object.

MBCdata = array2table(accumulatedData,VariableNames=["Trq","Id","Iq","n"]);

Calculate these derived quantities and assign them to variables in the MBCdata table object:

• Phase current, MBCdata.Is

• Slip speed, MBCdata.ws

• Synchronous speed, MBCdata.we

• Phase voltage, MBCdata.Vs

• Slip, MBCdata.slip

MBCdata.Is = sqrt(MBCdata.Id.^2 + MBCdata.Iq.^2);
MBCdata.ws = interpn(idAxis,iqAxis,ws,MBCdata.Id,MBCdata.Iq);
MBCdata.we = 2*pi.*MBCdata.n.*pp./60 + MBCdata.ws;

Fluxd = interpn(idAxis,iqAxis,Lambdad,MBCdata.Id,MBCdata.Iq);
Fluxq = interpn(idAxis,iqAxis,Lambdaq,MBCdata.Id,MBCdata.Iq);

Vd = Rs.*MBCdata.Id - (MBCdata.we).*Fluxq;  
Vq = Rs.*MBCdata.Iq + (MBCdata.we).*Fluxd; 

MBCdata.Vs = sqrt(Vd.^2 + Vq.^2);
MBCdata.slip = MBCdata.ws./MBCdata.we;

Set the units for each variable in the MBCdata table object.

MBCdata.Properties.VariableUnits = ...
    ["Nm","A","A","rpm","A","rad/s","rad/s","V" "slip"];

Induction Motor Calibration

Fit models with the preprocessed data, MBCData, and generate optimal current controller calibration tables using MBC functions.

Fit Models

Fit point-by-point models of Iq and Vs with local inputs of Iq and Trq, with speed operating points.

Create the MBC model project, proj.

proj = mbcmodel.CreateProject('IMmodels');

Create the data object, D, in the project object, proj, from the MBCdata table object.

D = CreateData(proj,MBCdata);

Edit the data object.

BeginEdit(D);

Sort the data points by speed, grouping those within a 0.5 tolerance of each other.

DefineTestGroups(D,{'n'},0.5,'',true);

Filter the data retaining current values that are within a 20% margin above the maximum value of the phase peak current, IsMax, and voltage values that are within a 20% margin above the maximum modulation voltage, VsMax.

AddFilter(D,sprintf('Is<=1.2*%g',IsMax));
AddFilter(D,sprintf('Vs<=1.2*%g',VsMax));

Apply the data object changes.

CommitEdit(D);

Create the test plan, T, and attach the data, D, to fit the models in the induction motor template, InductionMotor, provided with this example.

MBCModelTemplate = 'InductionMotor';
T = CreateTestplan(proj,MBCModelTemplate);
T.Name = 'Models';
AttachData(T,D,'UseDataRange',true,'Boundary',true);

Generate Optimal Calibration

Generate optimal calibrations with the optimizations defined in the induction motor CAGE template project, InductionMotor.cag, provided with this example.

CAGETemplate = 'InductionMotor.cag';

Create an mbcrunner object to the CAGE template project, CAGETemplate, to programmatically run CAGE operations.

runner = mbcrunner(CAGETemplate);

Define the lookup table breakpoints using the speed and torque percent breakpoints.

torquePercentBreakpoints = [1 4:4:100];
breakPoints = {speedBreakpoints, torquePercentBreakpoints};

Update the Id_Table and Iq_Table and TrqEnvelope breakpoints in the CAGE project with the look up table breakpoints, breakPoints. All tables changed are returned.

updateBreakpoints(runner,'Id_Table',breakPoints,true);

Update constraint bounds for the phase current, Is, and phase voltage, Vs, for the torque model and TPA optimizations in the CAGE project.

updateConstraints(runner,'TrqModel_Optimization','Is', ...
    IsMax,'Vs',VsMax)
updateConstraints(runner,'TPA_Optimization','Is_TrqGrid', ...
    IsMax,'Vs_TrqGrid',VsMax)

Update the models and rerun optimizations to fill the lookup tables.

importModels(runner,proj);

View Lookup Tables

Return the lookup tables from the mbcrunner object, runner.

optimalTables = runner.TableValues;

Plot the maximum Torque Envelope Table using speedBreakpoints and optimalTables.TrqEnvelope.

plot(speedBreakpoints, optimalTables.TrqEnvelope)
xlabel('Speed [rpm]')
ylabel('Torque [Nm]')
title('Torque Envelope Table')
grid

Plot the Id Table using able using torquePercentBreakpoints, speedBreakpoints, and optimalTables.Id_Table.

surf(torquePercentBreakpoints,speedBreakpoints, optimalTables.Id_Table)
ylabel('Speed [rpm]')
xlabel('Torque Percent [%]')
zlabel('Ids [A]')
title('Id Table')
grid

Plot the Iq Table using torquePercentBreakpoints, speedBreakpoints, and optimalTables.Iq_Table.

surf(torquePercentBreakpoints,speedBreakpoints, optimalTables.Iq_Table)
ylabel('Speed [rpm]')
xlabel('Torque Percent [%]')
zlabel('Iqs [A]')
title('Iq Table')

Generate an Excel® Data File for MBC

Alternatively, you can also import the MBCdata data to the MBC Model Fitting app to fit models and generate an optimal calibration using the MBC Optimization app.

If you want to write data to an Excel file to import to MBC, select the checkbox. Note that unit information will be lost.

if false
    writetable(MBCdata,"IM_MBCData.xlsx");
end

References

[1] Troncon, Diego, Matteo Carbonieri, Luigi Alberti, and Nicola Bianchi. “Measures and Simulations of Induction Machines Flux Linkage Characteristics Based on Rotor Field Orientation.” IEEE Transactions on Industry Applications 57, no. 5 (September 2021): 4686–93. https://doi.org/10.1109/TIA.2021.3089662.