psse_convert_xfmr

PSSE_CONVERT_XFMR Convert transformer data from PSS/E RAW to MATPOWER
   [XFMR, BUS, WARNINGS] = PSSE_CONVERT_XFMR(WARNINGS, TRANS2, TRANS3, ...
                                   VERBOSE, BASEMVA, BUS)
   [XFMR, BUS, WARNINGS, BUS_NAME] = PSSE_CONVERT_XFMR(WARNINGS, TRANS2, ...
                                   TRANS3, VERBOSE, BASEMVA, BUS, BUS_NAME)

   Convert all transformer data read from a PSS/E RAW data file
   into MATPOWER format. Returns a branch matrix corresponding to
   the transformers and an updated bus matrix, with additional buses
   added for the star points of three winding transformers.

   Inputs:
       WARNINGS :  cell array of strings containing accumulated
                   warning messages
       TRANS2  : matrix of raw two winding transformer data returned
                 by PSSE_READ in data.trans2.num
       TRANS3  : matrix of raw three winding transformer data returned
                 by PSSE_READ in data.trans3.num
       VERBOSE :   1 to display progress info, 0 (default) otherwise
       BASEMVA : system MVA base
       BUS     : MATPOWER bus matrix
       BUS_NAME: (optional) cell array of bus names

   Outputs:
       XFMR    : MATPOWER branch matrix of transformer data
       BUS     : updated MATPOWER bus matrix, with additional buses
                 added for star points of three winding transformers
       WARNINGS :  cell array of strings containing updated accumulated
                   warning messages
       BUS_NAME: (optional) updated cell array of bus names

   See also PSSE_CONVERT.

