fm_cymflow.m

电力系统的psat

function check = fm_cymflow(filename, pathname, flag)
% FM_CYMFLOW convert CYME power flow data files
%        into PSAT data format
%
%CHECK = FM_CYMFLOW(FILENAME,PATHNAME,FLAG)
%       FILENAME name of the file to be converted
%       PATHNAME path of the file to be converted
%       FLAG = 1 append original file as a comment
%       FLAG = 0 do not append original file
%
%       CHECK = 1 conversion completed
%       CHECK = 0 problem encountered (no data file created)
%
%Author:    Federico Milano
%Date:      15-May-2003
%Version:   1.0.0
%
%E-mail:    [e-mail]
%Web-site:  http://thunderbox.uwaterloo.ca/~fmilano
%
% Copyright (C) 2002-2005 Federico Milano
%
% This toolbox is free software; you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation; either version 2.0 of the License, or
% (at your option) any later version.
%
% This toolbox is distributed in the hope that it will be useful, but
% WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANDABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
% General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with this toolbox; if not, write to the Free Software
% Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307,
% USA.

global Settings

check = 1;
if ~strcmp(pathname(end),filesep)
  pathname = [pathname,filesep];
end

fm_disp
fm_disp('Conversion from CYME power flow data format...');
fm_disp(['Source data file "',pathname,filename,'"'])

fid = fopen([pathname,filename]);
b80 = blanks(80);
if fid == -1,
  fm_disp(['Can''t open file ',pathname,filename],2),
  check = 0;
  return,
end

% header, power base and frequency
% ----------------------------------------------------------------------

heading = fgetl(fid);
heading = [heading,b80];

if ~strcmp(heading(1:2),'91') & ~strcmp(heading(1:2),'81')
  fm_disp('Header Row is in a wrong format. Conversion Process Interrupted',2)
  check = 0;
  return
end

ftype = 1;
if strcmp(heading(1:2),'81')
  fm_disp('CYME data file is in free format',2)
  ftype = 0;
end

if ftype
  casetitle = heading(3:42);
  base = valrow(heading(45:50),100);
  vmax = valrow(heading(51:56),1.1);
  vmin = valrow(heading(57:62),0.9);
  init = valrow(heading(68));
  iter = valrow(heading(69:71),60);
  dtol = valrow(heading(72:76),1e-4);
  intc = valrow(heading(79));
  ovld = valrow(heading(80));
else
  [casetitle,base,vamx,vmin,init,iter,dtol,intc,ovld] = ...
      strread(heading(4:end),'%s%f%f%f%d%d%d%d%d','delimiter',',');
end

freq = 60;

% initialize data matrices
bus = [];
namebus = [];
swbus = [];
pvbus = [];
pqbus = [];
shunt = [];
linedata = [];

% read first data line
crow = fgetl(fid);
if crow == -1
  fm_disp('CYME data file is likely empty. Conversion Process Interrupted',2)
  fclose(fid);
end
crow = [crow,b80];
ccod = crow(1:2);
eof = 0;

