%
%  This script plots output from the simple hurricane model and allows one to
%    make .png files of graphs. 
%
%  For UNIX and LINUX systems, set 'sepchar' to '/'; for windows, set to '\'
%
sepchar='\';
%
file1=input('Enter output file (default = output)','s');
if isempty(file1)
   file1='output';
end   
file=strcat(file1,sepchar);
load(strcat(file,'time.out'))
load(strcat(file,'radius.out'))
rb=radius(:,1);
rbi=radius(:,2);
rt=radius(:,3);
time1=time(:,1);
load(strcat(file,'rgraph.out'))
load(strcat(file,'zgraph.out'))
r=rgraph;
z=zgraph;
as=size(rb);
aa=as(1,1);
itlast=0;
name={'Maximum_surface_wind','Radius_of_maximum_surface_wind','Minimum_V_at_tropopause','Radius_of_minimum V',...
      'Central_pressure','Updraft_mass_flux','Downdraft_mass_flux',...
      'Storm_center_ocean_mixed_layer_depth','Storm_center_SST_anomaly',...
      'Surface_azimuthal_velocity','Tropopause_azimuthal_velocity',...
      'Surface_radial_velocities','Vertical_mass_fluxes','Surface_pressure',...
      'Equivalent_potential_temperature','Ocean_mixed_layer_depth',...
      'Sea_surface_temperature_anomaly','Azimuthal_velocity','Radial_velocity',...
      'Vertical_velocity','Streamfunction','Perturbation_pressure',...
      'Perturbation_temperature','Equivalent_potential_temperature_(rz)',...
      'Rain_water','Updraft_mass_flux_(hov)','Azimuthal_velocity_(hov)',...
      'Radial_velocity_(hov)'};
i=1;
while i == 1
   ir=menu('Select from the following:','Time series','Radial plots',...
      'Contour plots','Hovmuller plots','Exit');
j=1;
while j == 1
if ir == 1,
   it=menu('Time series plots','Maximum surface wind','Radius of maximum surface wind','Minimum V at tropopause','Radius of minimum V',...
      'Central pressure','Updraft mass flux','Downdraft mass flux',...
      'Storm center ocean mixed layer depth','Storm center SST anomaly',...
      'Make .png plot','Exit');
   if it == 10
      fig=strcat(name(itlast),'.png');
      fig=char(fig);
      exportfig(gcf,fig,'format','png','color','rgb','bounds','loose','Renderer','Painters','Resolution',300);
      it=itlast;
   end
   if it == 11,
      break
   end  
elseif ir == 2
   it=menu('Radial plots','Surface azimuthal velocity','Tropopause azimuthal velocity',...
      'Radial velocities','Vertical mass fluxes','Surface pressure',...
      'Equivalent potential temperature','Ocean mixed layer depth',...
      'Sea surface temperature anomaly','Make .png plot','Exit');
   it=it+9;
   if it == 18
      fig=strcat(name(itlast),'.png');
      fig=char(fig);
      exportfig(gcf,fig,'format','png','color','rgb','bounds','loose','Renderer','Painters','Resolution',300);
      it=itlast;
   end
   if it == 19
      break
   end 
elseif ir == 3
   it=menu('Contour plots','Azimuthal velocity','Radial velocity',...
      'Vertical velocity','Streamfunction','Perturbation pressure',...
      'Perturbation temperature','Equivalent potential temperature',...
      'Rain water','Make .png plot','Exit');
   it=it+17;
   if it == 26
      fig=strcat(name(itlast),'.png');
      fig=char(fig);
      exportfig(gcf,fig,'format','png','color','rgb','bounds','loose','Renderer','Painters','Resolution',300);
      it=itlast;
   end
   if it == 27
      break
   end
elseif ir == 4
   it=menu('Hovmuller plots','Updraft mass flux','Azimuthal velocity',...
      'Radial velocity','Make .png plot','Exit');
   it=it+25;
   if it == 29
      fig=strcat(name(itlast),'.png');
      fig=char(fig);
      exportfig(gcf,fig,'format','png','color','rgb','bounds','loose','Renderer','Painters','Resolution',300);
      it=itlast;
   end
   if it == 30
      break
   end
