MATLAB 图形文件



我对轮廓图有问题,当我运行语法时,我收到如下错误:

Error using contourf (line 66)
Z must be size 2x2 or greater.
Error in example (line 660)
    contourf (x,y,c); colorbar;

这是我的语法:

%---------------------stack1--------------------------------------------
v = 0.5;                    % velocity
lambda = 0;                 % decay rate
Q = 1;                      % emission rate(s)
u = 3;                      % wind speed
P = 1;                      % Ambient Pressure
D = 6;                      % inside diameter
V = 22.4;                   % volumetric flow rate of stack gas
Ts = 400;                   % temperature of stack gas  
Ta = 283;                   % temperature of ambient air  
xstack = 0; ystack = 60;    % stack location(s)
xmin = 0.5; xmax = 3.5;     % x-axis interval
ymin = 0; ymax = 100;       % y-axis interval (used only for d>1))
h = 50;                     % physical stack height
z = 0;                      % height of observation (=0 for ground surface)
gplot = 1;                  % plot option (=1 yes; =0 no)
gcont = 2;                  % contour plot option (=2 filled; =1 yes; =0 none) 
%----------------------------------execution-------------------------------
[x,y] = meshgrid (linspace(xmin,xmax,100),linspace(ymin,ymax,100));
c = zeros (size(x)); e = ones(size(x));
Dy = linspace(1,100,100); %  in meters
Dz = Dy'; %  in meters
[Dy,Dz] = meshgrid(Dz,Dy);
for i=1:100
   for j=1:100
   %...Pasquill-gifford for Dy
      c = 24.167; d = 2.5334;%...Pasquill Stability Category is A
      if x(i,j)<0.10; % x in kilometers
         th = 0.017453293.*(c - d.*log(x(i,j)));
         Dy(i,j) = 465.11628.*x(i,j).*(tan(th));
      end
   c = 18.3330; d = 1.8096;%...Pasquill Stability Category is B
   if x(i,j)<0.20 % x in kilometers
      th = 0.017453293.*(c - d.*log(x(i,j)));
      Dy(i,j) = 465.11628.*x(i,j).*(tan(th)); 
   end

   c = 12.5000; d = 1.0857;%...Pasquill Stability Category is C
   if x(i,j) == 4; % x in kilometers 
      th = 0.017453293.*(c - d.*log(x(i,j)));
      Dy(i,j) = 465.11628.*x(i,j).*(tan(th)); 
   end
   c = 8.3330; d = 0.72382;%...Pasquill Stability Category is D
   if x(i,j)== 0.30 % x in kilometers
      th = 0.017453293.*(c - d.*log(x(i,j)));
      Dy(i,j) = 465.11628.*x(i,j).*(tan(th)); 
   end

   c = 6.2500; d = 0.54287; %...Pasquill Stability Category is E
   if x(i,j) < 0.10 % x in kilometers 
      th = 0.017453293.*(c - d.*log(x(i,j)));
      Dy(i,j) = 465.11628.*x(i,j).*(tan(th)); 
   end

   c = 4.1667; d = 0.36191; %...Pasquill Stability Category is F
   if x(i,j) < 0.20 % x in kilometers
      th = 0.017453293.*(c - d.*log(x(i,j)));
      Dy(i,j) = 465.11628.*x(i,j).*(tan(th));
   end

%...Pasquill-gifford for Dz

   if x(i,j)<0.10 % x in kilometers
      a = 122.8; b = 0.9447;%...stabilitas is A
      Dz(i,j)= a.*x(i,j).^b; 

      if Dz(i,j)>5000
         Dz(i,j)=5000;
      end
   elseif 0.100 <= x(i,j) && x(i,j) < 0.150
      a = 158.08; b = 1.0542;
      Dz(i,j)= a.*x(i,j).^b;

      if Dz(i,j)>5000
         Dz(i,j)=5000;
      end
   elseif 0.150 <= x(i,j) && x(i,j) < 0.200
      a= 170.22 ; b= 1.0932 ;
      Dz(i,j)= a.*x(i,j).^b;

      if Dz(i,j)>5000
         Dz(i,j)=5000;
      end
   elseif 0.200 <= x(i,j) && x(i,j) < 0.250
      a= 179.52 ; b= 1.1262 ;
      Dz(i,j)= a.*x(i,j).^b;

      if Dz(i,j)>5000
         Dz(i,j)=5000;
      end
   elseif 0.250 <= x(i,j) && x(i,j) < 0.300
      a= 217.41 ; b= 1.2644 ;
      Dz(i,j)= a.*x(i,j).^b;

      if Dz(i,j)>5000
         Dz(i,j)=5000;
      end
   elseif 0.300 <= x(i,j) && x(i,j) < 0.400