while 1

  if eof, break, end

  arow = fgetl(fid);
  if arow == -1
    fclose(fid);
    eof = 1;
    arow = '02';
  end
  arow = [arow,b80];
  acod = arow(1:2);

  switch ccod

   case '92' % comment

    rtype = 'comment';
    if ftype
      cmnt = deblank(crow(8:79));
    else
      cmnt = deblank(crow(6:end));
    end

   case '01' % area data

    rtype = 'area';
    if ftype
      area = deblank(crow(3:42));
      ano  = valrow(crow(45:48));
      ana  = crow(49:54);
      azo  = valrow(crow(55:56));
      sint = valrow(crow(57:62));
      tol  = valrow(crow(63:68));
      pmax = valrow(crow(69:74));
      pmin = valrow(crow(75:80));
    else
      ana = '';
      azo = 0;
      [area,ano,sint,tol,pmax,pmin] = ...
          strread(crow(4:end),'%s%d%f%f%f%f','delimiter',',');
    end

   case '02' % bus data

    rtype = 'bus';
    bvmax = vmax;
    bvmin = vmin;
    if ftype
      bno  = valrow(crow(3:6));
      bna  = crow(7:12);
      bzo  = valrow(crow(13:14));
      kvb  = valrow(crow(15:20),230);
      pl   = valrow(crow(21:26));
      ql   = valrow(crow(27:32));
      shur = valrow(crow(33:38));
      shux = valrow(crow(39:44));
      kvi  = valrow(crow(45:50),kvb);
      pgen = valrow(crow(51:56));
      qmax = valrow(crow(57:62));
      qmin = valrow(crow(63:69));
      type = valrow(crow(70),2);
      angl = valrow(crow(71:77));
    else
      [bno,bna,bzo,type,kvb,kvi,d1,d2,angl,pl,ql,shur,shux] = ...
          strread(crow(4:end),'%d%s%d%d%f%f%d%d%f%f%f%f%f','delimiter',',');
    end

   case '03' % ac line data

    rtype = 'acline';
    if ftype
      sno  = valrow(crow(3:6));
      sna  = crow(7:12);
      szo  = valrow(crow(13:14));
      rno  = valrow(crow(15:18));
      rna  = crow(19:24);
      rzo  = valrow(crow(25:26));
      cnum = valrow(crow(29),1);
      rl   = valrow(crow(30:35));
      xl   = valrow(crow(36:41));
      bl   = valrow(crow(42:47));
      ll   = valrow(crow(48:53));
      kvf  = valrow(crow(54:59));
      kvt  = valrow(crow(60:65),1);
      angl = valrow(crow(66:71));
    else
      sna = '';
      szo = 0;
      rna = '';
      rzo = 0;
      [sno,rno,cnum,rl,xl,bl,ll,kvf,kvt,angl] = ...
          strread(crow(4:end),'%d%d%d%f%f%f%f%f%f%f','delimiter',',');
    end

   case '04' % dc line data

    rtype = 'dcline';
    if ftype
      rno  = valrow(crow(3:6));
      rna  = crow(7:12);
      rzo  = valrow(crow(13:14));
      rnb  = valrow(crow(15:18));
      kvbr = valrow(crow(19:24));
      alfa = valrow(crow(37:42));
      rdc  = valrow(crow(43:48));
    else
      rna = '';
      rzo = 0;
      [rno,rnb,kvbr,alfa,rdc] = ...
          strread(crow(4:end),'%d%s%d%d%f%f%f','delimiter',',');
    end

  end

  if strcmp(acod,'00')
    switch ccod

     case '01' % area data
      if ftype
        zones = str2num(arow(3:80));
      else
        zones = sscanf(arow(4:end),'%d,');
      end

     case '02' % bus data
      if ftype
        cno  = valrow(arow(3:6));
        cna  = arow(7:12);
        czo  = valrow(arow(13:14));
        kvcb = valrow(arow(15:20));
        ib   = valrow(arow(21:24));
        r    = valrow(arow(25:32));
        x    = valrow(arow(33:40));
        bvmax = valrow(arow(41:46),vmax);
        bvmin = valrow(arow(47:52),vmin);
      else
        cna = '';
        czo = 0;
        [pgen,qmax,qmin,cno,kvcb,ib,r,x] = ...
            strread(arow(4:end),'%f%f%f%d%f%d%f%f','delimiter',',');
      end

     case '03' % ac line data
      if angl
        % phase shifting transformer
        rtype = 'phshifter';
        if ftype
          amax = valrow(arow(15:20));
          amin = valrow(arow(21:26));
          step = valrow(arow(27:29));
          pmax = valrow(arow(30:36));
          pmin = valrow(arow(37:43));
        else
          [amax,amin,taps,pmax,pmin] = ...
              strread(arow(5:end),'%d%f%f%d%f%f','delimiter',',');
        end

      elseif isempty(num2str(arow(3:6)))
        % reactive power regulating transformer
        rtype = 'qrtc';
        if ftype
          kvhf = valrow(arow(15:20));
          kvlf = valrow(arow(21:26));
          taps = valrow(arow(27:29));
          qmax = valrow(arow(30:36));
          qmin = valrow(arow(37:43));
        else
          [cno,kvhf,kvlf,taps,qmax,qmin] = ...
              strread(arow(4:end),'%d%f%f%d%f%f','delimiter',',');
        end

      else
        % voltage regulating transformer
        rtype = 'vrtc';
        if ftype
          cno  = valrow(arow(3:6));
          cna  = arow(7:12);
          czo  = valrow(arow(13:14));
          kvhf = valrow(arow(15:20));
          kvlf = valrow(arow(21:26));
          taps = valrow(arow(27:29));
          kvht = valrow(arow(30:36));
          kvlt = valrow(arow(37:43));
        else
          cna = '';
          czo = 0;
          [cno,kvhf,kvlf,taps,kvht,kvlt] = ...
              strread(arow(4:end),'%d%f%f%d%f%f','delimiter',',');
        end

      end

     case '04' % dc line data
      if ftype
        ino  = valrow(arow(3:6));
        ina  = arow(7:12);
        izo  = valrow(arow(13:14));
        inb  = valrow(arow(15:18));
        kvbi = valrow(arow(19:24));
        gama = valrow(arow(37:42));
        pdes = valrow(arow(43:48));
        kvdc = valrow(arow(49:54));
      else
        ina = '';
        izo = 0;
        [ino,inb,kvbr,gama,pdes,kvdc] = ...
            strread(arow(4:end),'%d%d%f%f%f%f','delimiter',',');
      end

    end

    crow = fgetl(fid);
    if crow == -1
      eof = 1;
    else
      crow = [crow,b80];
      ccod = crow(1:2);
    end

  else

    crow = arow;
    ccod = acod;

  end

  switch rtype

   case 'bus'
    v0 = kvi/kvb;
    if init, a0 = 0; else, a0 = angl*pi/180; end
    if init, vb = 1; else, vb = v0; end
    bus = [bus; bno, kvb, vb, a0];
    if isempty(deblank(bna))
      bna = ['Bus',num2str(bno),'   '];
      bna = bna(1:6);
    end
    namebus = [namebus;bna];
    switch type
     case 1
      swbus = [swbus; bno,base,kvb,v0,a0,qmax/base,qmin/base, ...
              bvmax,bvmin,pgen/base,1];
     case 2
      pvbus = [pvbus; bno,base,kvb,pgen/base,v0,qmax/base, ...
               qmin/base,bvmax,bvmin,0];
    end
    if pl | ql
      pqbus = [pqbus; bno,base,kvb,pl/base,ql/base,bvmax,bvmin,0];
    end
    if shur | shux
      shunt = [shunt; bno,base,kvb,freq,shur,shux];
    end

   case 'acline'
    linedata = [linedata; sno, rno, base, 0, freq, 0, kvf/kvt, ...
                rl, xl, bl, 1, angl, 0, 0, 0];
    if ovld
      linedata(end,15) = ll/base;
    else
      linedata(end,13) = ll;
    end
   case 'phshifter'
    fm_disp(['Phase shifter on line #', ...
             num2str(length(linedata(:,1)+1)),'will be neglected.'])
    linedata = [linedata; sno, rno, base, 0, freq, 0, kvf/kvt, ...
                rl, xl, bl, 1, angl, 0, 0, 0];
    if ovld
      linedata(end,15) = ll/base;
    else
      linedata(end,13) = ll;
    end
   case 'vrtc'

   case 'qrtc'

   case 'dcline'

  end