elseif ir == 5
   i=0;
   it=0;
   break
end
itlast=it;
if it == 0
   break
elseif it == 1
v=time(:,2);
vtheory=time(:,3);
h=plot(time1,v,'-',time1,vtheory,'--');
%h=plot(time1,vcon,'-',time1,vtheory,'--',time1,v,'-');
set(h,'linewidth',3)
set(gca,'fontweight','bold','fontsize',12)
xlabel('Time (days)','fontweight','bold')         
ylabel('Maximum surface wind speed (m/s)','fontweight','bold')
elseif it == 2
rmax=time(:,4);
h=plot(time1,rmax);
set(h,'linewidth',3)
set(gca,'fontweight','bold','fontsize',12)
xlabel('Time (days)','fontweight','bold')         
ylabel('Radius of maximum winds (km)','fontweight','bold')
elseif it == 3
vt=time(:,5);
h=plot(time1,vt);
set(h,'linewidth',3)
set(gca,'fontweight','bold','fontsize',12)
xlabel('Time (days)','fontweight','bold')         
ylabel('Minimum tropopause wind speed (m/s)','fontweight','bold')
elseif it == 4
rt2=time(:,6);
h=plot(time1,rt2);
set(h,'linewidth',3)
set(gca,'fontweight','bold','fontsize',12)
xlabel('Time (days)','fontweight','bold')         
ylabel('Radius of minimum tropopause winds (km)','fontweight','bold')
elseif it == 5
pc=time(:,7);
ptheory=time(:,8);
h=plot(time1,pc,'-',time1,ptheory,'--');
set(h,'linewidth',3)
set(gca,'fontweight','bold','fontsize',12)
xlabel('Time (days)','fontweight','bold')         
ylabel('Central surface pressure (mb)','fontweight','bold')
elseif it == 6
mu=time(:,9);
h=plot(time1,mu);
set(h,'linewidth',3)
set(gca,'fontweight','bold','fontsize',12)
xlabel('Time (days)','fontweight','bold')         
ylabel('Convective updraft mass flux (m/s)','fontweight','bold')
elseif it == 7
md=time(:,10);
h=plot(time1,md);
set(h,'linewidth',3)
set(gca,'fontweight','bold','fontsize',12)
xlabel('Time (days)','fontweight','bold')         
ylabel('Convective downdraft mass flux (m/s)','fontweight','bold')
elseif it == 8
hmc=time(:,11);
h=plot(time1,hmc);
set(h,'linewidth',3)
set(gca,'fontweight','bold','fontsize',12)
xlabel('Time (days)','fontweight','bold')
ylabel('Storm center ocean mixed layer depth (m)','fontweight','bold')
elseif it == 9
sstc=time(:,12);
h=plot(time1,sstc);
set(h,'linewidth',3)
set(gca,'fontweight','bold','fontsize',12)
xlabel('Time (days)','fontweight','bold')
ylabel('Storm center SST anomaly (K)','fontweight','bold')
elseif it == 10
v=radius(:,4);
h=plot(rb,v);
set(h,'linewidth',3)
set(gca,'fontweight','bold','fontsize',12)
xlabel('Radius (km)','fontweight','bold')
ylabel('Surface azimuthal velocity (m/s)','fontweight','bold')
elseif it == 11
vt=radius(:,5);
h=plot(rt,vt);
set(h,'linewidth',3)
set(gca,'fontweight','bold','fontsize',12)
xlabel('Radius (km)','fontweight','bold')
ylabel('Tropopause azimuthal velocity (m/s)','fontweight','bold')
elseif it == 14
p=radius(:,6);
h=plot(rbi,p);
set(h,'linewidth',3)
set(gca,'fontweight','bold','fontsize',12)
xlabel('Radius (km)','fontweight','bold')
ylabel('Surface pressure (mb)','fontweight','bold')
elseif it == 12
u0=radius(:,7);
u=radius(:,8);
h=plot(rb,u0,'-',rb,u,'--');
set(h,'linewidth',3)
set(gca,'fontweight','bold','fontsize',12)
xlabel('Radius (km)','fontweight','bold')
ylabel('Radial velocities (m/s)','fontweight','bold')
elseif it == 13
mu=radius(:,9);
we=radius(:,10);
h=plot(rbi,mu,'-',rbi,we,'--');
set(h,'linewidth',3)
set(gca,'fontweight','bold','fontsize',12)
xlabel('Radius (km)','fontweight','bold')
ylabel('Updraft and subsidence velocities (m/s)','fontweight','bold')
elseif it == 15
the=radius(:,11);
thes=radius(:,12);
them=radius(:,13);
h=plot(rbi,the,'-',rbi,thes,'--',rbi,them,'-.');
set(h,'linewidth',3)
set(gca,'fontweight','bold','fontsize',12)
xlabel('Radius (km)','fontweight','bold')
ylabel('PBL, saturation and mid-trop theta_e','fontweight','bold')
elseif it == 16
hm=radius(:,14);
h=plot(rbi,hm);
set(h,'linewidth',3)
set(gca,'fontweight','bold','fontsize',12)
xlabel('Radius (km)','fontweight','bold')
ylabel('Mixed layer depth (m)','fontweight','bold')
elseif it == 17
sst=radius(:,15);
h=plot(rbi,sst);
set(h,'linewidth',3)
set(gca,'fontweight','bold','fontsize',12)
xlabel('Radius (km)','fontweight','bold')
ylabel('Sea surface temperature anomaly (K)','fontweight','bold')
elseif it == 18
load(strcat(file,'vcon.out'))
b=max(max(vcon));
a=min(min(vcon));
pcolor(r,z,vcon)
colormap(jet)
%pcolor(r,z,vcon)
set(gca,'fontweight','bold','fontsize',12)
title(['Azimuthal velocity, from ',num2str(a),' to ',num2str(b), ' m/s'],'fontweight','bold')
xlabel('Radius (km)','fontweight','bold')
ylabel ('z (km)','fontweight','bold')
shading interp
g1=colorbar('horiz');
set(g1,'fontweight','bold')
elseif it == 19
load(strcat(file,'ucon.out'))
b=max(max(ucon));
a=min(min(ucon));
pcolor(r,z,ucon)
colormap(hsv)
shading interp
set(gca,'fontweight','bold','fontsize',12)
title(['Radial velocity, from ',num2str(a),' to ',num2str(b),' m/s'],'fontweight','bold')
xlabel('Radius (km)','fontweight','bold')
ylabel ('z (km)','fontweight','bold')
g1=colorbar('horiz');
set(g1,'fontweight','bold')
elseif it == 20
load(strcat(file,'wcon.out'))
b=max(max(wcon));
a=min(min(wcon));
pcolor(r,z,wcon)
colormap(jet)
shading interp
g1=colorbar('horiz');
set(g1,'fontweight','bold')
set(gca,'fontweight','bold','fontsize',12)
title(['Vertical velocity, from ',num2str(a),' to ',num2str(b),' m/s (- values X 100)'],'fontweight','bold')
xlabel('Radius (km)','fontweight','bold')
ylabel ('z (km)','fontweight','bold')
elseif it == 21
load(strcat(file,'psicon.out'))
b=max(max(psicon));
a=min(min(psicon));
pcolor(r,z,psicon)
colormap(jet)
shading interp
set(gca,'fontweight','bold','fontsize',12)
title(['Streamfunction, from ',num2str(a),' to ',num2str(b),' 10**8 Kg/s'],'fontweight','bold')
xlabel('Radius (km)','fontweight','bold')
ylabel ('z (km)','fontweight','bold')
hold on
contour(r,z,psicon,'k')
hold off
g1=colorbar('horiz');
set(g1,'fontweight','bold')
elseif it == 22
clf
load(strcat(file,'pcon.out'))
b=max(max(pcon));
a=min(min(pcon));
pcolor(r,z,pcon)
colormap(jet)
shading interp
set(gca,'fontweight','bold','fontsize',12)
title(['Perturbation pressure, from ',num2str(a),' to ',num2str(b),' mb'],'fontweight','bold')
xlabel('Radius (km)','fontweight','bold')
ylabel ('z (km)','fontweight','bold')
hold on
contour(r,z,pcon,'k')
hold off
g1=colorbar('horiz');
set(g1,'fontweight','bold')
elseif it == 23
load(strcat(file,'tcon.out'))
b=max(max(tcon));
a=min(min(tcon));
colormap(jet)
pcolor(r,z,tcon)
shading interp
set(gca,'fontweight','bold','fontsize',12)
title(['Perturbation temperature, from ',num2str(a),' to ',num2str(b),' K'],'fontweight','bold')
xlabel('Radius (km)','fontweight','bold')
ylabel ('z (km)','fontweight','bold')
g1=colorbar('horiz');
set(g1,'fontweight','bold')
elseif it == 24
load(strcat(file,'tecon.out'))
b=max(max(tecon));
a=min(min(tecon));
pcolor(r,z,tecon)
colormap(jet)
shading interp
set(gca,'fontweight','bold','fontsize',12)
title(['Equivalent potential temperature, from ',num2str(a),' to ',num2str(b),' K'],'fontweight','bold')
xlabel('Radius (km)','fontweight','bold')
ylabel ('z (km)','fontweight','bold')
g1=colorbar('horiz');
set(g1,'fontweight','bold')
elseif it == 25
load(strcat(file,'raincon.out'))
   lrain=log(raincon+1);
