!
!##################################################################
!##################################################################
!###### ######
!###### Advanced Regional Prediction System (ARPS) ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
!
!
PROGRAM PLT_RADAR_TILT,4
!
!-----------------------------------------------------------------------
!
! PURPOSE:
! Plot radar reflectivity and radial velocity in a tilt and
! vertical cross-sections
!
!-----------------------------------------------------------------------
!
! AUTHOR: Ming Xue
! 03/19/02.
!
! MODIFICATION HISTORY:
! 6/12/2005 Keith Brewster, CAPS
! Added reading of radar emulator volume files (radfmt=2)
!
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
!
IMPLICIT NONE
INTEGER, PARAMETER :: mxmapfile= 100
INTEGER, PARAMETER :: nref_tilts_max=14
INTEGER, PARAMETER :: nvel_tilts_max=14
INTEGER, PARAMETER :: n_elevation_max=14
INTEGER, PARAMETER :: n_azimuth_max=361
!
!-----------------------------------------------------------------------
!
! Namelist
!
!-----------------------------------------------------------------------
!
CHARACTER (LEN=132) :: filepath
CHARACTER (LEN=132) :: filetitle
INTEGER :: radfmt
INTEGER :: numfiles
CHARACTER (LEN=132) :: files(mxmapfile)
NAMELIST /input_data/ filepath,filetitle,radfmt,numfiles,files
CHARACTER(LEN=4) :: radar_symbol
REAL :: ctrlat
REAL :: ctrlon
NAMELIST /radar_index/ radar_symbol,ctrlat,ctrlon
INTEGER :: iout_ref_regular
INTEGER :: iout_vel_regular
CHARACTER (LEN=132) :: ref_regular_file
CHARACTER (LEN=132) :: vel_regular_file
CHARACTER (LEN=132) :: filetrailer
NAMELIST /output/ iout_ref_regular,iout_vel_regular,filetrailer
INTEGER :: mapproj
REAL :: trulat1
REAL :: trulat2
REAL :: trulon
REAL :: sclfct
NAMELIST /projection/ mapproj,trulat1,trulat2,trulon,sclfct
INTEGER :: overlayrefcontour
INTEGER :: n_elevation_ref
INTEGER :: n_elevation_vel
INTEGER :: n_azimuth_ref
INTEGER :: n_azimuth_vel
REAL :: elevations_ref(n_elevation_max)
REAL :: elevations_vel(n_elevation_max)
REAL :: azimuth_ref(n_azimuth_max)
REAL :: azimuth_vel(n_azimuth_max)
REAL :: h_range_max,h_range_max_e,h_range_max_w,h_range_max_s,h_range_max_n
REAL :: z_range_max
REAL :: elev_error_range
INTEGER :: coltab_ref
INTEGER :: ibgncol_ref
INTEGER :: iendcol_ref
CHARACTER(LEN=80) :: coltabfn_ref
INTEGER :: coltab_vel
INTEGER :: ibgncol_vel
INTEGER :: iendcol_vel
CHARACTER(LEN=80) :: coltabfn_vel
NAMELIST /plotting/ h_range_max_e,h_range_max_w,h_range_max_s, &
h_range_max_n,z_range_max,elev_error_range, &
overlayrefcontour, &
n_elevation_ref,elevations_ref, &
n_elevation_vel,elevations_vel, &
n_azimuth_ref,azimuth_ref, &
n_azimuth_vel,azimuth_vel, &
coltab_ref,ibgncol_ref,iendcol_ref,coltabfn_ref, &
coltab_vel,ibgncol_vel,iendcol_vel,coltabfn_vel
INTEGER :: ovrmap
INTEGER :: mapgrid
INTEGER :: nmapfile
CHARACTER (LEN=132) :: mapfile(mxmapfile)
REAL :: latgrid
REAL :: longrid
NAMELIST /map_plot/ ovrmap,mapgrid,latgrid,longrid,nmapfile,mapfile
INTEGER :: lmapfile
REAL :: xorig, yorig
!
!-----------------------------------------------------------------------
!
! plot variable
!
!-----------------------------------------------------------------------
!
REAL :: cl(100)
REAL, ALLOCATABLE :: x(:,:)
REAL, ALLOCATABLE :: y(:,:)
REAL, ALLOCATABLE :: zp(:,:)
REAL, ALLOCATABLE :: radialp(:,:)
REAL, ALLOCATABLE :: zm(:,:)
REAL, ALLOCATABLE :: radialm(:,:)
REAL, ALLOCATABLE :: xw(:)
REAL, ALLOCATABLE :: yw(:)
REAL, ALLOCATABLE :: ref_rz(:,:)
REAL, ALLOCATABLE :: vel_rz(:,:)
INTEGER, ALLOCATABLE :: iw(:,:)
INTEGER :: ncl
INTEGER :: mode
!
!-----------------------------------------------------------------------
!
! INPUT RADAR PARAMETER AND OBSERVATION
!
!-----------------------------------------------------------------------
!
REAL :: radar_lat
REAL :: radar_lon
REAL :: radar_elevation
INTEGER :: iyear,imon,iday,ihour,imin,isec
INTEGER, ALLOCATABLE :: iryear(:),irmon(:),irday(:),irhour(:),irmin(:),irsec(:)
INTEGER, ALLOCATABLE :: ivyear(:),ivmon(:),ivday(:),ivhour(:),ivmin(:),ivsec(:)
INTEGER :: nrefgate
INTEGER :: nrefray
REAL :: rfstgat
REAL :: refgatsp
INTEGER :: nvelgate
INTEGER :: nvelray
REAL :: vfstgat
REAL :: velgatsp
REAL :: azmlutint
CHARACTER(LEN=80) :: runname
INTEGER :: lfnkey,abstsec
INTEGER :: wrtuaref,wrtuavel,wrtvort,idummy,vcp
REAL :: beamwid,rmisval,rngfval,curtim
REAL, ALLOCATABLE :: razim(:)
REAL, ALLOCATABLE :: razim_tem(:)
REAL, ALLOCATABLE :: ref(:,:)
REAL, ALLOCATABLE :: ref_tem(:,:)
REAL, ALLOCATABLE :: ref_volume(:,:,:) ! Array containing one
! full volume scan of
! reflectivity intepolated to
! 1-degree spaced polar angles
REAL, ALLOCATABLE :: vel(:,:)
REAL, ALLOCATABLE :: vel_tem(:,:)
REAL, ALLOCATABLE :: vel_volume(:,:,:) ! Array containing one
! full volume scan of
! reflectivity intepolated to
! 1-degree spaced polar angles
REAL :: int_razim(361)
REAL :: ref_eleva(nref_tilts_max)
REAL :: vel_eleva(nvel_tilts_max)
INTEGER :: maxgatevel
INTEGER :: maxgateref
INTEGER :: maxazim
INTEGER :: maxelev
REAL, ALLOCATABLE :: rngvvol(:,:)
REAL, ALLOCATABLE :: azmvvol(:,:)
REAL, ALLOCATABLE :: elvvvol(:,:)
REAL, ALLOCATABLE :: velvol(:,:,:)
REAL, ALLOCATABLE :: rngrvol(:,:)
REAL, ALLOCATABLE :: azmrvol(:,:)
REAL, ALLOCATABLE :: elvrvol(:,:)
REAL, ALLOCATABLE :: refvol(:,:,:)
REAL, ALLOCATABLE :: vnyquist(:)
REAL, ALLOCATABLE :: uarefvol(:,:,:)
REAL, ALLOCATABLE :: uavelvol(:,:,:)
REAL, ALLOCATABLE :: vorvol(:,:,:)
!
!-----------------------------------------------------------------------
!
! Misc. local variables:
!
!-----------------------------------------------------------------------
!
INTEGER :: nref_tilts
INTEGER :: nvel_tilts
INTEGER :: ifiles
CHARACTER(LEN=100) :: string
CHARACTER(LEN=100) :: stringtime
real :: radmax ,rad
INTEGER :: j_razim_start
INTEGER :: iazimuth
real :: razim_start
REAL :: elev_to_plot, xscale, yscale
INTEGER :: kref,kvel,ktilt
INTEGER :: kref_tilts,kvel_tilts
REAL :: refmax, refmin
REAL :: velmax, velmin
REAL :: eleva, deg2rad
REAL :: xradar, yradar, tem0,tem1,tem2,tem3,tem4
LOGICAL :: fexist
INTEGER :: zxplot_initialized = 0
REAL :: missing_value = -999.0
INTEGER :: iunit
INTEGER :: i,j,k
INTEGER :: jj
CHARACTER(LEN=100) :: psfilename
INTEGER :: ilen
REAL :: factorx,factory
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
! Beginning of executable code...
!
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
!
!-----------------------------------------------------------------------
! Input data file name and grid
!-----------------------------------------------------------------------
READ(5,input_data,ERR=100)
WRITE(6,'(a)') 'Namelist input_data was sucessfully read in.'
! WRITE(*,input_data)
! WRITE(*,*)
READ(5,radar_index)
WRITE(6,'(a)') 'Namelist radar_index was sucessfully read in.'
! WRITE(*,radar_index)
! WRITE(*,*)
READ(5,output)
WRITE(6,'(a)') 'Namelist output was sucessfully read in.'
!-----------------------------------------------------------------------
! Input map projection parameters
!-----------------------------------------------------------------------
sclfct = 1.0
READ (5,projection,ERR=100)
WRITE(6,'(a)') 'Namelist projection was sucessfully read in.'
! WRITE(*,projection)
! WRITE(*,*)
READ(5,plotting)
WRITE(6,'(a)') 'Namelist plotting sucessfully read in.'
READ(5,map_plot,ERR=100)
WRITE(6,'(a)') 'Namelist map_plot sucessfully read in.'
!
!
h_range_max = (h_range_max_e+h_range_max_w)/2.0
!-----------------------------------------------------------------------
! Read in reflectivity data and interpolate into INTEGER azimuth grid
!-----------------------------------------------------------------------
DO j=1,361
int_razim(j)=(j-1)*1.0
ENDDO
DO ifiles=1,numfiles
!
! read in data on a tilt
!
INQUIRE(FILE=trim(filepath)//trim(files(ifiles)),EXIST = fexist )
IF( .not. fexist) THEN
write(6,'(a,a,a)') 'Radar data file ', &
trim(files(ifiles)),' not found.'
cycle
ENDIF
WRITE(6,'(a,i4)') ' Radar format = ',radfmt
IF( radfmt == 1 ) THEN
print*,'To read file ', trim(files(ifiles))
OPEN(27,file=trim(filepath)//trim(files(ifiles)),form='unformatted', &
status='old')
READ(27) maxgatevel,maxgateref,maxazim,maxelev
READ(27) radar_elevation,radar_lat, radar_lon
READ(27) iyear,imon,iday,ihour,imin,isec
! write(*,*) maxgatevel,maxgateref,maxazim,maxelev
! write(*,*) radar_elevation,radar_lat, radar_lon
! write(*,*) iyear,imon,iday,ihour,imin,isec
ALLOCATE( rngrvol(maxgateref,maxelev) )
ALLOCATE( azmrvol(maxazim,maxelev) )
ALLOCATE( elvrvol(maxazim,maxelev) )
ALLOCATE( refvol(maxgateref,maxazim,maxelev) )
ALLOCATE( rngvvol(maxgatevel,maxelev) )
ALLOCATE( azmvvol(maxazim,maxelev) )
ALLOCATE( elvvvol(maxazim,maxelev) )
ALLOCATE( velvol(maxgatevel,maxazim,maxelev) )
READ(27) rngrvol
READ(27) azmrvol
READ(27) elvrvol
READ(27) refvol
READ(27) rngvvol
READ(27) azmvvol
READ(27) elvvvol
READ(27) velvol
close(27)
psfilename=trim(files(ifiles))
ilen=len(trim(psfilename))
iday=ICHAR(psfilename(ilen-12:ilen-12))-ICHAR("0")
iday=iday*10+(ICHAR(psfilename(ilen-11:ilen-11))-ICHAR("0"))
ihour=ICHAR(psfilename(ilen-9:ilen-9))-ICHAR("0")
ihour=ihour*10+(ICHAR(psfilename(ilen-8:ilen-8))-ICHAR("0"))
imin=ICHAR(psfilename(ilen-7:ilen-7))-ICHAR("0")
imin=imin*10+(ICHAR(psfilename(ilen-6:ilen-6))-ICHAR("0"))
ELSE IF (radfmt == 2) THEN ! radar emulator binary file
iunit=27
WRITE(6,'(a,a)') ' Opening: ',TRIM(filepath)//TRIM(files(ifiles))
OPEN(iunit,FILE=trim(filepath)//trim(files(ifiles)),STATUS='old', &
FORM='unformatted')
print *, ' opened file OK'
READ(iunit) runname
!print *, ' runname: ',runname
lfnkey=80
CALL gtlfnkey
(runname,lfnkey)
! WRITE(6,'(a,a)') ' Runname:',runname(1:lfnkey)
READ(iunit) radar_symbol
! WRITE(6,'(a,a)') ' Radar name:',radar_symbol
READ(iunit) maxgatevel,maxgateref,maxazim,maxelev
! WRITE(6,'(a,i6,a,i6,a,i6)') &
! ' Ngate : ',maxgatevel,' Nazim: ',maxazim,' maxelev: ',maxelev
READ(iunit) iyear,imon,iday,ihour,imin,isec
! WRITE(6,'(2(a,i2.2),a,i4.4,3(a,i2.2))') &
! ' Date: ',imon,'/',iday,'/',iyear,' Time: ',ihour,':',imin,':',isec
READ(iunit) radar_elevation,radar_lat,radar_lon
! WRITE(6,'(a,f9.2,a,f9.2,a,f9.2)') &
! ' Latitude:',radar_lat,' Longitude:',radar_lon, &
! ' Elevation:',radar_elevation
READ(iunit) vcp,beamwid,rmisval,rngfval,curtim
! WRITE(6,'(a,i5,a,f9.2,a,f10.1,a,f10.1)') &
! ' VCP:',vcp,' Beamwidth:',beamwid, &
! ' MissVal: ',rmisval,' RngfVal :',rngfval
!
CALL ctim2abss(iyear,imon,iday,ihour,imin,isec,abstsec)
READ(iunit) wrtuaref,wrtuavel,wrtvort, &
idummy,idummy,idummy,idummy,idummy
WRITE(6,'(a)') ' Allocating volume memory '
ALLOCATE(rngrvol(maxgateref,maxelev))
ALLOCATE(azmrvol(maxazim,maxelev))
ALLOCATE(elvrvol(maxazim,maxelev))
ALLOCATE(refvol(maxgateref,maxazim,maxelev))
IF(wrtuaref /= 0) ALLOCATE(uarefvol(maxgateref,maxazim,maxelev))
ALLOCATE(vnyquist(maxelev))
ALLOCATE(rngvvol(maxgatevel,maxelev))
ALLOCATE(azmvvol(maxazim,maxelev))
ALLOCATE(elvvvol(maxazim,maxelev))
ALLOCATE(velvol(maxgatevel,maxazim,maxelev))
IF(wrtuavel /= 0) ALLOCATE(uavelvol(maxgatevel,maxazim,maxelev))
IF(wrtvort /= 0) ALLOCATE(vorvol(maxgatevel,maxazim,maxelev))
WRITE(6,'(a)') ' Reading Reflectiving Arrays'
READ(iunit) rngrvol
READ(iunit) azmrvol
READ(iunit) elvrvol
READ(iunit) refvol
IF( wrtuaref /= 0 ) READ(iunit) uarefvol
WRITE(6,'(a)') ' Reading Velocity Arrays'
READ(iunit) vnyquist
READ(iunit) rngvvol
READ(iunit) azmvvol
READ(iunit) elvvvol
READ(iunit) velvol
IF( wrtuavel /= 0 ) READ(iunit) uavelvol
IF( wrtvort /= 0 ) READ(iunit) vorvol
!
WRITE(6,'(a)') ' Reading successfully completed'
CLOSE(iunit)
ELSE IF (radfmt == 4) THEN ! CASA-formatted NetCDF file
WRITE(6,'(a,i6,a)') ' Radar format number ',radfmt,' not yet implemented'
STOP
ELSE
WRITE(6,'(a,i6,a)') ' Radar format number ',radfmt,' not supported.'
STOP
END IF
IF ( mapproj == 0 ) THEN
trulat1 = radar_lat
trulat2 = radar_lat
trulon = radar_lon
END IF
kref_tilts=0
DO k=1,maxelev
if(elvrvol(1,k) >= -100 ) THEN
kref_tilts=kref_tilts+1
ref_eleva(kref_tilts) = elvrvol(1,k)
ENDIF
ENDDO
nref_tilts = kref_tilts
! WRITE(6,'(a)') ' Radar elevations read: '
! write(6,'(10f8.1)') (ref_eleva(i),i=1,nref_tilts)
nrefgate=0
DO j=1,maxgateref
if(rngrvol(j,1) >= -700 ) THEN
nrefgate=nrefgate+1
ENDIF
ENDDO
rfstgat=rngrvol(1,1)
refgatsp=rngrvol(3,1)-rngrvol(2,1)
ALLOCATE( ref_volume(nrefgate,361,nref_tilts) )
ALLOCATE( razim(maxazim) )
ALLOCATE( razim_tem(maxazim) )
print *, ' Begin Reflectivity Interpolation'
DO ktilt=1,nref_tilts
nrefray=0
DO j=1,maxazim
if(azmrvol(j,ktilt) >= -100 ) THEN
nrefray=nrefray+1
razim(nrefray)=azmrvol(j,ktilt)
ENDIF
ENDDO
IF( nrefray == maxazim) nrefray=nrefray-1
ALLOCATE( ref(nrefgate,nrefray+1) )
ALLOCATE( ref_tem(nrefgate,nrefray+1) )
j_razim_start = 1
razim_start = razim(1)
DO j=2,nrefray
if( razim(j).lt.razim_start ) then
j_razim_start = j
exit
endif
ENDDO
DO j=1,nrefray
DO i=1,nrefgate
ref(i,j)= refvol(i,j,ktilt)
ENDDO
ENDDO
DO j=1,nrefray
DO i=1,nrefgate
IF(abs(ref(i,j)) >= 800.0 ) ref(i,j)= missing_value
ENDDO
ENDDO
!
! resort the data
!
DO j=j_razim_start,nrefray
razim_tem(j-j_razim_start+1) = razim(j)
ENDDO
DO j=1,j_razim_start-1
razim_tem(nrefray-j_razim_start+1+j) = razim(j)
ENDDO
razim_tem(nrefray+1) = razim_tem(1) + 360.0
DO i=1,nrefgate
DO j=j_razim_start,nrefray
ref_tem(i,j-j_razim_start+1) = ref(i,j)
ENDDO
DO j=1,j_razim_start-1
ref_tem(i,nrefray-j_razim_start+1+j)=ref(i,j)
ENDDO
ref_tem(i,nrefray+1)=ref_tem(i,1)
ENDDO
razim = razim_tem
ref = ref_tem
!
! interpolate into regular azimuth value and save in 3d space
!
DO i=1,nrefgate
DO j=1,361
tem0 = int_razim(j)
IF( tem0 == 0.0 .and. tem0.lt. razim(1) ) then
tem0 = 360.0
ENDIF
DO jj=max(1,nint(tem0-20.0)),nrefray
tem1 = razim(jj)
tem2 = razim(jj+1)
IF((tem0.ge.tem1).and.(tem0.lt.tem2)) then
tem3 = ref(i,jj)
tem4 = ref(i,jj+1)
IF( abs(tem3-missing_value) <=0.001 .and. &
abs(tem4-missing_value) <=0.001 ) THEN
ref_volume(i,j,ktilt) = &
tem3+(tem4-tem3)*(tem0-tem1)/(tem2-tem1)
ELSE
IF( abs(tem0-tem1).lt.abs(tem0-tem2)) THEN
ref_volume(i,j,ktilt) = tem3
ELSE
ref_volume(i,j,ktilt) = tem4
ENDIF
ENDIF
GOTO 300
ENDIF
ENDDO
ref_volume(i,j,ktilt) = missing_value
300 CONTINUE
ENDDO
ENDDO
DEALLOCATE( ref )
DEALLOCATE( ref_tem )
DO i=1,nrefgate
ref_volume(i,361,ktilt) = ref_volume(i,1,ktilt)
ENDDO
ENDDO ! ktilt
DEALLOCATE( razim )
DEALLOCATE( razim_tem )
! output data for further retrival
IF( iout_ref_regular == 1 ) THEN
write(ref_regular_file,'(a,a,a,a)') trim(filepath),trim(files(ifiles)), &
'.ref_regular', trim(filetrailer)
CALL OUTPUT_REGULAR(ref_regular_file,nref_tilts_max,nref_tilts, &
radar_symbol,radar_lat,radar_lon,radar_elevation, &
nrefgate,361,rfstgat,refgatsp,int_razim,ref_eleva, &
ref_volume )
ENDIF
!-----------------------------------------------------------------------
! interpolate radial velocity into integer azimuth grid
!-----------------------------------------------------------------------
print *, ' Begin velocity interpolation'
kvel_tilts=0
DO k=1,maxelev
if(elvvvol(1,k) >= -100 ) THEN
kvel_tilts=kvel_tilts+1
vel_eleva(kvel_tilts) = elvvvol(1,k)
ENDIF
ENDDO
nvel_tilts = kvel_tilts
nvelgate=0
DO j=1,maxgatevel
if(rngvvol(j,1) >= -700 ) THEN
nvelgate=nvelgate+1
ENDIF
ENDDO
vfstgat=rngvvol(1,1)
velgatsp=rngvvol(3,1)-rngvvol(2,1)
ALLOCATE( vel_volume(nvelgate,361,nvel_tilts) )
ALLOCATE( razim(maxazim) )
ALLOCATE( razim_tem(maxazim) )
DO ktilt=1,nvel_tilts
nrefray=0
DO j=1,maxazim
if(azmvvol(j,ktilt) >= -100 ) THEN
nrefray=nrefray+1
razim(nrefray)=azmvvol(j,ktilt)
ENDIF
ENDDO
nvelray=nrefray
if( nrefray == maxazim ) nvelray=nrefray-1
ALLOCATE( vel(nvelgate,nvelray+1) )
ALLOCATE( vel_tem(nvelgate,nvelray+1) )
j_razim_start = 1
razim_start = razim(1)
DO j=2,nvelray
if( razim(j).lt.razim_start ) then
j_razim_start = j
exit
endif
ENDDO
DO j=1,nvelray
DO i=1,nvelgate
vel(i,j)= velvol(i,j,ktilt)
ENDDO
ENDDO
! rearrange missing data value
DO j=1,nvelray
DO i=1,nvelgate
IF(abs(vel(i,j)) >= 700.0 ) vel(i,j)= missing_value
ENDDO
ENDDO
!
! resort the data
!
DO j=j_razim_start,nvelray
razim_tem(j-j_razim_start+1) = razim(j)
ENDDO
DO j=1,j_razim_start-1
razim_tem(nvelray-j_razim_start+1+j) = razim(j)
ENDDO
razim_tem(nvelray+1) = razim_tem(1) + 360.0
DO i=1,nvelgate
DO j=j_razim_start,nvelray
vel_tem(i,j-j_razim_start+1) = vel(i,j)
ENDDO
DO j=1,j_razim_start-1
vel_tem(i,nvelray-j_razim_start+1+j)=vel(i,j)
ENDDO
vel_tem(i,nvelray+1)=vel_tem(i,1)
ENDDO
razim = razim_tem
vel = vel_tem
!
! interpolate into regular azimuth value and save in 3d space
!
DO i=1,nvelgate
DO j=1,361
tem0 = int_razim(j)
IF( tem0 == 0.0 .and. tem0.lt. razim(1) ) then
tem0 = 360.0
ENDIF
DO jj=max(1,nint(tem0-20.0)),nvelray
tem1 = razim(jj)
tem2 = razim(jj+1)
IF((tem0.ge.tem1).and.(tem0.lt.tem2)) then
tem3 = vel(i,jj)
tem4 = vel(i,jj+1)
IF( abs(tem3-missing_value) <=0.001 .and. &
abs(tem4-missing_value) <=0.001 ) THEN
vel_volume(i,j,ktilt) = &
tem3+(tem4-tem3)*(tem0-tem1)/(tem2-tem1)
ELSE
IF( abs(tem0-tem1).lt.abs(tem0-tem2)) THEN
vel_volume(i,j,ktilt) = tem3
ELSE
vel_volume(i,j,ktilt) = tem4
ENDIF
ENDIF
GOTO 200
ENDIF
ENDDO
vel_volume(i,j,ktilt) = missing_value
200 CONTINUE
ENDDO
ENDDO
DEALLOCATE( vel )
DEALLOCATE( vel_tem )
DO i=1,nvelgate
vel_volume(i,361,ktilt) = vel_volume(i,1,ktilt)
ENDDO
ENDDO ! ktilt
DEALLOCATE( razim )
DEALLOCATE( razim_tem )
! output data for further retrival
IF( iout_vel_regular == 1 ) THEN
write(vel_regular_file,'(a,a,a,a)') trim(filepath),trim(files(ifiles)), &
'.vel_regular', trim(filetrailer)
CALL OUTPUT_REGULAR(vel_regular_file,nvel_tilts_max,nvel_tilts, &
radar_symbol,radar_lat,radar_lon,radar_elevation, &
nvelgate,361,vfstgat,velgatsp,int_razim,vel_eleva, &
vel_volume )
ENDIF
!-----------------------------------------------------------------------
! plot
!-----------------------------------------------------------------------
print *, ' Begin plotting'
factorx=0.3
factory=0.3
if(h_range_max_w + h_range_max_e > h_range_max_s + h_range_max_n ) then
h_range_max=(h_range_max_e+h_range_max_w)/2
factory=factorx*(h_range_max_s + h_range_max_n)/(h_range_max_w + h_range_max_e)
else
factorx=factory*(h_range_max_w + h_range_max_e)/(h_range_max_s + h_range_max_n)
h_range_max=(h_range_max_n+h_range_max_s)/2
endif
!-----------------------------------------------------------------------
! plot reflectivity on a tilt
!-----------------------------------------------------------------------
deg2rad = atan(1.0)/45.0
IF( n_elevation_ref > 0 .or. n_azimuth_ref > 0 .or. &
n_elevation_vel > 0 .or. n_azimuth_vel > 0 ) THEN
IF( zxplot_initialized == 0 ) then
CALL xpsfn('zxout'//'.ps', 5)
! PS output will be called zxout.ps
CALL xpaprlnth( 1.0 )
CALL xdevic ! beigin graphic
zxplot_initialized = 1
ENDIF
ENDIF
print *, ' nref_tilts, ref_elva: ',nref_tilts,ref_eleva
IF(n_elevation_ref > 0) THEN
IF( nref_tilts <= 0 ) THEN
WRITE(*,*) ' There are no data to be plotted ???'
STOP 119
ENDIF
ALLOCATE( x(nrefgate,361))
ALLOCATE( y(nrefgate,361))
DO ktilt = 1,n_elevation_ref
elev_to_plot = elevations_ref(ktilt)
kref = 0
DO k = 1,nref_tilts
IF( abs( ref_eleva(k)-elev_to_plot ) &
< elev_error_range ) THEN
kref = k
EXIT
ENDIF
ENDDO
IF(kref == 0 ) CYCLE ! no elevation angle
! within error range was found
refmin = 1000.0
refmax =-1000.0
velmin = 1000.0
refmax =-1000.0
eleva=ref_eleva(kref)
DO j=1,361
DO i=1,nrefgate
rad = 0.001 * rngrvol(i,kref)
x(i,j)= rad * cos( (90-int_razim(j))*deg2rad )
y(i,j)= rad * sin( (90-int_razim(j))*deg2rad )
tem0 = ref_volume(i,j,kref)
IF(tem0 /= missing_value ) then
refmax=max(tem0,refmax)
refmin=min(tem0,refmin)
ENDIF
END DO
END DO
! print*,'refmin, refmax = ', refmin, refmax
!
!-----------------------------------------------------------------------
! Define plotting space.
!
CALL xpspac(0.5-factorx,0.5+factorx,0.5-factory,0.5+factory)
CALL xmap(-h_range_max_w, h_range_max_e, -h_range_max_s, h_range_max_n)
CALL xlabmask(1)
IF( coltab_ref == -1 ) CALL xstctfn(coltabfn_ref)
CALL xsetclrs(coltab_ref) ! Use internally defined color map number 6
CALL xcolor(1) ! Set color to black
CALL xcfont(3)
!
call xaxfmt('(I4)')
!
! Set tentative color interval and limit on the min/max number of contours
!
cl(1)=0.0
cl(2)=5.0 ! Set tentative contour interval
CALL xnctrs(2,100) ! Set lower and upper limits
! on the number of contours allowed
! It's used by both xcolfil and xclimit
CALL xnctrs(5,20) ! Set lower and upper limits
mode=2
call xctrbadv(1) ! Turn on bad value checking
! in contouring routine
call xbadval( missing_value ) ! Specify bad value flag
!
! plot a color filled contour field for positive values
!
CALL xctrclr(ibgncol_ref,iendcol_ref) ! Use colors between 10 and 39 inclusive
CALL xctrlim(5.0, 75.0) ! Only contours above zero are plotted
ALLOCATE( iw(nrefgate,361) )
ALLOCATE( xw(8*nrefgate) )
ALLOCATE( yw(8*nrefgate) )
CALL xwindw(-h_range_max_w, h_range_max_e, -h_range_max_s, h_range_max_n)
Print*,'To plot reflectivity field at elevation angle = ', eleva
CALL xcolfil(ref_volume(1,1,kref),x,y,iw,xw,yw, &
nrefgate,nrefgate,361, cl, ncl,mode)
if ( overlayrefcontour == 1 ) THEN
! call XCLFMT('I2')
call XCLTYP(0)
CALL xctrclr(14,15) ! Use colors between 10 and 39 inclusive
mode=3
ncl=1
cl(1)=15.0
CALL XCONTA(ref_volume(1,1,kref),x,y,IW,nrefgate,nrefgate,361, CL, NCL,MODE)
cl(1)=30.0
CALL XCONTA(ref_volume(1,1,kref),x,y,IW,nrefgate,nrefgate,361, CL, NCL,MODE)
cl(1)=45.0
CALL XCONTA(ref_volume(1,1,kref),x,y,IW,nrefgate,nrefgate,361, CL, NCL,MODE)
cl(1)=55.0
CALL XCONTA(ref_volume(1,1,kref),x,y,IW,nrefgate,nrefgate,361, CL, NCL,MODE)
endif
CALL xwdwof
CALL xaxsca(-h_range_max_w, h_range_max_e, 0.1*h_range_max, &
-h_range_max_s, h_range_max_n, 0.1*h_range_max)
DEALLOCATE( iw )
DEALLOCATE( xw )
DEALLOCATE( yw )
WRITE(string,'(a,f5.2)') 'Reflectivity at Elevation Angle=',eleva
CALL xchsiz(0.025*2*h_range_max ) ! Set character size
CALL xcharc((-h_range_max_w+h_range_max_e)*0.5, &
h_range_max_n+2*h_range_max*0.04, trim(string) )
WRITE(stringtime,'(a,i2.2,a,i2.2,a,i4,a,i2.2,a,i2.2)') 'Scan Time:', &
imon,'/',iday,'/',iyear, ' ',ihour,':',imin
CALL xchsiz(0.025*2*h_range_max ) ! Set character size
CALL xcharc((-h_range_max_w+h_range_max_e)*0.5, &
h_range_max_n+2*h_range_max*0.01, trim(stringtime) )
CALL xchsiz(0.02* 2*h_range_max ) ! Set character size
CALL xcpalet(1) ! Plot horizontal color palette
xscale = h_range_max
yscale = h_range_max
xorig = -(h_range_max_e+h_range_max_w)/2.0
yorig = -(h_range_max_s+h_range_max_n)/2.0
xradar = 0.0
yradar = 0.0
CALL xthick(3)
CALL xpenup( xradar-xscale*0.02, yradar )
CALL xpendn( xradar+xscale*0.02, yradar )
CALL xpenup( xradar, yradar-yscale*0.02 )
CALL xpendn( xradar, yradar+yscale*0.02 )
call xthick(1)
CALL xchsiz(0.02* 2*h_range_max ) ! Set character size
CALL xcharc(xradar, yradar-0.05*yscale, radar_symbol)
xscale = h_range_max_e+h_range_max_w
yscale = h_range_max_s+h_range_max_n
CALL XSTPJGRD(mapproj,trulat1,trulat2,trulon,ctrlat,ctrlon, &
xscale,yscale,xorig,yorig)
CALL xwindw(-h_range_max_w, h_range_max_e, -h_range_max_s, h_range_max_n)
IF( ovrmap /= 0 ) THEN
DO i=1,MIN(mxmapfile,nmapfile)
lmapfile=len_trim(mapfile(i))
WRITE(6,'(1x,a,a)') 'Input was ',trim(mapfile(i))
INQUIRE(FILE=trim(mapfile(i)), EXIST = fexist )
IF( .NOT.fexist) THEN
WRITE(6,'(a)') 'Map file '//mapfile(i)(1:lmapfile)// &
' not found. Corresponding map no plotted.'
ELSE
CALL xdrawmap(10,mapfile(i)(1:lmapfile),latgrid,longrid)
END IF
END DO
END IF
CALL xwdwof
CALL xframe
ENDDO ! ktilt
DEALLOCATE( x )
DEALLOCATE( y )
ENDIF ! n_elevation_ref>0
!
!-----------------------------------------------------------------------
! To plot vertical cross-sections of reflectivity
!-----------------------------------------------------------------------
!
IF(n_azimuth_ref > 0) THEN
IF( nref_tilts <= 0 ) THEN
WRITE(*,*) ' There are no data to be plotted ???'
STOP 119
ENDIF
ALLOCATE( zp(nrefgate,nref_tilts))
ALLOCATE( zm(nrefgate,nref_tilts))
ALLOCATE( radialp(nrefgate,nref_tilts))
ALLOCATE( radialm(nrefgate,nref_tilts))
DO k=1,nref_tilts
DO i=1,nrefgate
rad = 0.001*(rfstgat + (i-1)*refgatsp)
radialp(i,k)=rad*cos( ref_eleva(k)*deg2rad )
zp(i,k) =rad*sin( ref_eleva(k)*deg2rad )
END DO
DO i=1,nrefgate
rad = 0.001*(rfstgat + (i-1)*refgatsp)
radialm(i,k)=-rad*cos( ref_eleva(k)*deg2rad )
zm(i,k) = rad*sin( ref_eleva(k)*deg2rad )
END DO
END DO
DO iazimuth= 1,n_azimuth_ref
Print*,'To plot vertical cross-section at azimuth angle = ', &
azimuth_ref(iazimuth)
CALL xpspac(0.1,0.9,0.2,0.6)
z_range_max = 12.0
! CALL xmap(-h_range_max, h_range_max, 0.0, z_range_max)
CALL xmap(0, h_range_max, 0.0, z_range_max)
CALL xlabmask(1)
IF( coltab_ref == -1 ) CALL xstctfn(coltabfn_ref)
CALL xsetclrs(coltab_ref) ! Use internally defined color map number 6
CALL xcolor(1) ! Set color to black
CALL xcfont(3)
!
!
! Set tentative color interval and limit on the min/max number of contours
!
cl(1)=0.0
cl(2)=5.0 ! Set tentative contour interval
CALL xnctrs(2,100) ! Set lower and upper limits
! on the number of contours allowed
! It's used by both xcolfil and xclimit
CALL xnctrs(5,20) ! Set lower and upper limits
mode=2
CALL xctrbadv(1) ! Turn on bad value checking
! in contouring routine
CALL xbadval( missing_value ) ! Specify bad value flag
!
! plot a color filled contour field for positive values
!
CALL xctrclr(ibgncol_ref,iendcol_ref) ! Use colors between 10 and 39 inclusive
CALL xctrlim(5.0, 75.0) ! Only contours above zero are plotted
CALL xwindw(-h_range_max, h_range_max, 0.0, z_range_max)
ALLOCATE( iw(nrefgate,nref_tilts) )
ALLOCATE( xw(8*nrefgate) )
ALLOCATE( yw(8*nrefgate) )
ALLOCATE( ref_rz(nrefgate,nref_tilts) )
jj = min(361,max(1,nint(azimuth_ref(iazimuth))+1 ))
DO k=1,nref_tilts
DO i=1,nrefgate
ref_rz(i,k) = ref_volume(i,jj,k)
ENDDO
ENDDO
CALL xcolfil(ref_rz,radialp,zp,iw,xw,yw,nrefgate,nrefgate, &
nref_tilts,cl,ncl,mode)
! tem0 = azimuth_ref(iazimuth)+180.0
101 CONTINUE
! IF( tem0 > 360.0 ) then
! tem0 = tem0 - 360.0
! GOTO 101
! ENDIF
!
! jj = nint(tem0)+1
! DO k=1,nref_tilts
! DO i=1,nrefgate
! ref_rz(i,k) = ref_volume(i,jj,k)
! ENDDO
! ENDDO
! CALL xcolfil(ref_rz,radialm,zm,iw,xw,yw,nrefgate,nrefgate, &
! nref_tilts,cl,ncl,mode)
DEALLOCATE( iw )
DEALLOCATE( xw )
DEALLOCATE( yw )
DEALLOCATE( ref_rz )
CALL xwdwof
! CALL xaxsca(-h_range_max, h_range_max, 0.1*h_range_max, &
! 0.0, z_range_max, 0.5 )
CALL xaxsca(0, h_range_max, 10.0,0.0, z_range_max, 1.0)
WRITE(string,'(a,f6.2)') &
' Reflectivity at Azimuth Angle=',int_razim(jj)
CALL xchsiz(0.06*z_range_max ) ! Set character size
CALL xcharc(h_range_max/2.0, z_range_max*1.10, trim(string) )
WRITE(stringtime,'(a,i2.2,a,i2.2,a,i4,a,i2.2,a,i2.2)') 'First Scan Time:',&
imon,'/',iday,'/',iyear, &
' ',ihour,':',imin
CALL xchsiz(0.05*z_range_max ) ! Set character size
CALL xcharc(h_range_max/2.0, z_range_max*1.04, trim(stringtime) )
CALL xchsiz(0.04*z_range_max ) ! Set character size
CALL xcpalet(1) ! Plot horizontal color palette
CALL xframe
ENDDO ! iazimuth
DEALLOCATE(zp)
DEALLOCATE(zm)
DEALLOCATE(radialp)
DEALLOCATE(radialm)
ENDIF ! n_azimuth_ref
!-----------------------------------------------------------------------
! plot radial velocity on a tilt
!-----------------------------------------------------------------------
IF(n_elevation_vel > 0) THEN
IF( nvel_tilts <= 0 ) THEN
WRITE(*,*) ' There are no radial velocity data to be plotted ???'
STOP 119
ENDIF
ALLOCATE( x(nvelgate,361))
ALLOCATE( y(nvelgate,361))
DO ktilt = 1,n_elevation_vel
elev_to_plot = elevations_vel(ktilt)
kvel = 0
DO k = 1,nvel_tilts
IF( abs( vel_eleva(k)-elev_to_plot ) &
< elev_error_range ) THEN
kvel = k
EXIT
ENDIF
ENDDO
IF(kvel == 0 ) CYCLE ! no elevation angle
! within error range was found
eleva=vel_eleva(kvel)
velmin = 1000.0
velmax =-1000.0
DO j=1,361
DO i=1,nvelgate
rad = rngvvol(i,kvel)
rad = 0.001 * rad
x(i,j)= rad * cos( (90-int_razim(j))*deg2rad )
y(i,j)= rad * sin( (90-int_razim(j))*deg2rad )
tem0 = vel_volume(i,j,kvel)
IF(tem0 /= missing_value ) then
velmax=max(tem0,velmax)
velmin=min(tem0,velmin)
ENDIF
END DO
END DO
print*,'velmin, velmax = ', velmin, velmax
!
! Define plotting space.
!
CALL xpspac(0.5-factorx,0.5+factorx,0.5-factory,0.5+factory)
CALL xmap(-h_range_max_w, h_range_max_e, -h_range_max_s, h_range_max_n)
CALL xlabmask(1)
IF( coltab_vel == -1 ) CALL xstctfn(coltabfn_vel)
CALL xsetclrs(coltab_vel) ! Use internally defined color map number 6
CALL xcolor(1) ! Set color to black
CALL xcfont(3)
!
call xaxfmt('(I4)')
!
! Set tentative color interval and limit on the min/max number of contours
!
cl(1)=0.0
cl(2)=4.0 ! Set tentative contour interval
CALL xnctrs(2,100) ! Set lower and upper limits
! on the number of contours allowed
! It's used by both xcolfil and xclimit
CALL xnctrs(2,100) ! Set lower and upper limits
mode=2
call xctrbadv(1) ! Turn on bad value checking
! in contouring routine
call xbadval( missing_value ) ! Specify bad value flag
!
! plot a color filled contour field for positive values
!
CALL xctrclr(ibgncol_vel,iendcol_vel) ! Use colors between 10 and 39 inclusive
CALL xctrlim(-40.0, 40.0) ! Only contours above zero are plotted
! CALL xctrlim(-9999.0, -9999.0) ! Only contours above zero are plotted
ALLOCATE( iw(nvelgate,361) )
ALLOCATE( xw(8*nvelgate) )
ALLOCATE( yw(8*nvelgate) )
CALL xwindw(-h_range_max_w, h_range_max_e, -h_range_max_s, h_range_max_n)
Print*,'To plot radial velocity field at elevation angle =', eleva
CALL xcolfil(vel_volume(1,1,kvel),x,y,iw,xw,yw, &
nvelgate,nvelgate,361, cl, ncl,mode)
DEALLOCATE( iw )
DEALLOCATE( xw )
DEALLOCATE( yw )
CALL xwdwof
CALL xaxsca(-h_range_max_w, h_range_max_e, 0.1*h_range_max, &
-h_range_max_s, h_range_max_n, 0.1*h_range_max)
WRITE(string,'(a,f5.2)') 'Radial velocity at Elevation Angle=',eleva
CALL xchsiz(0.025*2*h_range_max ) ! Set character size
CALL xcharc((-h_range_max_w+h_range_max_e)*0.5, &
h_range_max_n+2*h_range_max*0.04, trim(string) )
WRITE(stringtime,'(a,i2.2,a,i2.2,a,i4.2,a,i2.2,a,i2.2)') 'Scan Time:', &
imon,'/',iday,'/',iyear, &
' ',ihour,':',imin
CALL xchsiz(0.025*2*h_range_max ) ! Set character size
CALL xcharc((-h_range_max_w+h_range_max_e)*0.5, &
h_range_max_n+2*h_range_max*0.01, trim(stringtime) )
CALL xchsiz(0.02* 2*h_range_max ) ! Set character size
CALL xcpalet(1) ! Plot horizontal color palette
xscale = h_range_max
yscale = h_range_max
xorig = -(h_range_max_e+h_range_max_w)/2.0
yorig = -(h_range_max_s+h_range_max_n)/2.0
xradar = 0.0
yradar = 0.0
CALL xthick(3)
CALL xpenup( xradar-xscale*0.02, yradar )
CALL xpendn( xradar+xscale*0.02, yradar )
CALL xpenup( xradar, yradar-yscale*0.02 )
CALL xpendn( xradar, yradar+yscale*0.02 )
call xthick(1)
CALL xchsiz(0.02* 2*h_range_max ) ! Set character size
CALL xcharc(xradar, yradar-0.05*yscale, radar_symbol)
xscale = h_range_max_e+h_range_max_w
yscale = h_range_max_s+h_range_max_n
CALL XSTPJGRD(mapproj,trulat1,trulat2,trulon,ctrlat,ctrlon, &
xscale,yscale,xorig,yorig)
CALL xwindw(-h_range_max_w, h_range_max_e, -h_range_max_s, h_range_max_n)
IF( ovrmap /= 0 ) THEN
DO i=1,MIN(mxmapfile,nmapfile)
lmapfile=len_trim(mapfile(i))
WRITE(6,'(1x,a,a)') 'Input was ',trim(mapfile(i))
INQUIRE(FILE=trim(mapfile(i)), EXIST = fexist )
IF( .NOT.fexist) THEN
WRITE(6,'(a)') 'Map file '//mapfile(i)(1:lmapfile)// &
' not found. Corresponding map no plotted.'
ELSE
CALL xdrawmap(10,mapfile(i)(1:lmapfile),latgrid,longrid)
END IF
END DO
END IF
CALL xwdwof
CALL xframe
ENDDO ! ktilt
DEALLOCATE( x )
DEALLOCATE( y )
ENDIF ! n_elevation_vel>0
!
!-----------------------------------------------------------------------
! To plot vertical cross-sections of radial velocity
!-----------------------------------------------------------------------
!
IF(n_azimuth_vel > 0) THEN
IF( nvel_tilts <= 0 ) THEN
WRITE(*,*) ' There are no radial velocity data to be plotted ???'
STOP 119
ENDIF
ALLOCATE( zp(nvelgate,nvel_tilts))
ALLOCATE( zm(nvelgate,nvel_tilts))
ALLOCATE( radialp(nvelgate,nvel_tilts))
ALLOCATE( radialm(nvelgate,nvel_tilts))
DO k=1,nvel_tilts
DO i=1,nvelgate
rad = 0.001*(vfstgat + (i-1)*velgatsp)
radialp(i,k)=rad*cos( vel_eleva(k)*deg2rad )
zp(i,k) =rad*sin( vel_eleva(k)*deg2rad )
END DO
DO i=1,nvelgate
rad = 0.001*(vfstgat + (i-1)*velgatsp)
radialm(i,k)=-rad*cos( vel_eleva(k)*deg2rad )
zm(i,k) = rad*sin( vel_eleva(k)*deg2rad )
END DO
END DO
DO iazimuth= 1,n_azimuth_vel
Print*,'To plot vertical cross-section at azimuth angle = ', &
azimuth_vel(iazimuth)
CALL xpspac(0.1,0.9,0.2,0.6)
z_range_max = 12.0
! CALL xmap(-h_range_max, h_range_max, 0.0, z_range_max)
CALL xmap(0.0, h_range_max, 0.0, z_range_max)
CALL xlabmask(1)
IF( coltab_vel == -1 ) CALL xstctfn(coltabfn_vel)
CALL xsetclrs(coltab_vel) ! Use internally defined color map number 6
CALL xcolor(1) ! Set color to black
CALL xcfont(3)
!
!
! Set tentative color interval and limit on the min/max number of contours
!
cl(1)=0.0
cl(2)=4.0 ! Set tentative contour interval
CALL xnctrs(2,100) ! Set lower and upper limits
! on the number of contours allowed
! It's used by both xcolfil and xclimit
CALL xnctrs(5,20) ! Set lower and upper limits
mode=2
CALL xctrbadv(1) ! Turn on bad value checking
! in contouring routine
CALL xbadval( missing_value ) ! Specify bad value flag
!
! plot a color filled contour field for positive values
!
CALL xctrclr(ibgncol_vel,iendcol_vel) ! Use colors between 10 and 39 inclusive
CALL xctrlim(-40.0, 40.0) ! Only contours above zero are plotted
CALL xwindw(-h_range_max, h_range_max, 0.0, z_range_max)
ALLOCATE( iw(nvelgate,nvel_tilts) )
ALLOCATE( xw(8*nvelgate) )
ALLOCATE( yw(8*nvelgate) )
ALLOCATE( vel_rz(nvelgate,nvel_tilts) )
jj = min(361,max(1,nint(azimuth_vel(iazimuth))+1 ))
DO k=1,nvel_tilts
DO i=1,nvelgate
vel_rz(i,k) = vel_volume(i,jj,k)
ENDDO
ENDDO
CALL xcolfil(vel_rz,radialp,zp,iw,xw,yw,nvelgate,nvelgate, &
nvel_tilts,cl,ncl,mode)
! tem0 = azimuth_vel(iazimuth)+180.0
!
!102 CONTINUE
! IF( tem0 > 360.0 ) then
! tem0 = tem0 - 360.0
! GOTO 102
! ENDIF
!
! jj = nint(tem0)+1
! DO k=1,nvel_tilts
! DO i=1,nvelgate
! vel_rz(i,k) = vel_volume(i,jj,k)
! ENDDO
! ENDDO
! CALL xcolfil(vel_rz,radialm,zm,iw,xw,yw,nvelgate,nvelgate, &
! nvel_tilts,cl,ncl,mode)
DEALLOCATE( iw )
DEALLOCATE( xw )
DEALLOCATE( yw )
DEALLOCATE( vel_rz )
CALL xwdwof
! CALL xaxsca(-h_range_max, h_range_max, 0.1*h_range_max, &
! 0.0, z_range_max, 0.5 )
CALL xaxsca(0.0, h_range_max, 10.0,0.0, z_range_max, 1.0)
WRITE(string,'(a,f6.2)') &
' Radial velocity at Azimuth Angle=',int_razim(jj)
CALL xchsiz(0.06*z_range_max ) ! Set character size
CALL xcharc(h_range_max/2.0, z_range_max*1.10, trim(string) )
WRITE(stringtime,'(a,i2.2,a,i2.2,a,i4,a,i2.2,a,i2.2)') 'First Scan Time:',&
imon,'/',iday,'/',iyear, &
' ',ihour,':',imin
CALL xchsiz(0.05*z_range_max ) ! Set character size
CALL xcharc(h_range_max/2.0, z_range_max*1.04, trim(stringtime) )
CALL xchsiz(0.04*z_range_max ) ! Set character size
CALL xcpalet(1) ! Plot horizontal color palette
CALL xframe
ENDDO ! iazimuth
DEALLOCATE(zp)
DEALLOCATE(zm)
DEALLOCATE(radialp)
DEALLOCATE(radialm)
ENDIF ! n_azimuth_vel
!
!-----------------------------------------------------------------------
! END of PLOT
!-----------------------------------------------------------------------
!
DEALLOCATE( rngvvol)
DEALLOCATE( azmvvol)
DEALLOCATE( elvvvol)
DEALLOCATE( velvol)
DEALLOCATE( rngrvol)
DEALLOCATE( azmrvol)
DEALLOCATE( elvrvol)
DEALLOCATE( refvol)
DEALLOCATE( ref_volume )
DEALLOCATE( vel_volume )
ENDDO ! ifiles
IF( zxplot_initialized == 1 ) CALL xgrend ! End graphics
STOP 999
100 Continue
Write(6,'(a)') 'Error reading namelist input file.'
Write(6,'(a)') 'Program aborted.'
STOP 123
END PROGRAM PLT_RADAR_TILT
!
!##################################################################
!##################################################################
!###### ######
!###### Advanced Regional Prediction System (ARPS) ######
!###### ######
!###### Developed by ######
!###### Center for Analysis and Prediction of Storms ######
!###### University of Oklahoma ######
!###### ######
!##################################################################
!##################################################################
!
!
SUBROUTINE OUTPUT_REGULAR(regular_file, ntilts_max,ntilts, & 2
radar_symbol,radar_lat,radar_lon,radar_elevation, &
ngate,nray,fstgat,gatsp,int_razim,eleva, &
volume )
!
!-----------------------------------------------------------------------
!
! PURPOSE:
! Out put 3d regular polar grid data for further retrival.
!
!-----------------------------------------------------------------------
!
!
! AUTHOR: Ming Hu
! 03/19/02.
!
! MODIFICATION HISTORY:
!
! INPUT:
!
! regular_file file name that data will be saved into
! ntilts_max maximum number of tilts
! ntilts actual number of tilts in this observation data
! radar_symbol radar Name
! radar_lat radar location: latitude
! radar_lon radar location: longitude
! radar_elevation radar location: elavation
!
! ngate number of gates of observation data per radial
! nray number of radials of data( here = 361 )
! fstgat distance from radar of first gate of data(meters)
! gatsp gate spacing of data (meters)
!
! volume(ngate,nray,ntilts) Array containing one full volume scan of
! observation(reflectivity or radial velocity
! which has been intepolated to 1-degree
! regular spaced polar angles
! int_razim(nray) azimuth angle for each radial (degrees)
! eleva(ntilts_max) elevation angle for each tilt
!
!
!
! OUTPUT:
!
!-----------------------------------------------------------------------
!
! Variable Declarations.
!
!-----------------------------------------------------------------------
!
IMPLICIT NONE
!-----------------------------------------------------------------------
!
! INPUT RADAR PARAMETER AND OBSERVATION
!
!-----------------------------------------------------------------------
!
CHARACTER (LEN=132) :: regular_file
INTEGER :: ntilts_max
INTEGER :: ntilts
CHARACTER(LEN=4) :: radar_symbol
REAL :: radar_lat
REAL :: radar_lon
REAL :: radar_elevation
INTEGER :: ngate
INTEGER :: nray
REAL :: fstgat
REAL :: gatsp
REAL :: volume(ngate,nray,ntilts) ! Array containing one
! full volume scan of
! reflectivity intepolated to
! 1-degree spaced polar angles
REAL :: int_razim(nray)
REAL :: eleva(ntilts_max)
!
!-----------------------------------------------------------------------
!
! WRITE(*,*) regular_file
! WRITE(*,*) ntilts_max,ntilts,ngate,nray,ntilts
! WRITE(*,*) radar_symbol,radar_lat,radar_lon,radar_elevation
! WRITE(*,*) ngate,nray,fstgat,gatsp
! WRITE(*,*) int_razim
! WRITE(*,*) eleva
! WRITE(*,*) volume
!
OPEN(27,file=trim(regular_file),status='unknown',form='unformatted')
WRITE(27) ntilts_max,ntilts
WRITE(27) radar_symbol,radar_lat,radar_lon,radar_elevation
WRITE(27) ngate,nray,fstgat,gatsp
WRITE(27) int_razim
WRITE(27) eleva
WRITE(27) volume
CLOSE(27)
WRITE(*,*) ' The file:',trim(regular_file), &
'has been written out successfully'
END SUBROUTINE OUTPUT_REGULAR