ps format file does not display all plots

Hi there,

I am trying to plot vertical profiles of some quantities. For example
vertical profile of potential temperature for the first 15 levels (from
ground); i want these vertical profiles to be plotted for 3 hours, each
hour separately (for example vertical profile of potential temperature
for 49, 50 and 51 hour). I am using a 'do' loop in my script for this
purpose. The script is working fine for x11 windows; that is, it is
plotting vertical profiles of various quantities separately for every
hour. However when i generate post script 'ps' file, the .ps file only
creates the vertical profile for the last hour. It seems postscript
format does not like the loop. On the other hand if i use 'ncgm' option
for output, it creates several 'ncgm' files, one file for every hour. I
want to generate only one file that contains all the plots. I
appreciate, if anybody could help me how to handle this problem. My
script is attached and I am using ncl version 5.0.

Regards,

-- 
---
Basit A. Khan
PhD Candidate
Centre for Atmospheric Research
Department of Geography
University of Canterbury  Te Whare Wananga o Waitaha Private Bag 4800 
Christchurch 8020, New Zealand
Tel: (643) 364 2987 Ext. 7912
Fax: (643) 364 2907
bak41_at_student.canterbury.ac.nz 

load "$NCARG_ROOT/lib/ncarg/nclex/gsun/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRF_contributed.ncl"

begin
  
  cdf_file = addfile("wrfout_d04_2006-11-01_12:00:00.nc","r") ; %%%%%%%%%%%%%%%%%% USER INPUT REQUIRED %%%%%%%%%%%%%%
  
  
  
  times = wrf_user_list_times(cdf_file)
  ntimes = dimsizes(times)
  v_levels=(/0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35,36,37,38,39,40/)

  
  hour=new(ntimes,integer) ;creating an array with number of subscripts equal to number of times in the wrf file
  
  do counter=0,ntimes-1,1 ; running loop to fill 'hour' with the times in the wrf file
     hour(counter)=counter
  end do

  SL = 0 ;starting vertical level %%%%%%%%%%%%%%%%%% USER INPPUT REQUIRED %%%%%%%%%%%%%%%%
  EL = 3 ;Endig vertical level %%%%%%%%%%%%%%%%% USER INPPUT REQUIRED %%%%%%%%%%%%%%%%%
  ST = 48 ;Start time %%%%%%%%%%%%%%% USER INPPUT REQUIRED %%%%%%%%%%%%%%%%%%%
  ET = 50 ;End Time %%%%%%%%%%%%%%%%%% USER INPPUT REQUIRED %%%%%%%%%%%%%%%%
  xcord = 49 ;X-coordinate %%%%%%%%%%%%%%%%%%% USER INPPUT REQUIRED %%%%%%%%%%%%%%%
  ycord = 49 ;Y-coordinate %%%%%%%%%%%%%%%%%%%% USER INPPUT REQUIRED %%%%%%%%%%%%%%
  Time = hour(ST:ET)

  sub_ntimes=dimsizes(Time)
  VProfile= v_levels(SL:EL)
 
 ; XLAT = cdf_file->XLAT ; latitude
 ; XLONG = cdf_file->XLONG ; longitude
 
 do it =0, sub_ntimes -1, 1

; !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
; 3 D QUANTITIES - GET THE VARIABLES FROM NETCDF FILE

  U = cdf_file->U(ST,SL:EL,xcord,ycord) ; x-wind component from 48th hour to 71 hour x=2 and y=49
  V = cdf_file->V(ST,SL:EL,xcord,ycord) ; y-wind component
  W = cdf_file->W(ST,SL:EL,xcord,ycord) ; z-wind component
  T = cdf_file->T(ST,SL:EL,xcord,ycord) ; perturbation potential temperature (theta-t0)
  P = cdf_file->P(ST,SL:EL,xcord,ycord) ; perturbation pressure
  PH = cdf_file->PH(ST,SL:EL,xcord,ycord) ; perturbation geopotential
  PB = cdf_file->PB(ST,SL:EL,xcord,ycord) ; BASE STATE PRESSURE
  PHB = cdf_file->PHB(ST,SL:EL,xcord,ycord) ; base-state geopotential
  QVAPOR = cdf_file->QVAPOR(ST,SL:EL,xcord,ycord) ; Water vapor mixing ratio (Kg/kg)
 
  if (EL .le.3) then
    TSLB = cdf_file->TSLB(ST,SL:EL,xcord,ycord) ; SOIL TEMPERATURE
    SMOIS = cdf_file->SMOIS(ST,SL:EL,xcord,ycord) ; SOIL MOISTURE
  end if
  