b=max(max(raincon));
a=min(min(raincon));
pcolor(r,z,lrain)
colormap(jet)
shading interp
set(gca,'fontweight','bold','fontsize',12)
title(['Log of Rain water content, from ',num2str(a),' to ',num2str(b),' g/Kg'],'fontweight','bold')
xlabel('Radius (km)','fontweight','bold')
ylabel ('z (km)','fontweight','bold')
g1=colorbar('horiz');
set(g1,'fontweight','bold')
elseif it == 26
load(strcat(file,'whov.out'))
b=max(max(whov));
a=min(min(whov));
pcolor(r,time1,whov)
colormap(jet)
shading interp
set(gca,'fontweight','bold','fontsize',12)
title(['Cumulus mass flux, from ',num2str(a),' to ',num2str(b),' m/s'],'fontweight','bold')
xlabel('Radius (km)','fontweight','bold')
ylabel ('time (days)','fontweight','bold')
g1=colorbar('horiz');
set(g1,'fontweight','bold')
elseif it == 27
load(strcat(file,'vhov.out'))
b=max(max(vhov));
a=min(min(vhov));
pcolor(r,time1,vhov)
colormap(jet)
shading interp
set(gca,'fontweight','bold','fontsize',12)
title(['Azimuthal velocity, from ',num2str(a),' to ',num2str(b),' m/s'],'fontweight','bold')
xlabel('Radius (km)','fontweight','bold')
ylabel ('time (days)','fontweight','bold')
g1=colorbar('horiz');
set(g1,'fontweight','bold')
elseif it == 28
load(strcat(file,'uhov.out'))
b=max(max(uhov));
a=min(min(uhov));
pcolor(r,time1,uhov)
colormap(jet)
shading interp
set(gca,'fontweight','bold','fontsize',12)
title(['Radial velocity, from ',num2str(a),' to ',num2str(b),' m/s'],'fontweight','bold')
xlabel('Radius (km)','fontweight','bold')
ylabel ('time (days)','fontweight','bold')
g1=colorbar('horiz');
set(g1,'fontweight','bold')
end
end
end