a= 258.89 ; b= 1.4094 ;
Dz(i,j)= a.*x(i,j).^b; 

if Dz(i,j)>5000
Dz(i,j)=5000;
end
elseif 0.400 <= x(i,j) && x(i,j) < 0.500
a= 346.75 ; b= 1.7283 ;
Dz(i,j)= a.*x(i,j).^b;

if Dz(i,j)>5000
Dz(i,j)=5000;
end
elseif 0.500 <= x(i,j) && x(i,j) < 3.110
a= 453.85 ; b= 2.1166 ;
Dz(i,j)= a.*x(i,j).^b;

if Dz(i,j)>5000
Dz(i,j)=5000;
end
elseif x(i,j) >= 3.110
a= 453.85 ; b= 2.1166 ;
Dz(i,j)= a.*x(i,j).^b;  

if Dz(i,j)>5000
Dz(i,j)=5000;
end
if x(i,j)<0.20 % x in kilometers
a = 90.673; b = 0.93198;%...stabilitas is B
Dz(i,j)= a.*x(i,j).^b;  
end

if Dz(i,j)>5000
Dz(i,j)=5000;
end
elseif 0.21 <= x(i,j) && x(i,j) < 0.40
a = 98.483; b = 0.98332;
Dz(i,j)= a.*x(i,j).^b; 

if Dz(i,j)>5000
Dz = 5000;
end
elseif x(i,j)>= 0.40
a = 109.400; b = 1.09710;
Dz(i,j)= a.*x(i,j).^b;  
if Dz(i,j)>5000
Dz = 5000;
end

if x(i,j) == all % x in kilometers
a = 61.141; b = 0.91465;%...stabilitas is C
Dz(i,j)= a.*x(i,j).^b;  

if Dz(i,j)>5000
Dz = 5000;
end


if x(i,j)== 0.30 % x in kilometers
a = 34.459; b = 0.86974;%...stabilitas is D
Dz(i,j)= a.*x(i,j).^b;  
end

if 0.31 <= x(i,j) && x(i,j) < 1.00
a = 32.093; b = 0.81066;
Dz(i,j)= a.*x(i,j).^b;  
end

if 1.01 <= x(i,j) && x(i,j) < 3.00
a = 32.093; b = 0.64403;
Dz(i,j)= a.*x(i,j).^b;  
end

if 3.01 <= x(i,j) && x(i,j) < 10.00
a = 33.504; b = 0.60486;
Dz(i,j)= a.*x(i,j).^b;  
end
if 10.01 <= x(i,j) && x(i,j) < 30.00
a = 36.650; b = 0.56589;
Dz(i,j)= a.*x(i,j).^b;  
end
if x(i,j)>= 30.00
a = 44.053; b =0.51179;
Dz(i,j)= a.*x(i,j).^b;  
end

if x(i,j) < 0.10 % x in kilometers
a = 24.260; b = 0.83660; %...stabilitas is E  
Dz(i,j)= a.*x(i,j).^b; 
end

if 0.10 <= x(i,j) && x(i,j) < 30.00
a = 23.331; b = 0.81956;
Dz(i,j)= a.*x(i,j).^b;  
end

if 0.31 <= x(i,j) && x(i,j) < 1.00
a = 21.628; b = 0.75660;
Dz(i,j)= a.*x(i,j).^b;  
end

if 1.01 <= x(i,j) && x(i,j) < 2.00
a = 21.628; b = 0.63077;
Dz(i,j)= a.*x(i,j).^b;  
end
if 2.01 <= x(i,j) && x(i,j) < 4.00
a = 22.534; b = 0.57154;
Dz(i,j)= a.*x(i,j).^b;  
end
if 4.01 <= x(i,j) && x(i,j) < 10.00
a = 24.703; b = 0.50527;
Dz(i,j)= a.*x(i,j).^b;  
end
if 10.01 <= x(i,j) && x(i,j) < 20.00
a = 26.970; b = 0.46713;
Dz(i,j)= a.*x(i,j).^b;  
end
if 20.01 <= x(i,j) && x(i,j) < 40.00
a = 35.420; b = 0.37615;
Dz(i,j)= a.*x(i,j).^b;  
end
if x(i,j) >=40.00
a = 44.053; b = 0.51179;
Dz(i,j)= a.*x(i,j).^b;  
end