end

for i = 1:length(linedata(:,1))
  idxf = find(bus(:,1)==linedata(i,1));
  idxt = find(bus(:,1)==linedata(i,2));
  vbf = bus(idxf,2);
  vbt = bus(idxt,2);
  linedata(i,4) = vbf;
  if vbf == vbt
    linedata(i,7) = 0;
  end
  if linedata(i,13)
    linedata(i,13) = linedata(i,13)*sqrt(3)*vbf/1000/base;
  end
end

% load original file into cell string
if flag == 1,
  orig_file = textread([pathname,filename],'%s','delimiter','\n', ...
                       'whitespace','');
end

% definition of file name for PSAT data file
extension = findstr(filename,'.');
newfile = ['d_',fm_filenum(filename(1:extension(end)-1),1),'.m'];
% open *.m file for writing data
fid = fopen([pathname,newfile], 'wt');

% Bus data: Bus.con
% ----------------------------------------------------------------------

count = fprintf(fid,['%% ',datestr(now,2), ...
                    ' File data originated from CYME data file\n']);
count = fprintf(fid,['%% \n']);
count = fprintf(fid, ['%%  ', casetitle, '\n\n']);
count = fprintf(fid, 'Bus.con = [ ...\n');
for i = 1:length(bus(:,1))-1
  count = fprintf(fid,'%4d %8.4g %8.4g %8.4g; ', bus(i,:));
  if rem(i,5) == 0; count = fprintf(fid,'\n'); end
