cdf2matp

CDF2MATP Converts data from IEEE Common Data Format to MATPOWER format.
   CDF2MATP                         : prompts for input & output file names
   CDF2MATP(CDF_FILE, MATP_FILE)    : uses specified file names
   CDF2MATP(CDF_FILE, MATP_FILE, v) : and MATPOWER case version
                                      v = '1' or '2', default is '2'

   Optional arguments are the names of the input IEEE CDF file and the
   output MATPOWER case file. If these are not given you will be prompted
   to enter them.

   The IEEE CDF does not include some data need to run an optimal power
   flow. This script creates default values for some of this data as
   follows:

       Bus data:
           Vmin = 0.94 p.u.
           Vmax = 1.06 p.u.
       Gen data:
           Pmin = 0 MW
           Pmax = Pg + baseMVA
       Gen cost data:
           Quadratic costs with:
               c2 = 10 / Pg, c1 = 20, c0 = 0, if Pg is non-zero, and
               c2 = 0.01,    c1 = 40, c0 = 0, if Pg is zero
           This should yield an OPF solution "close" to the
           existing solution (assuming it is a solved case)
           with lambdas near $40/MWh. See 'help caseformat'
           for details on the cost curve format.

   CDF2MATP may modify some of the data which are "infeasible" for
   running optimal power flow. If so, warning information will be
   printed out on screen.

   Note: Since our code can not handle transformers with variable tap,
   you may not expect to get exactly the same power flow solution
   using converted data. This is the case when we converted ieee300.cdf.