if x(i,j) < 0.20 % x in kilometers
a = 15.209; b = 0.81558;%...stabilitas is F
Dz(i,j)= a.*x(i,j).^b;  
end
if 0.21 <= x(i,j) && x(i,j) < 0.70
a = 14.457; b = 0.78407;
Dz(i,j)= a.*x(i,j).^b; 
end
if 0.71 <= x(i,j) && x(i,j) < 1.00
a = 13.953; b = 0.68465;
Dz(i,j)= a.*x(i,j).^b; 
end
if 1.01 <= x(i,j) && x(i,j) < 2.00
a = 13.953; b = 0.63227;
Dz(i,j)= a.*x(i,j).^b; 
end
if 2.01 <= x(i,j) && x(i,j) < 3.00
a = 14.823; b = 0.54503;
Dz(i,j)= a.*x(i,j).^b;  
end
if 3.01 <= x(i,j) && x(i,j) < 7.00
a = 16.187; b = 0.46490;
Dz(i,j)= a.*x(i,j).^b;  
end
if 7.01 <= x(i,j) && x(i,j) < 15.00
a = 17.836; b = 0.41507;
Dz(i,j)= a.*x(i,j).^b;  
end

if 15.01 <= x(i,j) && x(i,j) < 30.00
a = 22.651; b = 0.32681;
Dz(i,j)= a.*x(i,j).^b;  
end
if 30.01 <= x(i,j) && x(i,j) < 60.00
a = 27.074; b = 0.27436;
Dz(i,j)= a.*x(i,j).^b;  
end

if x(i,j)>= 60.00
a = 34.219; b = 0.21716;
Dz(i,j)= a.*x(i,j).^b;
end
end
end

end
end
for i = 1:size(Q,2)
    xx = x - xstack(i);
    yy = y - ystack(i); 
end

    deltah(i,j) = V.*D./u.*1.5 + 2.68.*(10^(-3)).*P.*(Ts-Ta).*D./Ts;
    H = h + deltah;
 c1 = c + Q(i).*e./(4.*pi.*xx.*sqrt(Dy.*Dz)).*exp(-v.*yy.*yy./(4.*Dy.*xx)).*... 
(exp(-v.*(z-H(i)).*(z-H(i)).*e./(4.*Dz.*xx))+exp(-v.*(z+H(i)).*(z+H(i)).*e./(4.*Dz.*xx)))...
.*exp(-lambda.*xx./v);
%--------------------------stack2---------------------------------------------
v = 0.5;                    % velocity
lambda = 0;                 % decay rate
Q = 1;                      % emission rate(s)
u = 3;                      % wind speed
P = 1;                      % Ambient Pressure
D = 6;                      % inside diameter
V = 22.4;                   % volumetric flow rate of stack gas
Ts = 400;                   % temperature of stack gas  
Ta = 283;                   % temperature of ambient air  
xstack = 0; ystack = 40;    % stack location(s)
xmin = 0.5; xmax = 3.5;     % x-axis interval
ymin = 0; ymax = 100;       % y-axis interval (used only for d>1))
h = 50;                     % physical stack height
z = 0;                      % height of observation (=0 for ground surface)
gplot = 1;                  % plot option (=1 yes; =0 no)
gcont = 2;                  % contour plot option (=2 filled; =1 yes; =0 none) 
%----------------------------------execution-------------------------------
[x,y] = meshgrid (linspace(xmin,xmax,100),linspace(ymin,ymax,100));
c = zeros (size(x)); e = ones(size(x));
Dy = linspace(1,100,100); %  in meters
Dz = Dy'; %  in meters
[Dy,Dz] = meshgrid(Dz,Dy);
for i=1:100
for j=1:100
%...Pasquill-gifford for Dy
c = 24.167; d = 2.5334;%...Pasquill Stability Category is A
if x(i,j)<0.10; % x in kilometers
th = 0.017453293.*(c - d.*log(x(i,j)));
Dy(i,j) = 465.11628.*x(i,j).*(tan(th));
end
c = 18.3330; d = 1.8096;%...Pasquill Stability Category is B
if x(i,j)<0.20 % x in kilometers
th = 0.017453293.*(c - d.*log(x(i,j)));
Dy(i,j) = 465.11628.*x(i,j).*(tan(th)); 
end

c = 12.5000; d = 1.0857;%...Pasquill Stability Category is C
if x(i,j) == 4; % x in kilometers 
th = 0.017453293.*(c - d.*log(x(i,j)));
Dy(i,j) = 465.11628.*x(i,j).*(tan(th)); 
end
c = 8.3330; d = 0.72382;%...Pasquill Stability Category is D
if x(i,j)== 0.30 % x in kilometers
th = 0.017453293.*(c - d.*log(x(i,j)));
Dy(i,j) = 465.11628.*x(i,j).*(tan(th)); 
end