end
count = fprintf(fid, '%4d %8.4g %8.4g %8.4g];\n\n\n', bus(end,:));

% Swing Generator data: SW.con
% ----------------------------------------------------------------------

if ~isempty(swbus)
  count = fprintf(fid, 'SW.con = [ ...\n');
  format = '%4d %8.4g %8.4g %8.4g %8.4g %8.4g %8.4g %8.4g %8.4g %8.4g %2u';
  count = fprintf(fid,[format,';\n'],swbus(1:end-1,:)');
  count = fprintf(fid,[format,'];\n\n\n'],swbus(end,:));
end

% PV Generator data: PV.con
% ----------------------------------------------------------------------

if ~isempty(pvbus)
  count = fprintf(fid, 'PV.con = [ ...\n');
  format = '%4d %8.4g %8.4g %8.4g %8.4g %8.4g %8.4g %8.4g %8.4g %2u';
  count = fprintf(fid,[format,';\n'],pvbus(1:end-1,:)');
  count = fprintf(fid,[format,'];\n\n\n'],pvbus(end,:));
end

% Constant Power Load data: PQ.con
% ----------------------------------------------------------------------

if ~isempty(pqbus)
  count = fprintf(fid, 'PQ.con = [ ...\n');
  format = '%4d %8.4g %8.4g %8.4g %8.4g %8.4g %8.4g %2u';
  count = fprintf(fid,[format,';\n'],pqbus(1:end-1,:)');
  count = fprintf(fid,[format,'];\n\n\n'],pqbus(end,:));
end

% Shunt Impedance data: Shunt.con
% ----------------------------------------------------------------------

if ~isempty(shunt)
  count = fprintf(fid, 'Shunt.con = [ ...\n');
  count = fprintf(fid,'%4d %8.4g %8.4g %8.4g %8.4g %8.4g;\n',shunt(1:end-1,:)');
  count = fprintf(fid,'%4d %8.4g %8.4g %8.4g %8.4g %8.4g];\n\n\n',shunt(end,:));
end

% Branch data: Line.con
% ----------------------------------------------------------------------

count = fprintf(fid, 'Line.con = [ ...\n');
format = ['%4d %4d %8.4g %8.4g %8.4g %8.4g %8.4g %8.4g %8.4g %8.4g ' ...
          '%8.4g %8.4g %8.4g %8.4g %8.4g'];
count = fprintf(fid,[format,';\n'],linedata(1:end-1,:)');
count = fprintf(fid,[format,'];\n\n\n'],linedata(end,:));

% Bus Names: Varname.bus
% ----------------------------------------------------------------------

count = fprintf(fid, 'Varname.bus = {...\n      ');
for i = 1:length(namebus(:,1))-1
  count = fprintf(fid, ['''', deblank(namebus(i,:)),'''; ']);
  if rem(i,5) == 0; count = fprintf(fid,'\n      '); end
end
count = fprintf(fid, ['''', deblank(namebus(end,:)),'''};\n\n']);

% append original data file
if flag == 1
  count = fprintf(fid, '\n\n%% ORIGINAL DATA FILE: \n\n');
  for i = 1:length(orig_file); count = fprintf(fid,'%% %s \n\n',orig_file{i,1}); end
end

% end of operations
fm_disp(['Conversion into data file "',pathname,newfile,'" completed.'])
if Settings.beep, beep, end
count = fclose(fid);

% ------------------------------------
function output = valrow(varargin)
% ------------------------------------

output = str2num(varargin{1});
if isempty(output)
  if nargin == 1
    output = 0;
  elseif nargin == 2
    output = varargin{2};
  end
end