fm_plot.m

电力系统的psat

function fm_plot(flag)
% FM_PLOT plot results of Continuation Power Flow,
%         Optimal Power Flow and Tiam Domain
%         Simulations.
%
% FM_PLOT(FLAG)
%    FLAG  0 -> create variable list
%          1 -> plot selected variables
%          2 -> save graph
%          3 -> set layout
%
%Author:    Federico Milano
%Date:      11-Nov-2002
%Update:    25-Feb-2003
%Update:    26-Jan-2005
%Version:   1.0.2
%
%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 DAE Bus Line Syn Settings Fig Snapshot Hdl
global Varout Varname Path OPF Theme Exc Oxl GAMS

if isempty(Varout.t) & isempty(OPF.varout)
  fm_disp('Plotting Utilities: No data available for plotting.')
  return
end

%allbaby = get(Fig.plot,'Children');
%hdlfig = allbaby(end-1);
%hdlfig = findobj(Fig.plot,'Tag','Axes1');
%hdlfig2 = allbaby(end);
%hdlfig2 = findobj(Fig.plot,'Tag','Axes2');

hdlfig = Hdl.axesplot;
hdlfig2 = Hdl.axeslogo;

Hdl_grid = findobj(Fig.plot,'Tag','Checkbox1');
Hdl_legend = findobj(Fig.plot,'Tag','Checkbox2');
Hdl_listvar = findobj(Fig.plot,'Tag','Listbox1');
Hdl_listplot = findobj(Fig.plot,'Tag','Listbox2');
Hdl_tipoplot = findobj(Fig.plot,'Tag','PopupMenu1');
Hdl_angref = findobj(Fig.plot,'Tag','PopupMenu2');
Hdl_snap = findobj(Fig.plot,'Tag','Radiobutton1');
hdl_zoom1 = findobj(Fig.plot,'Tag','Pushbutton12');
hdl_zoom2 = findobj(Fig.plot,'Tag','Pushbutton11');
hdl_zoom3 = findobj(Fig.plot,'Tag','Pushbutton4');
hdl_x = findobj(Fig.plot,'Tag','Pushbutton9');
hdl_y = findobj(Fig.plot,'Tag','Pushbutton5');
hdl_xy = findobj(Fig.plot,'Tag','Pushbutton10');

if flag == 99,  % output data as plain text file 
  flag = 1;
  out_matlab = 0;
  out_mtv = 0;
  out_text = 1;
elseif flag == 98, % output data as MTV plot file
  flag = 1;
  out_matlab = 0;
  out_mtv = 1;
  out_text = 0;
elseif flag == 97, % output data as Matlab script file
  flag = 1;
  out_matlab = 1;
  out_mtv = 0;
  out_text = 0;
else
  out_matlab = 0;
  out_mtv = 0;
  out_text = 0;
end