c = 6.2500; d = 0.54287; %...Pasquill Stability Category is E
if x(i,j) < 0.10 % x in kilometers 
th = 0.017453293.*(c - d.*log(x(i,j)));
Dy(i,j) = 465.11628.*x(i,j).*(tan(th)); 
end

c = 4.1667; d = 0.36191; %...Pasquill Stability Category is F
if x(i,j) < 0.20 % x in kilometers
th = 0.017453293.*(c - d.*log(x(i,j)));
Dy(i,j) = 465.11628.*x(i,j).*(tan(th));
end

%...Pasquill-gifford for Dz

if x(i,j)<0.10 % x in kilometers
a = 122.8; b = 0.9447;%...stabilitas is A
Dz(i,j)= a.*x(i,j).^b; 

if Dz(i,j)>5000
Dz(i,j)=5000;
end
elseif 0.100 <= x(i,j) && x(i,j) < 0.150
a = 158.08; b = 1.0542;
Dz(i,j)= a.*x(i,j).^b;

if Dz(i,j)>5000
Dz(i,j)=5000;
end
elseif 0.150 <= x(i,j) && x(i,j) < 0.200
a= 170.22 ; b= 1.0932 ;
Dz(i,j)= a.*x(i,j).^b;

if Dz(i,j)>5000
Dz(i,j)=5000;
end
elseif 0.200 <= x(i,j) && x(i,j) < 0.250
a= 179.52 ; b= 1.1262 ;
Dz(i,j)= a.*x(i,j).^b;

if Dz(i,j)>5000
Dz(i,j)=5000;
end
elseif 0.250 <= x(i,j) && x(i,j) < 0.300
a= 217.41 ; b= 1.2644 ;
Dz(i,j)= a.*x(i,j).^b;

if Dz(i,j)>5000
Dz(i,j)=5000;
end
elseif 0.300 <= x(i,j) && x(i,j) < 0.400
a= 258.89 ; b= 1.4094 ;
Dz(i,j)= a.*x(i,j).^b; 

if Dz(i,j)>5000
Dz(i,j)=5000;
end
elseif 0.400 <= x(i,j) && x(i,j) < 0.500
a= 346.75 ; b= 1.7283 ;
Dz(i,j)= a.*x(i,j).^b;

if Dz(i,j)>5000
Dz(i,j)=5000;
end
elseif 0.500 <= x(i,j) && x(i,j) < 3.110
a= 453.85 ; b= 2.1166 ;
Dz(i,j)= a.*x(i,j).^b;

if Dz(i,j)>5000
Dz(i,j)=5000;
end
elseif x(i,j) >= 3.110
a= 453.85 ; b= 2.1166 ;
Dz(i,j)= a.*x(i,j).^b;  

if Dz(i,j)>5000
Dz(i,j)=5000;
end
if x(i,j)<0.20 % x in kilometers
a = 90.673; b = 0.93198;%...stabilitas is B
Dz(i,j)= a.*x(i,j).^b;  
end

if Dz(i,j)>5000
Dz(i,j)=5000;
end
elseif 0.21 <= x(i,j) && x(i,j) < 0.40
a = 98.483; b = 0.98332;
Dz(i,j)= a.*x(i,j).^b; 

if Dz(i,j)>5000
Dz = 5000;
end
elseif x(i,j)>= 0.40
a = 109.400; b = 1.09710;
Dz(i,j)= a.*x(i,j).^b;  
if Dz(i,j)>5000
Dz = 5000;
end

if x(i,j) == all % x in kilometers
a = 61.141; b = 0.91465;%...stabilitas is C
Dz(i,j)= a.*x(i,j).^b;  

if Dz(i,j)>5000
Dz = 5000;
end


if x(i,j)== 0.30 % x in kilometers
a = 34.459; b = 0.86974;%...stabilitas is D
Dz(i,j)= a.*x(i,j).^b;  
end

if 0.31 <= x(i,j) && x(i,j) < 1.00
a = 32.093; b = 0.81066;
Dz(i,j)= a.*x(i,j).^b;  
end

if 1.01 <= x(i,j) && x(i,j) < 3.00
a = 32.093; b = 0.64403;
Dz(i,j)= a.*x(i,j).^b;  
end

if 3.01 <= x(i,j) && x(i,j) < 10.00
a = 33.504; b = 0.60486;
Dz(i,j)= a.*x(i,j).^b;  
end
if 10.01 <= x(i,j) && x(i,j) < 30.00
a = 36.650; b = 0.56589;
Dz(i,j)= a.*x(i,j).^b;  
end
if x(i,j)>= 30.00
a = 44.053; b =0.51179;
Dz(i,j)= a.*x(i,j).^b;  
end