; !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
; 3D DIAGNOSTIC QUANTITIES - GET THE VARIABLES

   THk = T + 300
   PRES = (P + PB)/100
   Tk = (T+300)*((P+PB)/100000)^0.286
   UV = ((U^2)+(V^2))^(1/2) ; horizontal wind speed
   UVW = ((U^2)+(V^2)+(W^2))^(1/2) ; Three dimensional wind speed
   RH = wrf_rh(QVAPOR, PRES, Tk)
   Td = wrf_td(PRES, QVAPOR)
   Z = (PH+PHB)/9.81
            dim = dimsizes(Z)
   GHT = 0.5*(Z(0:dim(0)-2)+Z(1:dim(0)-1))
   
   ST = ST+1

    Time_at_long_name = "Time (NZST)"
    U_at_long_name = "U m/s"
    V_at_long_name ="V m/s"
    W_at_long_name ="W m/s"
    QVAPOR_at_long_name = "QVAPOR (kg/kg)"
    Tk_at_long_name = "Temp degK"
    Td_at_long_name = "Dew Point degC"
    RH_at_long_name = "RH %"
    THk_at_long_name = "TH degK"
    UV_at_long_name = "UV m/s"
    UVW_at_long_name = "UVW m/s"
    PRES_at_long_name = "PRES (hPa)"
    GHT_at_long_name = "GHT (m)"
    Z_at_long_name = "Z (m)"
  
 if (EL .le. 3) then
   TSLB_at_long_name = "TSLB (K)"
   SMOIS_at_long_name = "SMOIS (m3/m3)"
 end if
    
  wks = gsn_open_wks("ps","vp_new") ; %%%%%%%%%%%%%%%%%%%%%%%%%%%% Open an X11 workstation.
  
  res = True
  res_at_tmXBLabelStride= 1
  res_at_tiMainString = "Vertical Profile; Hour: " + ST
  res_at_tiYAxisString = "Sigma Lvls"
  
 if (EL .le.3) then
    plot = gsn_xy(wks,TSLB,VProfile,res)
    plot = gsn_xy(wks,SMOIS,VProfile,res)
 end if
   
  plot = gsn_xy(wks,W,VProfile,res)
  plot = gsn_xy(wks,Tk,VProfile,res)
  plot = gsn_xy(wks,Td,VProfile,res)
  plot = gsn_xy(wks,THk,VProfile,res)
  plot = gsn_xy(wks,PRES,VProfile,res)
  plot = gsn_xy(wks,Z,VProfile,res)
  plot = gsn_xy(wks,GHT,VProfile,res)
  plot = gsn_xy(wks,RH,VProfile,res)
  plot = gsn_xy(wks,QVAPOR,VProfile,res)

; !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
; Plotting U and V verses time series

  data_UV = new((/2,dimsizes(U)/),float)
  data_UV(0,:) = U
  data_UV(1,:) = V

;;;;;;;;;;;; Resources for U and V plot
  res = True
  res_at_tiMainString = "Vertical Profile; Hour: " + ST
  res_at_tiYAxisString = "Sigma Lvls"
  res_at_xyLineColors = (/"black","red"/)
  res_at_xyLabelMode = "Custom"
  res_at_xyExplicitLabels= (/"U","V"/)
  res_at_xyLineLabelFontHeightF = 0.01
  res_at_xyLineLabelFontColor =2
    plot = gsn_xy(wks,data_UV,VProfile,res)
  delete(plot)
  delete(res)

; !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
; Plotting UV(horizontal windspeed) and UVW(3D wind speed) verses time series

  data_UVUVW = new((/2,dimsizes(UV)/),float)
  data_UVUVW(0,:) = UV
  data_UVUVW(1,:) = UVW

;;;;;;;;;;;; Resources for UV and UVW plot
  res = True
  res_at_tiMainString = "Vertical Profile; Hour: " + ST
  res_at_tiYAxisString = "Sigma Lvls"
  res_at_xyLineColors = (/"blue","red"/)
  res_at_xyLabelMode = "Custom"
  res_at_xyExplicitLabels= (/"UV","UVW"/)
  res_at_xyLineLabelFontHeightF = 0.01
  res_at_xyLineLabelFontColor =2
    plot = gsn_xy(wks,data_UVUVW,VProfile,res)
  delete(plot)
  delete(res)
 
 end do

end

_______________________________________________
ncl-talk mailing list
ncl-talk_at_ucar.edu
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Sun Jun 22 2008 - 16:41:11 MDT