switch flag
 case 0
  if ~isempty(OPF.varout)
    set(gcbo,'String',enum(OPF.uname));
    Varname.changes = sparse(length(OPF.uname),1);
    set(gcf,'UserData',OPF.fname);
  elseif ~Settings.vs
    set(gcbo,'String',enum([Varname.ux; Varname.uV; Varname.uang; ...
                        Varname.uP; Varname.uQ; Varname.uPflow; ...
                        Varname.uQflow; Varname.uPm; Varname.uVf; ...
                        Varname.uVref; Varname.uIf]));
    Varname.changes = sparse(DAE.n+4*Bus.n+4*Line.n+ ...
                             2*Syn.n+Exc.n+Oxl.n,1);
    set(gcf,'UserData',[Varname.fx; Varname.fV; Varname.fang; ...
                        Varname.fP; Varname.fQ; Varname.fPflow; ...
                        Varname.fQflow; Varname.fPm; Varname.fVf; ...
                        Varname.fVref; Varname.fIf]);
  elseif Settings.vs
    pflf =  cell(Bus.n*Bus.n,1);
    pflfd = cell(Bus.n*Bus.n,1);
    count = 1;
    for i = 1:Bus.n
      for j = 1:Bus.n
        pflf{count} = Varname.upflf{j,i};
        pflfd{count} = Varname.upflfd{j,i};
        count = count + 1;
      end
    end
    set(gcbo,'String',enum([Varname.ux; Varname.uV; Varname.uang; ...
                        Varname.uP; Varname.uQ; Varname.uPflow; ...
                        Varname.uQflow; Varname.uPm; Varname.uVf; ...
                        Varname.uVref; Varname.uIf; Varname.udet; ...
                        Varname.uautostate; Varname.uautojlfr; ...
                        Varname.uautojlfdr; pflf; pflfd]));
    Varname.changes = sparse(2*DAE.n+4*Line.n+2*Syn.n+3+ ...
                             (6+2*Bus.n)*Bus.n+Exc.n+Oxl.n,1);
    set(gcf,'UserData',[Varname.fx; Varname.fV; Varname.fang; ...
                        Varname.fP; Varname.fQ; Varname.fPflow; ...
                        Varname.fQflow; Varname.fPm; Varname.fVf; ...
                        Varname.fVref; Varname.fIf; Varname.fdet; ...
                        Varname.fautostate; Varname.fautojlfr; ...
                        Varname.fautojlfdr; pflf; pflfd]);
  end
  set(gcf, ...
      'DefaultAxesColorOrder',Settings.color, ...
      'DefaultAxesLineStyle','-');

 case 1
  
  nB = Bus.n;
  nL = Line.n;
  nS = Syn.n;
  nE = Exc.n;
  nO = Oxl.n;
  
  Value = get(Hdl_listvar,'Value');
  if isempty(Value),
    return,
  end
  AxesFont = get(hdlfig,'FontName');
  AxesColor = get(hdlfig,'Color');
  AxesWeight = get(hdlfig,'FontWeight');
  AxesAngle = get(hdlfig,'FontAngle');
  AxesSize = get(hdlfig,'FontSize');
  AxesUnits = get(hdlfig,'FontUnits');
  
  if ~isempty(Varout.t)
    ang_ref = get(Hdl_angref,'Value');
    switch ang_ref
     case 1
      angolo = zeros(length(Varout.t),1);
     otherwise
      if ang_ref <= length(Syn.delta)+1
        val_ang = Syn.delta(ang_ref-1);
        angolo = Varout.x(:,val_ang);
      elseif ang_ref <= length(Syn.delta)+nB+1
        val_ang = DAE.n+ang_ref-1-length(Syn.delta)+nB;
        angolo = Varout.ang(:,val_ang-DAE.n-nB);
      end
    end
  end
  
  plot_snap = get(Hdl_snap,'Value');
  if plot_snap
    snap_idx = [];
    for i = 1:length(Snapshot);
      a = find(Varout.t == Snapshot(i).time);
      if isempty(a)
        fm_disp(['Plotting utilities: Snapshots do not match ' ...
                 'current simulation data'],2)
        Hdl_rad1 = findobj(gcf,'Tag','Radiobutton1');
        set(Hdl_rad1,'Value',0);
        plot_snap = 0;
        return
      else
        snap_idx(i) = a;
      end
    end
  else
    snap_idx = 1:length(Varout.t);
  end
  
  Stringa = get(gcf,'UserData');
  legenda = Stringa(Value(Varname.po));
  leg_value = get(Hdl_legend,'Value');
  numero = length(Value);
  punti = length(Varout.t);
  if punti == 0, punti = length(OPF.omega); end
  Y = zeros(punti,numero);
  nD = DAE.n;
  if isempty(Varname.ux), nD = 0; end
  
  if ~isempty(Varout.t)
    for i = 1:numero
      val = Value(i);
      if val <= nD
        didx = Syn.delta;
        if isempty(didx),
          didx = 0;
        end
        if isempty(find(didx == val)) | Varname.changes(val)
          Y(:,i) = Varout.x(:,val);
        else
          Y(:,i) = Varout.x(:,val) - angolo;
          j = find(Varname.po == i);
          if ang_ref > 1,
            legenda{j} = [legenda{j},' - ',Stringa{val_ang}];
          end
        end
      elseif val <= nD+nB
        Y(:,i) = Varout.V(:,val-nD);
      elseif val <= nD + 2*nB
        Y(:,i) =  Varout.ang(:,val-nD-nB) - angolo;
        j = find(Varname.po == i);
        if ang_ref > 1 & ~Varname.changes(val),
          legenda{j} = [legenda{j},' - ',Stringa{val_ang}];
        end
      elseif val <= nD + 3*nB
        Y(:,i) = Varout.p(:,val-nD-2*nB);
      elseif val <= nD + 4*nB
        Y(:,i) = Varout.q(:,val-nD-3*nB);
      elseif val <= nD + 4*nB + nL
        Y(:,i) = Varout.psflow(:,val-nD-4*nB);
      elseif val <= nD + 4*nB + 2*nL
        Y(:,i) = Varout.prflow(:,val-nD-4*nB-nL);
      elseif val <= nD+4*nB+3*nL
        Y(:,i) = Varout.qsflow(:,val-nD-4*nB-2*nL);
      elseif val <= nD+4*nB+4*nL
        Y(:,i) = Varout.qrflow(:,val-nD-4*nB-3*nL);
      elseif val <= nD+4*nB+4*nL+nS
        Y(:,i) = Varout.Pm(:,val-nD-4*nB-4*nL);
      elseif val <= nD+4*nB+4*nL+2*nS
        Y(:,i) = Varout.Vf(:,val-nD-4*nB-4*nL-nS);
      elseif val <= nD+4*nB+4*nL+2*nS+nE
        Y(:,i) = Varout.Vref(:,val-nD-4*nB-4*nL-2*nS);
      elseif val <= nD+4*nB+4*nL+2*nS+nE+nO
        Y(:,i) = Varout.If(:,val-nD-4*nB-4*nL-2*nS-nE);
      elseif val == nD+4*nB+4*nL+2*nS+nE+nO+1
        Y(:,i) =  Varout.detJlf;
      elseif val == nD+4*nB+4*nL+2*nS+nE+nO+2
        Y(:,i) = Varout.detJlfv;
      elseif val == nD+4*nB+4*nL+2*nS+nE+nO+3
        Y(:,i) = Varout.detJlfd;
      elseif val <= 2*nD+4*nB+4*nL+2*nS+nE+nO+3
        Y(:,i) = real(Varout.AutoState(:,val-(nD+4*nB+4*nL+2*nS+nE+nO+3)));
      elseif val <= 2*nD+5*nB+4*nL+2*nS+nE+nO+3
        Y(:,i) = Varout.AutoJlfr(:,val-(2*nD+4*nB+4*nL+2*nS+nE+nO+3));
      elseif val <= 2*nD+6*nB+4*nL+2*nS+nE+nO+3
        Y(:,i) = Varout.AutoJlfdr(:,val-(2*nD+5*nB+4*nL+2*nS+nE+nO+3));
      elseif val <= 2*nD+4*nL+2*nS+nE+nO+3+(6+nB)*nB
        a = rem(val - (2*nD+6*nB+4*nL+2*nS+nE+nO+4),nB);
        b = round((val - (2*nD+6*nB+4*nL+2*nS+nE+nO+4)-a)/nB);
        Y(:,i) = Varout.pf_lf(:,a+1,b+1);
      elseif val <= 2*nD+4*nL+2*nS+nE+nO+3+(6+2*nB)*nB
        a = rem(val-(2*nD+4*nL+2*nS+nE+nO+4)-(6+nB)*nB,nB);
        b = round((val - (2*nD+4*nL+2*nS+nE+nO+4)-(6+nB)*nB-a)/nB);
        Y(:,i) = Varout.pf_lfd(:,a+1,b+1);
      end
    end
  else
    Y = OPF.varout(:,Value);
  end
  
  % X axis variable
  hdlab = findobj(Fig.plot,'Tag','PopupMenu3');
  AbValue = get(hdlab,'Value');
  if ~isempty(Varout.t)
    if AbValue == 1
      X = Varout.t;
    else
      val = AbValue-1;
      if val <= nD
        didx = Syn.delta;
        if isempty(didx), didx = 0; end
        if isempty(find(didx == val)) | Varname.changes(val)
          X = Varout.x(:,val);
        else
          X = Varout.x(:,val)-angolo;
        end
      elseif val <= nD+nB
        X = Varout.V(:,val-nD);
      elseif val <= nD + 2*nB
        X =  Varout.ang(:,val-nD-nB)-angolo;
      elseif val <= nD + 3*nB
        X = Varout.p(:,val-nD-2*nB);
      elseif val <= nD + 4*nB
        X = Varout.q(:,val-nD-3*nB);
      elseif val <= nD + 4*nB + nL
        X = Varout.psflow(:,val-nD-4*nB);
      elseif val <= nD + 4*nB + 2*nL
        X = Varout.prflow(:, val-nD-4*nB-nL);
      elseif val <= nD+4*nB+3*nL
        X = Varout.qsflow(:,val-nD-4*nB-2*nL);
      elseif val <= nD+4*nB+4*nL
        X = Varout.qrflow(:,val-nD-4*nB-3*nL);
      elseif val <= nD+4*nB+4*nL+nS
        X = Varout.Pm(:,val-nD-4*nB-4*nL);
      elseif val <= nD+4*nB+4*nL+2*nS
        X = Varout.Vf(:,val-nD-4*nB-4*nL-nS);
      elseif val <= nD+4*nB+4*nL+2*nS+nE
        X = Varout.Vref(:,val-nD-4*nB-4*nL-2*nS);
      elseif val <= nD+4*nB+4*nL+2*nS+nE+nO
        X = Varout.If(:,val-nD-4*nB-4*nL-2*nS-nE);
      elseif val == nD+4*nB+4*nL+2*nS+nE+nO+1
        X =  Varout.detJlf;
      elseif val == nD+4*nB+4*nL+2*nS+nE+nO+2
        X = Varout.detJlfv;
      elseif val == nD+4*nB+4*nL+2*nS+nE+nO+3
        X = Varout.detJlfd;
      elseif val <= 2*nD+4*nB+4*nL+2*nS+nE+nO+3
        X = real(Varout.AutoState(:,val-(nD+4*nB+4*nL+2*nS+nE+nO+3)));
      elseif val <= 2*nD+5*nB+4*nL+2*nS+nE+nO+3
        X = Varout.AutoJlfr(:,val-(2*nD+4*nB+4*nL+2*nS+nE+nO+3));
      elseif val <= 2*nD+6*nB+4*nL+2*nS+nE+nO+3
        X =  Varout.AutoJlfdr(:,val-(2*nD+5*nB+4*nL+2*nS+nE+nO+3));
      elseif val <= 2*nD+4*nL+2*nS+nE+nO+3+(6+nB)*nB
        a = rem(val-(2*nD+6*nB+4*nL+2*nS+nE+nO+4),nB);
        b = round((val-(2*nD+6*nB+4*nL+2*nS+nE+nO+4)-a)/nB);
        X = Varout.pf_lf(:,a+1,b+1);
      elseif val <= 2*nD+4*nL+2*nS+nE+nO+3+(6+2*nB)*nB
        a = rem(val-(2*nD+4*nL+2*nS+nE+nO+4)-(6+nB)*nB,nB);
        b = round((val-(2*nD+4*nL+2*nS+nE+nO+4)-(6+nB)*nB-a)/nB);
        X = Varout.pf_lfd(:,a+1,b+1);
      end
    end
  else
    if AbValue == 1
      if ~isempty(GAMS.hours)
        X = GAMS.hours;
      else
        X = OPF.wp;
      end
    else
      X = OPF.varout(:,AbValue-1);
    end
  end

  if isempty(Y), return, end
  tipoplot = get(Hdl_tipoplot,'Value');
  
  if out_text
    plainfile = fm_filenum('txt');
    fid = fopen([Path.data,plainfile,'.txt'],'wt');
    if fid == -1
      fm_disp('Cannot open file. Data not saved.')
      return
    end
    fprintf(fid,'C Legend:\n');
    fprintf(fid,'C %s,  ',Settings.xlabel);
    for i = 1:size(Y,2)
      fprintf(fid,'%s,  ',legenda{i});
    end
    fprintf(fid,'\nC Data:\n');
    fprintf(fid,[repmat('%8.5f  ',1,1+size(Y,2)),'\n'],[X,Y]');
    fclose(fid); 
    fm_disp(['Data exported to plain text file "',plainfile,'.txt"'])
  end
  
  if out_mtv
    plainfile = fm_filenum('mtv');
    fid = fopen([Path.data,plainfile,'.mtv'],'wt');
    if fid == -1
      fm_disp('Cannot open file. Data not saved.')
      return
    end
    %fprintf(fid,'$ DATA=CURVE2D\n');
    fprintf(fid,'%% xlabel = "%s"\n',Settings.xlabel);
    if min(X) < max(X)
      fprintf(fid,'%% xmin = %8.5f\n',min(X));
      fprintf(fid,'%% xmax = %8.5f\n',max(X));