function mcb_getCharacteristics(pmsm,inverter)
% mcb_getCharacteristics Computes the torque-speed characteristics for PMSM
% 
%   mcb_getCharacteristics calculates the torque and speed characteristics 
%   of the motor. This is helpful in control algorithm development for the motor. 
% 
%   mcb_getCharacteristics(pmsm,inverter) returns the following
%   characteristics the given PMSM:
%   - Torque vs Speed
%   - Power vs Speed 
%   - Idq vs Speed
%   - Ipeak vs Speed
%   
%   The function accepts a pmsm object and an inverter object as inputs. 
%   pmsm        = obj1;
%   inverter    = obj2;
%   I_max       = obj3;
% 
%   See HELP for additional information for this utility.

% Function      :   mcb_getCharacteristics
% Description   :   Compute Torque-Speed & Other curves for PMSM
% File name     :   mcb_getCharacteristics.m
% Inputs        :   pmsm object & inverter object
% Authors       :   Darshan Pandit

% Copyright 2020 The MathWorks, Inc.

%% Get Torque & Speed Characteristics
P_max = 2*pi*pmsm.N_base*pmsm.T_rated/60;
I_max = inverter.ISenseMax;

N = 1:pmsm.N_max;
T_max  = P_max./((2*pi/60).*N);
T_max(1:pmsm.N_base) = pmsm.T_rated*ones(1,pmsm.N_base);
P =  2*pi.*N.*T_max/60;

subplot(2,2,1);
plot(N,T_max,'b');
title(['Torque-Speed Characteristics for ' pmsm.model]);
xlabel('Speed (rpm)');
ylabel('Torque (Nm)');
grid on;

%% Get Power & Speed Characteristics

subplot(2,2,2)
plot(N,P,'g');
title(['Power-Speed Characteristics for ' pmsm.model]);
xlabel('Speed (rpm)');
ylabel('Power (Watts)');
grid on;

%% Get Id & Iq Characteristics

% Get Id & Iq Values
Omega_e     = pmsm.p*2*pi*N/60;
Omega_e_base= pmsm.p*2*pi*pmsm.N_base/60;

Id          = zeros(1,pmsm.N_max);  % Define size of variable
Iq          = ones(1,pmsm.N_max);   % Define size of variable

% Get Id for FWC
Id          = (((Omega_e_base - Omega_e).*pmsm.FluxPM)./(Omega_e.*pmsm.Ld));
Id(1:pmsm.N_base) = zeros(1,pmsm.N_base);

% Get Iq for FWC
index       = length(find(abs(Id) <= I_max));

max_speed_I    = N(index);  
Iq          = sqrt(I_max^2 - Id.^2);
Iq_limits   = T_max./((3/2)*pmsm.p*pmsm.FluxPM);

Rs          = pmsm.Rs;

Vd          = Id(1:index).*Rs - Omega_e(1:index).*pmsm.Lq.*Iq(1:index);
Vq          = Iq(1:index).*Rs + Omega_e(1:index).*pmsm.Ld.*Id(1:index) + Omega_e(1:index)*pmsm.FluxPM;

V_max       = sqrt(Vd.^2 + Vq.^2);
Vpu         = V_max/(inverter.V_dc/sqrt(3));
idx         = find(Vpu<=1);
max_speed_V  = N(length(idx));
I_peak       = sqrt(Id.^2 + Iq_limits.^2);

% Plot the currents Id Iq
subplot(2,2,3)
plot(N(1:index),Id(1:index),'g',N(1:index),Iq_limits(1:index),'b');
title('I_{dq} Vs Speed');
xlabel('Speed (rpm)');
ylabel('I_{dq} (Amps)');
legend('I_d','I_q');
grid on;
scale = axis;
axis([scale(1:2) -I_max I_max]);

% Plot the currents I_peak
subplot(2,2,4)
plot(N(1:index),I_peak(1:index),'k',N(1:index),pmsm.I_rated*ones(1,index),'r-.');
title('I_{peak} Vs Speed');
xlabel('Speed (rpm)');
ylabel('I_{peak} (Amps)');
legend('I_{peak}','I_{rated}');
grid on;
scale = axis;
axis([scale(1:2) -I_max I_max]);

% Scale axis
subplot(2,2,1);
scale2 = axis;
axis([scale(1:2) 0 1.1*scale2(4)]);

subplot(2,2,2);
scale2 = axis;
axis([scale(1:2) 0 1.1*scale2(4)]);

% Return the computed values
% Id_ref = Id(1:index);
% Iq_ref = Iq(1:index);
% Speed_ref = N(1:index);

end