0001 function [xfmr, bus, warns, bus_name] = psse_convert_xfmr(warns, trans2, trans3, verbose, baseMVA, bus, bus_name)
0002 %PSSE_CONVERT_XFMR Convert transformer data from PSS/E RAW to MATPOWER
0003 %   [XFMR, BUS, WARNINGS] = PSSE_CONVERT_XFMR(WARNINGS, TRANS2, TRANS3, ...
0004 %                                   VERBOSE, BASEMVA, BUS)
0005 %   [XFMR, BUS, WARNINGS, BUS_NAME] = PSSE_CONVERT_XFMR(WARNINGS, TRANS2, ...
0006 %                                   TRANS3, VERBOSE, BASEMVA, BUS, BUS_NAME)
0007 %
0008 %   Convert all transformer data read from a PSS/E RAW data file
0009 %   into MATPOWER format. Returns a branch matrix corresponding to
0010 %   the transformers and an updated bus matrix, with additional buses
0011 %   added for the star points of three winding transformers.
0012 %
0013 %   Inputs:
0014 %       WARNINGS :  cell array of strings containing accumulated
0015 %                   warning messages
0016 %       TRANS2  : matrix of raw two winding transformer data returned
0017 %                 by PSSE_READ in data.trans2.num
0018 %       TRANS3  : matrix of raw three winding transformer data returned
0019 %                 by PSSE_READ in data.trans3.num
0020 %       VERBOSE :   1 to display progress info, 0 (default) otherwise
0021 %       BASEMVA : system MVA base
0022 %       BUS     : MATPOWER bus matrix
0023 %       BUS_NAME: (optional) cell array of bus names
0024 %
0025 %   Outputs:
0026 %       XFMR    : MATPOWER branch matrix of transformer data
0027 %       BUS     : updated MATPOWER bus matrix, with additional buses
0028 %                 added for star points of three winding transformers
0029 %       WARNINGS :  cell array of strings containing updated accumulated
0030 %                   warning messages
0031 %       BUS_NAME: (optional) updated cell array of bus names
0032 %
0033 %   See also PSSE_CONVERT.
0034 
0035 %   MATPOWER
0036 %   Copyright (c) 2014-2016, Power Systems Engineering Research Center (PSERC)
0037 %   by Yujia Zhu, PSERC ASU
0038 %   and Ray Zimmerman, PSERC Cornell
0039 %   Based on mptransin.m and mptransficbus.m, written by:
0040 %       Yujia Zhu, Jan 2014, yzhu54@asu.edu.
0041 %
0042 %   This file is part of MATPOWER.
0043 %   Covered by the 3-clause BSD License (see LICENSE file for details).
0044 %   See https://matpower.org for more info.
0045 
0046 %% define named indices into bus, gen, branch matrices
0047 [PQ, PV, REF, NONE, BUS_I, BUS_TYPE, PD, QD, GS, BS, BUS_AREA, VM, ...
0048     VA, BASE_KV, ZONE, VMAX, VMIN, LAM_P, LAM_Q, MU_VMAX, MU_VMIN] = idx_bus;
0049 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
0050     TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
0051     ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
0052 
0053 %% compatibility
0054 use_winding_baseV = 1;  %% If true, will use winding base V for any conversions
0055                         %% instead of bus base V if they are different.
0056                         %% Turn off for compatibility with, apparently buggy,
0057                         %% PSS/E implementation (checked by YZ both for
0058                         %% PSS/E internal model and conversion from RAW to
0059                         %% IEEE format)
0060 
0061 %% sizes
0062 nb = size(bus, 1);
0063 nt2 = size(trans2, 1);
0064 nt3 = size(trans3, 1);
0065 
0066 %%-----  create fictitious buses for star point of 3 winding transformers  -----
0067 starbus = zeros(nt3, VMIN);     %% initialize additional bus matrix
0068 if nt3 > 0
0069     mb =  max(bus(:, BUS_I));       %% find maximum bus number
0070     b = 10^ceil(log10(mb+1));       %% start new numbers at next order of magnitude
0071     wind1bus = trans3(:,1);         %% winding1 bus number
0072     e2i = sparse(bus(:, BUS_I), ones(nb, 1), 1:nb, max(bus(:, BUS_I)), 1);
0073     starbus(:, BUS_I) = (b+1:b+nt3)';   %% bus numbers follow originals
0074     starbus(:, BUS_TYPE) = PQ;
0075     starbus(trans3(:,12)==0, BUS_TYPE) = NONE;  %% isolated if transformer is off-line
0076     starbus(:, VA) = trans3(:,31);  %% VA = star point voltage angle
0077     starbus(:, VM) = trans3(:,30);  %% VM = star point voltage magnitude (PU)
0078     starbus(:, [BUS_AREA, ZONE]) = bus(e2i(wind1bus), [7,11]); %% wind1 bus area, zone
0079     starbus(:, BASE_KV) = 1;        %% baseKV = 1 kV (RayZ: why?)
0080     starbus(:, VMAX) = 1.1;
0081     starbus(:, VMIN) = 0.9;
0082 end
0083 
0084 %% two winding transformer data
0085 [tf,fbus] = ismember(trans2(:,1), bus(:, BUS_I));      %% I
0086 [tf,tbus] = ismember(trans2(:,2), bus(:, BUS_I));      %% J
0087 %% check for bad bus numbers
0088 k = find(fbus == 0 | tbus == 0);
0089 if ~isempty(k)
0090     warns{end+1} = sprintf('Ignoring %d two-winding transformers with bad bus numbers', length(k));
0091     if verbose
0092         fprintf('WARNING: Ignoring %d two-winding transformers with bad bus numbers', length(k));
0093     end
0094     fbus(k) = [];
0095     tbus(k) = [];
0096     trans2(k, :) = [];
0097     nt2 = nt2 - length(k);
0098 end
0099 if use_winding_baseV
0100     Zbs = bus(:, BASE_KV).^2 / baseMVA;     %% system impedance base
0101 end
0102 if nt2 > 0
0103     cw2 = find(trans2(:,5) == 2);   %% CW = 2
0104     cw3 = find(trans2(:,5) == 3);   %% CW = 3
0105     cw23 = [cw2;cw3];               %% CW = 2 or 3
0106     cz2 = find(trans2(:,6) == 2);   %% CZ = 2
0107     cz3 = find(trans2(:,6) == 3);   %% CZ = 3
0108     cz23 = [cz2;cz3];               %% CZ = 2 or 3
0109 
0110     R = trans2(:,21);
0111     X = trans2(:,22);
0112     if use_winding_baseV
0113         Zb = ones(nt2, 1);
0114         Zb(cz23) = trans2(cz23,25).^2 ./ trans2(cz23,23);
0115     end
0116     R(cz3) = 1e-6 * R(cz3, 1) ./ trans2(cz3,23);
0117     X(cz3) = sqrt(X(cz3).^2 - R(cz3).^2);   %% R, X for cz3, pu on winding bases
0118     if use_winding_baseV
0119         R(cz23) = R(cz23) .* Zb(cz23) ./ Zbs(fbus(cz23));
0120         X(cz23) = X(cz23) .* Zb(cz23) ./ Zbs(fbus(cz23));
0121     else    %% use bus base V (even if winding base V is different)
0122         R(cz23) = baseMVA * R(cz23, 1) ./ trans2(cz23,23);
0123         X(cz23) = baseMVA * X(cz23, 1) ./ trans2(cz23,23);
0124     end
0125     tap = trans2(:,24) ./ trans2(:,40);                 %% WINDV1/WINDV2
0126     tap(cw23) = tap(cw23, 1) .* bus(tbus(cw23), BASE_KV)./bus(fbus(cw23), BASE_KV);
0127     tap(cw3)  = tap(cw3, 1)  .* trans2(cw3,25)./trans2(cw3,41);
0128     shift = trans2(:, 26);
0129 end
0130 
0131 %% three winding transformer data
0132 [tf,ind1] = ismember(trans3(:,1), bus(:, BUS_I));
0133 [tf,ind2] = ismember(trans3(:,2), bus(:, BUS_I));
0134 [tf,ind3] = ismember(trans3(:,3), bus(:, BUS_I));
0135 %% check for bad bus numbers
0136 k = find(ind1 == 0 | ind2 == 0 | ind3 == 0);
0137 if ~isempty(k)
0138     warns{end+1} = sprintf('Ignoring %d three-winding transformers with bad bus numbers', length(k));
0139     if verbose
0140         fprintf('WARNING: Ignoring %d three-winding transformers with bad bus numbers', length(k));
0141     end
0142     ind1(k) = [];
0143     ind2(k) = [];
0144     ind3(k) = [];
0145     trans3(k, :) = [];
0146     starbus(k, :) = [];
0147     nt3 = nt3 - length(k);
0148 end
0149 %% Each three winding transformer will be converted into 3 branches:
0150 %% The branches will be in the order of
0151 %% # winding1 -> # winding2
0152 %% # winding2 -> # winding3
0153 %% # winding3 -> # winding1
0154 if nt3 > 0
0155     cw2 = find(trans3(:,5) == 2);   %% CW = 2
0156     cw3 = find(trans3(:,5) == 3);   %% CW = 3
0157     cw23 = [cw2;cw3];               %% CW = 2 or 3
0158     cz2 = find(trans3(:,6) == 2);   %% CZ = 2
0159     cz3 = find(trans3(:,6) == 3);   %% CZ = 3
0160     cz23 = [cz2;cz3];               %% CZ = 2 or 3
0161 
0162     tap1 = trans3(:, 32);   %% off nominal tap ratio of branch 1
0163     tap2 = trans3(:, 48);   %% off nominal tap ratio of branch 2
0164     tap3 = trans3(:, 64);   %% off nominal tap ratio of branch 3
0165     tap1(cw23) = tap1(cw23, 1) ./ bus(ind1(cw23), BASE_KV);
0166     tap2(cw23) = tap2(cw23, 1) ./ bus(ind2(cw23), BASE_KV);
0167     tap3(cw23) = tap3(cw23, 1) ./ bus(ind3(cw23), BASE_KV);
0168     tap1(cw3)  = tap1(cw3, 1)  .* trans3(cw3, 33);
0169     tap2(cw3)  = tap2(cw3, 1)  .* trans3(cw3, 49);
0170     tap3(cw3)  = tap3(cw3, 1)  .* trans3(cw3, 65);
0171     shift1 = trans3(:, 34);
0172     shift2 = trans3(:, 50);
0173     shift3 = trans3(:, 66);
0174 
0175     %% replace winding base voltage with bus base voltage
0176     % commented out: Yujia thinks this is wrong
0177     % trans3(cz3, 33) = bus(ind1(cz3), BASE_KV);
0178     % trans3(cz3, 49) = bus(ind1(cz3), BASE_KV);
0179     % trans3(cz3, 65) = bus(ind1(cz3), BASE_KV);
0180 
0181     R12 = trans3(:, 21);
0182     X12 = trans3(:, 22);
0183     R23 = trans3(:, 24);
0184     X23 = trans3(:, 25);
0185     R31 = trans3(:, 27);
0186     X31 = trans3(:, 28);
0187     Zb1 = trans3(:, 33).^2 ./ trans3(:, 23);
0188     Zb2 = trans3(:, 49).^2 ./ trans3(:, 26);
0189     Zb3 = trans3(:, 65).^2 ./ trans3(:, 29);
0190 
0191     R12(cz3) = 1e-6 * R12(cz3, 1) ./ trans3(cz3,23);
0192     X12(cz3) = sqrt(X12(cz3).^2 - R12(cz3).^2);
0193     R23(cz3) = 1e-6 * R23(cz3, 1) ./ trans3(cz3,26);
0194     X23(cz3) = sqrt(X23(cz3).^2 - R23(cz3).^2);
0195     R31(cz3) = 1e-6 * R31(cz3, 1) ./ trans3(cz3,29);
0196     X31(cz3) = sqrt(X31(cz3).^2 - R31(cz3).^2);
0197 
0198     if use_winding_baseV
0199         R12(cz23) = R12(cz23) .* Zb1(cz23) ./ Zbs(ind1(cz23));
0200         X12(cz23) = X12(cz23) .* Zb1(cz23) ./ Zbs(ind1(cz23));
0201         R23(cz23) = R23(cz23) .* Zb2(cz23) ./ Zbs(ind2(cz23));
0202         X23(cz23) = X23(cz23) .* Zb2(cz23) ./ Zbs(ind2(cz23));
0203         R31(cz23) = R31(cz23) .* Zb3(cz23) ./ Zbs(ind3(cz23));
0204         X31(cz23) = X31(cz23) .* Zb3(cz23) ./ Zbs(ind3(cz23));
0205     else    %% use bus base V (even if winding base V is different)
0206         R12(cz23) = baseMVA * R12(cz23, 1) ./ trans3(cz23,23);
0207         X12(cz23) = baseMVA * X12(cz23, 1) ./ trans3(cz23,23);
0208         R23(cz23) = baseMVA * R23(cz23, 1) ./ trans3(cz23,26);
0209         X23(cz23) = baseMVA * X23(cz23, 1) ./ trans3(cz23,26);
0210         R31(cz23) = baseMVA * R31(cz23, 1) ./ trans3(cz23,29);
0211         X31(cz23) = baseMVA * X31(cz23, 1) ./ trans3(cz23,29);
0212     end
0213 
0214     R1 = (R12+R31-R23) ./ 2;
0215     R2 = (R12+R23-R31) ./ 2;
0216     R3 = (R31+R23-R12) ./ 2;
0217     X1 = (X12+X31-X23) ./ 2;
0218     X2 = (X12+X23-X31) ./ 2;
0219     X3 = (X31+X23-X12) ./ 2;
0220 end
0221 
0222 %%-----  assemble transformer data into MATPOWER branch format  -----
0223 %%% two winding transformers %%%
0224 xfmr2 = zeros(nt2, ANGMAX);
0225 if nt2 > 0
0226     xfmr2(:, [F_BUS T_BUS]) = trans2(:,[1,2]);
0227     xfmr2(:, [BR_R BR_X TAP SHIFT]) = [R X tap shift];
0228     xfmr2(:, [RATE_A RATE_B RATE_C]) = trans2(:,[27,28,29]);
0229     xfmr2(:, BR_STATUS) = trans2(:,12);
0230     xfmr2(:, ANGMIN) = -360;
0231     xfmr2(:, ANGMAX) = 360;
0232 end
0233 
0234 %%% three winding transformers %%%
0235 xfmr3 = zeros(3*nt3, ANGMAX);
0236 if nt3 > 0
0237     idx1 = (1:3:3*nt3)';        %% indices of winding 1
0238     idx2 = idx1+1;              %% indices of winding 2
0239     idx3 = idx1+2;              %% indices of winding 3
0240     %% bus numbers
0241     xfmr3(idx1, [F_BUS T_BUS]) = [trans3(:,1), starbus(:,1)];
0242     xfmr3(idx2, [F_BUS T_BUS]) = [trans3(:,2), starbus(:,1)];
0243     xfmr3(idx3, [F_BUS T_BUS]) = [trans3(:,3), starbus(:,1)];
0244     %% impedances, tap ratios & phase shifts
0245     xfmr3(idx1, [BR_R BR_X TAP SHIFT]) = [R1 X1 tap1 shift1];
0246     xfmr3(idx2, [BR_R BR_X TAP SHIFT]) = [R2 X2 tap2 shift2];
0247     xfmr3(idx3, [BR_R BR_X TAP SHIFT]) = [R3 X3 tap3 shift3];
0248     %% ratings
0249     xfmr3(idx1, [RATE_A RATE_B RATE_C]) = trans3(:,[35,36,37]);
0250     xfmr3(idx2, [RATE_A RATE_B RATE_C]) = trans3(:,[51,52,53]);
0251     xfmr3(idx3, [RATE_A RATE_B RATE_C]) = trans3(:,[67,68,69]);
0252     xfmr3(:, ANGMIN) = -360;        %% angle limits
0253     xfmr3(:, ANGMAX) =  360;
0254     %% winding status
0255     xfmr3(:, BR_STATUS) = 1;        %% initialize to all in-service
0256     status = trans3(:, 12);
0257     k1 = find(status == 0 | status == 4);   %% winding 1 out-of-service
0258     k2 = find(status == 0 | status == 2);   %% winding 2 out-of-service
0259     k3 = find(status == 0 | status == 3);   %% winding 3 out-of-service
0260     if ~isempty(k1)
0261         xfmr3(idx1(k1), BR_STATUS) = 0;
0262     end
0263     if ~isempty(k2)
0264         xfmr3(idx2(k2), BR_STATUS) = 0;
0265     end
0266     if ~isempty(k3)
0267         xfmr3(idx3(k3), BR_STATUS) = 0;
0268     end
0269 end
0270 
0271 %% combine 2-winding and 3-winding transformer data
0272 xfmr = [xfmr2; xfmr3];
0273 
0274 % %% delete out-of-service windings
0275 % k = find(xfmr(:, BR_STATUS) == 0);
0276 % xfmr(k, :) = [];
0277 % k = find(trans3(:,12) <= 0);
0278 % starbus(k, :) = [];
0279 
0280 %% finish adding the star point bus
0281 bus = [bus; starbus];
0282 if nt3 > 0
0283     warns{end+1} = sprintf('Added buses %d-%d as star-points for 3-winding transformers.', ...
0284         starbus(1, BUS_I), starbus(end, BUS_I));
0285     if verbose
0286         fprintf('Added buses %d-%d as star-points for 3-winding transformers.\n', ...
0287             starbus(1, BUS_I), starbus(end, BUS_I));
0288     end
0289 end
0290 if nargin > 6 && nargout > 3
0291     starbus_name = cell(nt3, 1);
0292     for k = 1:nt3
0293         starbus_name{k} = sprintf('STR_PT_XF_%-2d', k);
0294     end
0295     bus_name = [bus_name; starbus_name];
0296 end