0001 function cdf2matp(cdf_file, matp_file, mpcver)
0002 %CDF2MATP Converts data from IEEE Common Data Format to MATPOWER format.
0003 %   CDF2MATP                         : prompts for input & output file names
0004 %   CDF2MATP(CDF_FILE, MATP_FILE)    : uses specified file names
0005 %   CDF2MATP(CDF_FILE, MATP_FILE, v) : and MATPOWER case version
0006 %                                      v = '1' or '2', default is '2'
0007 %
0008 %   Optional arguments are the names of the input IEEE CDF file and the
0009 %   output MATPOWER case file. If these are not given you will be prompted
0010 %   to enter them.
0011 %
0012 %   The IEEE CDF does not include some data need to run an optimal power
0013 %   flow. This script creates default values for some of this data as
0014 %   follows:
0015 %
0016 %       Bus data:
0017 %           Vmin = 0.94 p.u.
0018 %           Vmax = 1.06 p.u.
0019 %       Gen data:
0020 %           Pmin = 0 MW
0021 %           Pmax = Pg + baseMVA
0022 %       Gen cost data:
0023 %           Quadratic costs with:
0024 %               c2 = 10 / Pg, c1 = 20, c0 = 0, if Pg is non-zero, and
0025 %               c2 = 0.01,    c1 = 40, c0 = 0, if Pg is zero
0026 %           This should yield an OPF solution "close" to the
0027 %           existing solution (assuming it is a solved case)
0028 %           with lambdas near $40/MWh. See 'help caseformat'
0029 %           for details on the cost curve format.
0030 %
0031 %   CDF2MATP may modify some of the data which are "infeasible" for
0032 %   running optimal power flow. If so, warning information will be
0033 %   printed out on screen.
0034 %
0035 %   Note: Since our code can not handle transformers with variable tap,
0036 %   you may not expect to get exactly the same power flow solution
0037 %   using converted data. This is the case when we converted ieee300.cdf.
0038 
0039 %   MATPOWER
0040 %   $Id: cdf2matp.m,v 1.25 2011/12/01 19:06:47 cvs Exp $
0041 %   by Deqiang (David) Gan, PSERC Cornell & Zhejiang University
0042 %   Copyright (c) 1996-2010 by Power System Engineering Research Center (PSERC)
0043 %
0044 %   This file is part of MATPOWER.
0045 %   See http://www.pserc.cornell.edu/matpower/ for more info.
0046 %
0047 %   MATPOWER is free software: you can redistribute it and/or modify
0048 %   it under the terms of the GNU General Public License as published
0049 %   by the Free Software Foundation, either version 3 of the License,
0050 %   or (at your option) any later version.
0051 %
0052 %   MATPOWER is distributed in the hope that it will be useful,
0053 %   but WITHOUT ANY WARRANTY; without even the implied warranty of
0054 %   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
0055 %   GNU General Public License for more details.
0056 %
0057 %   You should have received a copy of the GNU General Public License
0058 %   along with MATPOWER. If not, see <http://www.gnu.org/licenses/>.
0059 %
0060 %   Additional permission under GNU GPL version 3 section 7
0061 %
0062 %   If you modify MATPOWER, or any covered work, to interface with
0063 %   other modules (such as MATLAB code and MEX-files) available in a
0064 %   MATLAB(R) or comparable environment containing parts covered
0065 %   under other licensing terms, the licensors of MATPOWER grant
0066 %   you additional permission to convey the resulting work.
0067 
0068 %% define named indices into bus, gen, branch matrices
0069 [PQ, PV, REF, NONE, BUS_I, BUS_TYPE, PD, QD, GS, BS, BUS_AREA, VM, ...
0070     VA, BASE_KV, ZONE, VMAX, VMIN, LAM_P, LAM_Q, MU_VMAX, MU_VMIN] = idx_bus;
0071 [GEN_BUS, PG, QG, QMAX, QMIN, VG, MBASE, GEN_STATUS, PMAX, PMIN, ...
0072     MU_PMAX, MU_PMIN, MU_QMAX, MU_QMIN, PC1, PC2, QC1MIN, QC1MAX, ...
0073     QC2MIN, QC2MAX, RAMP_AGC, RAMP_10, RAMP_30, RAMP_Q, APF] = idx_gen;
0074 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
0075     TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
0076     ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
0077 [PW_LINEAR, POLYNOMIAL, MODEL, STARTUP, SHUTDOWN, NCOST, COST] = idx_cost;
0078 
0079 rev = '$Revision: 1.25 $';
0080 
0081 % ----- get the original load flow data -----
0082 if nargin < 1
0083     cdf_file =  input('Please enter IEEE Common Data Format file name : ', 's');
0084 end
0085 if nargin < 2
0086     matp_file = input('Please enter MATPOWER file name                : ', 's');
0087 end
0088 if nargin < 3
0089     mpcver = '2';
0090 end
0091 
0092 %% verify valid input filename
0093 [cdf_path cdf_name cdf_ext] = fileparts(cdf_file);
0094 if isempty(cdf_ext)
0095     cdf_ext = '.cdf';
0096     cdf_file = strcat(cdf_name , cdf_ext);
0097 end
0098 
0099 %% open input file
0100 [fid, msg] = fopen(cdf_file, 'r');
0101 if fid < 0
0102     disp(msg);
0103     error('cdf2matp: Can not read the input file:  %s', cdf_file )
0104 end
0105 
0106 %% verify valid output filename
0107 if isempty(matp_file) 
0108     matp_file = strcat(cdf_name, '.m');
0109     rootname = cdf_name;
0110 else 
0111     [matp_path matp_name matp_ext] = fileparts(matp_file);
0112     rootname = matp_name;
0113     matp_file = strcat(matp_name, matp_ext);
0114 end
0115 
0116 % set up some comments
0117 comments = { sprintf('%s', upper(rootname)), ...
0118     '  Please see ''help caseformat'' for details on the case file format.', ...
0119     '  This data was converted from IEEE Common Data Format', ...
0120     sprintf('  (%s) on %s by cdf2matp, rev. %s', cdf_file, date, rev(12:end-2)), ...
0121     '  See end of file for warnings generated during conversion.', ...
0122     '' };
0123 warnings = {};
0124 
0125 % get baseMVA
0126 title_cdf = fgetl(fid);
0127 if isnumeric(str2num(title_cdf(32:37))) && not(isempty(str2num(title_cdf(32:37))))
0128     baseMVA = str2num(title_cdf(32:37));
0129     if length(findstr(title_cdf(2:9), '/')) == 2    %% date in the file
0130         warnings{end+1} = sprintf('***** check the title format in the first line of the cdf file.');
0131     end
0132 else
0133     error('cdf2matp:argChk','Error getting the Base MVA, check the title format in the first line of the file.')
0134 end
0135     
0136 % adding the cdf title
0137 comments{end+1} = title_cdf;
0138 
0139 % find string 'BUS DATA FOLLOWS'
0140 while 1
0141     line = fgetl(fid);
0142     if line(1:16) == 'BUS DATA FOLLOWS', break, end
0143 end
0144 
0145 % ----- get bus data, feed them into matrix bus, gen, gencost
0146 ibus = 0;
0147 igen = 0;
0148 iarea = 0;
0149 
0150 while 1
0151     line = fgetl(fid);
0152     if line(1:4) == '-999', break, end
0153 
0154     % feed bus data
0155     ibus = ibus + 1;
0156     bus(ibus, BUS_I) = str2num(line(1:4));          % bus number
0157     busnames(ibus,:) = line(6:17);                  % bus names
0158     bus(ibus, BUS_TYPE) = str2num(line(25:26));
0159     if bus(ibus, BUS_TYPE) == 0                     % bus type
0160         bus(ibus, BUS_TYPE) = 1;
0161     end
0162     if (bus(ibus, BUS_TYPE) < 2)                    % Pd
0163         bus(ibus, PD) = str2num(line(41:49)) - str2num(line(60:67));
0164     elseif (bus(ibus, BUS_TYPE) >= 2)
0165         bus(ibus, PD) = str2num(line(41:49));
0166     end
0167     bus(ibus, QD) = str2num(line(50:59));           % Qd
0168     bus(ibus, GS) = baseMVA*str2num(line(107:114)); % Gs
0169     bus(ibus, BS) = baseMVA*str2num(line(115:122)); % Bs
0170     bus(ibus, BUS_AREA) = str2num(line(19:20));     % area
0171     bus(ibus, VM) = str2num(line(28:33));           % Vm
0172     bus(ibus, VA) = str2num(line(34:40));           % Va
0173     bus(ibus, BASE_KV) = str2num(line(77:83));      % baseKV
0174     bus(ibus, ZONE) = str2num(line(21:23));         % zone
0175     bus(ibus, VMAX) = 1.06;                         % default voltage upper limit
0176     bus(ibus, VMIN) = 0.94;                         % default voltage lower limit
0177 
0178     % feed gen and gencost
0179     Pg = str2num(line(60:67));
0180     Qg = str2num(line(68:75));
0181     Qmax = str2num(line(91:98));
0182     Qmin = str2num(line(99:106));
0183     if bus(ibus, BUS_TYPE) >= 2
0184         igen = igen + 1;
0185         if bus(ibus, BUS_TYPE) == 3, refgen = igen; end
0186         gen(igen, GEN_BUS) = bus(ibus, BUS_I);      % bus number
0187         gen(igen, PG) = Pg;                         % Pg
0188         if gen(igen, PG) < 0                        % negative Pg is transformed as load
0189             bus(ibus, PD) = bus(ibus, PD) - gen(igen, PG);
0190             warnings{end+1} = sprintf('***** negative Pg at bus %g treated as Pd', bus(ibus, BUS_I));
0191             fprintf('\n %s', warnings{end});
0192             gen(igen, PG) = 0;
0193         end
0194         gen(igen, QG)   = Qg;                       % Qg
0195         gen(igen, QMAX) = Qmax;                     % Qmax
0196         gen(igen, QMIN) = Qmin;                     % Qmin
0197         if Qmax - Qmin < 0.01                       % Qmax is modified
0198             gen(igen, QMAX) = Qmin + 0.1 * baseMVA;
0199             warnings{end+1} = sprintf('***** Qmax = Qmin at generator at bus %4i (Qmax set to Qmin + %g)', bus(ibus, BUS_I), baseMVA/10);
0200             fprintf('\n %s', warnings{end});
0201         end
0202         gen(igen, VG)    = str2num(line(85:90));    % specified voltage
0203         gen(igen, MBASE) = baseMVA;                 % baseMVA
0204         gen(igen, GEN_STATUS) = 1;                  % default status is 'on'
0205         gen(igen, PMAX)  = gen(igen, 2) + baseMVA;  % Pmax
0206         gen(igen, PMIN)  = 0;                       % Pmin = 0 by default
0207 
0208         gencost(igen, MODEL)    = POLYNOMIAL;       % by default, sets the model as polynomial
0209         gencost(igen, STARTUP)  = 0;                % start up cost is zero by default
0210         gencost(igen, SHUTDOWN) = 0;                % shut down cost is zero by default
0211         gencost(igen, NCOST)    = 3;                % number of coefficients in polynomial cost
0212 %       gencost(igen, COST)     = 0.01;             % default c2
0213 %       gencost(igen, COST+1)   = 0.3;              % default c1
0214 %       gencost(igen, COST+2)   = 0.2;              % default c0
0215     end
0216 end
0217 
0218 totload = sum(bus(:, PD));
0219 totgen = sum(gen(:, PG));
0220 if totgen < 1.04 * totload
0221     gen(refgen, PMAX) = gen(refgen, PG) + 1.1 * totload - totgen;   % Pg at slack bus is modified
0222     warnings{end+1} = sprintf('***** Insufficient generation, setting Pmax at slack bus (bus %d) to %g', gen(refgen, [GEN_BUS, PMAX]));
0223     fprintf('\n %s', warnings{end});
0224 end
0225 
0226 % ----- set up the cost coefficients of generators
0227 ng = size(gen, 1);
0228 % gencost(:, COST)    = zeros(ng, 1);
0229 % gencost(:, COST+1)  = 100*ones(ng, 1) ./ (gen(:, PG) + 10*ones(ng, 1));
0230 % gencost(:, COST+2)  = 100*ones(ng, 1) ./ (gen(:, PG) + 10*ones(ng, 1));
0231 zg  = find(gen(:, PG) == 0);                %% for Pg = 0
0232 gencost(zg, COST)  = 0.01 * ones(size(zg));
0233 gencost(zg, COST+1) = 40 * ones(size(zg));
0234 nzg = find(gen(:, PG) ~= 0);                %% Pg non-zero
0235 gencost(nzg, COST) = 10 * ones(size(nzg)) ./ gen(nzg, PG);
0236 gencost(nzg, COST+1) = 20 * ones(size(nzg));
0237 gencost(:, COST+2) = zeros(ng, 1);
0238 
0239 % find string 'BRANCH DATA FOLLOWS'
0240 while 1
0241     line = fgetl(fid);
0242     if line(1:19) == 'BRANCH DATA FOLLOWS', break, end
0243 end
0244 
0245 % ----- get branch data, feed them into matrix branch
0246 k = 0;
0247 while 1
0248     line = fgetl(fid);
0249     if line(1:4) == '-999', break, end
0250 
0251     k = k + 1;
0252     branch(k, F_BUS)  = str2num(line(1:4));     % fbus (also the tap bus)
0253     branch(k, T_BUS)  = str2num(line(6:9));     % tbus
0254     branch(k, BR_R)   = str2num(line(20:29));   % R
0255     branch(k, BR_X)   = str2num(line(30:40));   % X
0256     branch(k, BR_B)   = str2num(line(41:50));   % B
0257     branch(k, RATE_A) = str2num(line(51:55));   % RATE A
0258     if branch(k, RATE_A) < 0.000001
0259         branch(k, RATE_A) = 99 * baseMVA;       % RATE A is modified
0260         warnings{end+1} = sprintf('***** MVA limit of branch %d - %d not given, set to %g', branch(k, [1, 2, 6]));
0261         fprintf('\n %s', warnings{end});
0262     end
0263     branch(k, RATE_B) = str2num(line(57:61));   % RATE B
0264     branch(k, RATE_C) = str2num(line(63:67));   % RATE C
0265     branch(k, TAP)    = str2num(line(77:82));   % transformer turns ratio
0266     branch(k, SHIFT)  = 0;                      % phase shifter can not be modelled
0267     branch(k, BR_STATUS) = 1;                   % by default, branch is on
0268 end
0269 fprintf('\n');
0270 fclose(fid);
0271 
0272 % put in struct
0273 mpc.baseMVA = baseMVA;
0274 mpc.bus     = bus;
0275 mpc.branch  = branch;
0276 mpc.gen     = gen;
0277 mpc.gencost = gencost;
0278 mpc = loadcase(mpc);    %% convert to internal (e.g. v. '2') case format
0279 
0280 % ----- write data in MATPOWER format -----
0281 matp_file = savecase(matp_file, comments, mpc, mpcver);
0282 
0283 % print the name buses
0284 busnames = strtrim(busnames);       %% remove leading and trailing blanks
0285 if length(matp_file) < 4 || ~strcmp(matp_file(end-3:end), '.mat')
0286     fid = fopen(matp_file, 'at');   %% open for append
0287     fprintf(fid, '\n%% bus names\nmpc.busnames = [\n');
0288     for i = 1:size(busnames, 1)
0289         fprintf(fid, '\t''%s''\n', busnames(i, :));
0290     end
0291     fprintf(fid, '];\n');
0292     fclose(fid);
0293 end
0294 
0295 % print conversion warnings as comments
0296 if ~isempty(warnings)
0297     if length(matp_file) < 4 || ~strcmp(matp_file(end-3:end), '.mat')
0298         fid = fopen(matp_file, 'at');       % open for append
0299         fprintf(fid, '\n%% Warnings from cdf2matp conversion:\n%%\n');
0300         for i=1:length(warnings)
0301             fprintf(fid, '%% %s\n', warnings{i});
0302         end
0303         fclose(fid);
0304     end
0305 end
0306 
0307 % to do:
0308 %   PTI format can not be handled