if x(i,j) < 0.10 % x in kilometers
a = 24.260; b = 0.83660; %...stabilitas is E  
Dz(i,j)= a.*x(i,j).^b; 
end

if 0.10 <= x(i,j) && x(i,j) < 30.00
a = 23.331; b = 0.81956;
Dz(i,j)= a.*x(i,j).^b;  
end

if 0.31 <= x(i,j) && x(i,j) < 1.00
a = 21.628; b = 0.75660;
Dz(i,j)= a.*x(i,j).^b;  
end

if 1.01 <= x(i,j) && x(i,j) < 2.00
a = 21.628; b = 0.63077;
Dz(i,j)= a.*x(i,j).^b;  
end
if 2.01 <= x(i,j) && x(i,j) < 4.00
a = 22.534; b = 0.57154;
Dz(i,j)= a.*x(i,j).^b;  
end
if 4.01 <= x(i,j) && x(i,j) < 10.00
a = 24.703; b = 0.50527;
Dz(i,j)= a.*x(i,j).^b;  
end
if 10.01 <= x(i,j) && x(i,j) < 20.00
a = 26.970; b = 0.46713;
Dz(i,j)= a.*x(i,j).^b;  
end
if 20.01 <= x(i,j) && x(i,j) < 40.00
a = 35.420; b = 0.37615;
Dz(i,j)= a.*x(i,j).^b;  
end
if x(i,j) >=40.00
a = 44.053; b = 0.51179;
Dz(i,j)= a.*x(i,j).^b;  
end

if x(i,j) < 0.20 % x in kilometers
a = 15.209; b = 0.81558;%...stabilitas is F
Dz(i,j)= a.*x(i,j).^b;  
end
if 0.21 <= x(i,j) && x(i,j) < 0.70
a = 14.457; b = 0.78407;
Dz(i,j)= a.*x(i,j).^b; 
end
if 0.71 <= x(i,j) && x(i,j) < 1.00
a = 13.953; b = 0.68465;
Dz(i,j)= a.*x(i,j).^b; 
end
if 1.01 <= x(i,j) && x(i,j) < 2.00
a = 13.953; b = 0.63227;
Dz(i,j)= a.*x(i,j).^b; 
end
if 2.01 <= x(i,j) && x(i,j) < 3.00
a = 14.823; b = 0.54503;
Dz(i,j)= a.*x(i,j).^b;  
end
if 3.01 <= x(i,j) && x(i,j) < 7.00
a = 16.187; b = 0.46490;
Dz(i,j)= a.*x(i,j).^b;  
end
if 7.01 <= x(i,j) && x(i,j) < 15.00
a = 17.836; b = 0.41507;
Dz(i,j)= a.*x(i,j).^b;  
end

if 15.01 <= x(i,j) && x(i,j) < 30.00
a = 22.651; b = 0.32681;
Dz(i,j)= a.*x(i,j).^b;  
end
if 30.01 <= x(i,j) && x(i,j) < 60.00
a = 27.074; b = 0.27436;
Dz(i,j)= a.*x(i,j).^b;  
end

if x(i,j)>= 60.00
a = 34.219; b = 0.21716;
Dz(i,j)= a.*x(i,j).^b;
end
end
end

end
end
for i = 1:size(Q,2)
    xx = x - xstack(i);
    yy = y - ystack(i); 
end

    deltah(i,j) = V.*D./u.*1.5 + 2.68.*(10^(-3)).*P.*(Ts-Ta).*D./Ts;
    H = h + deltah;
 c2 = c + Q(i).*e./(4.*pi.*xx.*sqrt(Dy.*Dz)).*exp(-v.*yy.*yy./(4.*Dy.*xx)).*... 
(exp(-v.*(z-H(i)).*(z-H(i)).*e./(4.*Dz.*xx))+exp(-v.*(z+H(i)).*(z+H(i)).*e./(4.*Dz.*xx)))...
.*exp(-lambda.*xx./v);

%----------------------------------output----------------------------------
if gplot
    for i = 10:10:100
        plot (c1(:,i)); hold on;
        plot (c2(:,i)); hold on;
    end
end
if gcont
    figure;
    if gcont > 1
        contourf (x,y,c); colorbar;
    else
        contour (x,y,c); 
    end
end

有人可以帮我吗?!

当我运行您发布的代码时,我发现"c"变量是单个值,而"x"和"y"是矩阵。所以你的问题在于你试图一起绘制的变量的维度。

若要调试代码,请尝试在发生错误的代码中插入断点。然后,当您运行代码时,它将在断点处停止,您可以检查变量的值和维度。

最新更新