program stats_loc ! !******** compute statistics from an I/O API file ! compute max of 8hr moving average ! output goes to I/O API file ! modified from stats.F by George Pouliot ! by Dongming Hwang implicit none include 'parms3.ext' ! I/O API constants include 'fdesc3.ext' ! I/O API file description data structure include 'iodecl3.ext' ! I/O API function declarations INTEGER, EXTERNAL :: ENVINT ! ENVINT INTEGER STATUS ! ENVINT status integer :: maxwords integer :: maxwordlength, maxlinelength integer :: nwords parameter (maxwords = 100,maxlinelength=100,maxwordlength=16) character(len=maxwordlength) words(maxwords) character(len=maxlinelength) line character(len=maxwordlength) keyword character(len=maxwordlength) , allocatable :: spec_names(:) character(len=maxwordlength) , allocatable :: species(:) character(len=maxwordlength) , allocatable :: statschars(:) character(len=maxlinelength) :: species2op, default character(len=maxlinelength) :: stats_list integer, allocatable :: spec_lookup(:) integer, allocatable :: stats_lookup(:) integer logdev character(len=8) :: buf1, buf2, buf3 character(len=8) :: buf4, buf5 integer :: osteps, o_day ! !******* file header information for input & output file ! integer :: sdate, stime,tstep,mxrec,nvars integer :: ncols, nrows,nlays, nsteps integer :: nthik, gdtype character(len=16) :: gdname ! !******* file header information for input file ! integer vtypei( MXVARS3 ) ! variable type: M3(INT|REAL|DBLE) character(len=16) vnamei( MXVARS3 ) ! variable names (length MXDLEN3=80) character(len=16) unitsi( MXVARS3 ) ! " units or 'none' (MXDLEN3=80) character(len=80) vdesci( MXVARS3 ) ! " descriptions (MXDLEN3=80) ! !******* file header information for output file ! integer vtypeo( MXVARS3 ) ! variable type: M3(INT|REAL|DBLE) character(len=16) vnameo( MXVARS3 ) ! variable names (length MXDLEN3=80) character(len=16) unitso( MXVARS3 ) ! " units or 'none' (MXDLEN3=80) character(len=80) vdesco( MXVARS3 ) ! " descriptions (MXDLEN3=80) ! !******* environment variables ! integer :: infil_id character(2) infil_ch character(len=16) input, output, this_program,sainput,saoutput character(len=16) species_list, stats_env character(len=60) zone1, zone2, zone3, zone4 ! !******* data storage ! real, allocatable :: data( : , : , : ) real, allocatable :: data_buf( : , : , : , : ) real, allocatable :: sa_data( : , : , : ) real, allocatable :: sa_data_buf( : , : , : , : ) real, allocatable :: eighthr_max( : , : , : ) real, allocatable :: eighthr_sum( : , : , : ) real, allocatable :: eighthr_dat( : , : , : , : ) real, allocatable :: sa_eighthr_max( : , : , : ) real, allocatable :: sa_eighthr_sum( : , : , : ) real, allocatable :: sa_eighthr_dat( : , : , : , : ) integer, allocatable :: zone_col( : , : ) integer, allocatable :: zone_row( : , : ) integer :: zone_ptr, hr_id integer :: num_zn(4), i integer :: z1, z2, z3, z4 integer :: jdate, jtime integer :: tdate, ttime integer :: num_spe, num_stats integer :: rdate, rtime, wdate, wtime ! !******* indexing ! integer :: oindex, iindex ! !******* loop variables ! integer :: v_loop, t_loop, i_loop, j_loop integer :: l_loop integer :: k_loop, s_loop ! !******* error handling ! integer istatus character(len=80) mesg ! !******** external functions ! integer trimlen external trimlen ! ! read input file name ... ! input = 'CMAQ_FILE' infil_id = 1 write(infil_ch, '(I2.2)') infil_id write(*,*)'infil_ch ', infil_ch, 'infil_id ', infil_id input(10:11) = infil_ch write(*,*)"input ", input sainput = 'SA_FILE' ! infil_id = 1 ! write(infil_ch, '(I2.2)') infil_id ! write(*,*)'infil_ch ', infil_ch, 'infil_id ', infil_id sainput(8:9) = infil_ch write(*,*)"sainput ", sainput output = 'STATS_FILE' saoutput = 'SA_STATS_FILE' this_program = 'IOAPI_STATS' species_list = 'SPECIES' stats_env = 'STATS' logdev = init3() ! initialization returns unit # for log ! !******** open input file and get header information ! if ( .not. open3( input, FSREAD3, this_program ) ) & then MESG = 'Could not open file "' // & input( 1: TRIMLEN( input) ) // & '" for input' CALL M3EXIT( this_program, 0, 0, MESG, 2 ) end if if ( .not. desc3(input)) then mesg = 'Could not get description info for file "' // & input( 1: TRIMLEN( input) ) //'"' call m3exit( this_program, 0, 0, mesg, 2 ) end if ncols = NCOLS3D nrows = NROWS3D nvars = NVARS3D nsteps = MXREC3D sdate = SDATE3D stime = STIME3D tstep = TSTEP3D mxrec = MXREC3D nlays = NLAYS3D nthik = NTHIK3D gdtype = GDTYP3D gdname = GDNAM3D !*!*!*! Open SA file if ( .not. open3( sainput, FSREAD3, this_program ) ) & then MESG = 'Could not open file "' // & input( 1: TRIMLEN( sainput) ) // & '" for input' CALL M3EXIT( this_program, 0, 0, MESG, 2 ) end if if ( .not. desc3(sainput)) then mesg = 'Could not get description info for file "' // & input( 1: TRIMLEN( sainput) ) //'"' call m3exit( this_program, 0, 0, mesg, 2 ) end if ! !******* store variable information from input; construct default variable list ! default(1:1) = ' ' do v_loop = 1,nvars vtypei(v_loop) = VTYPE3D(v_loop) vnamei(v_loop) = VNAME3D(v_loop) unitsi(v_loop) = UNITS3D(v_loop) vdesci(v_loop) = VDESC3D(v_loop) if (v_loop .eq. 1) then default = & vnamei(v_loop)(1:trimlen(vnamei(v_loop))) else default = & default(1:trimlen(default))// & vnamei(v_loop)(1:trimlen(vnamei(v_loop))) endif enddo ! !***** get variable list from environment (use default list if not set by user) ! call envstr(species_list,'species to output', & default,species2op,istatus) ! !***** parse the variable list ! call parseline(maxlinelength,maxwords, nwords,maxwordlength, & species2op,words,keyword) num_spe = nwords write (logdev,*) & 'Number of species to put in output:',num_spe write (logdev,*) & 'Number of species to put in output:',num_spe ! !********* allocate storage for variable list for output ! if (num_spe .gt. 0) then allocate (species(num_spe)) allocate (spec_lookup(num_spe)) do i_loop = 1,num_spe species(i_loop) = & words(i_loop)(1:trimlen(words(i_loop))) enddo endif ! !********** create cross reference indexing for output variable list !********** from input variable list ! do i_loop = 1,num_spe do v_loop = 1,nvars if (vnamei(v_loop) .eq. species(i_loop)) then spec_lookup(i_loop) = v_loop endif enddo enddo ! !******** get list of statistics from environtment !******** (use all stats if not supplied by user) ! call envstr(stats_env,'stats to output', & '8HRMX',stats_list,istatus) ! !******** parse the stats character string ! call parseline(maxlinelength,maxwords, nwords,maxwordlength, & stats_list,words,keyword) num_stats = nwords write (logdev,*) & 'Number of stats to put in output:',num_stats write (logdev,*) & 'Number of stats to put in output:',num_stats ! !******** allocate storage to store the list of stats to compute ! if (num_stats .gt. 0) then allocate (statschars(num_stats)) allocate (stats_lookup(num_stats)) do i_loop = 1,num_stats statschars(i_loop) = & words(i_loop)(1:trimlen(words(i_loop))) enddo endif ! !********* create cross-reference index for stats list ! do i_loop = 1,num_stats if ( statschars(i_loop) .eq. 'HRMAX') then stats_lookup(i_loop) = 1 elseif (statschars(i_loop) .eq. '8HRMX') then stats_lookup(i_loop) = 2 elseif (statschars(i_loop) .eq. 'HRAVG') then stats_lookup(i_loop) = 5 elseif (statschars(i_loop) .eq. 'DYSUM') then stats_lookup(i_loop) = 6 endif enddo ! !********** setup of header information for output file ! NCOLS3D = ncols NROWS3D = nrows NVARS3D = num_stats*num_spe c... MXREC3D = (nsteps-offset)/daylength TSTEP3D = tstep NLAYS3D = nlays GDTYP3D = gdtype GDNAM3D = gdname NTHIK3D = nthik ! !******** set up output variable name descriptions and units ! do v_loop = 1,num_spe ! loop over variable names oindex = num_stats*(v_loop-1) iindex = spec_lookup(v_loop) do j_loop = 1,num_stats ! loop for each type of statistic VNAME3D(oindex+j_loop) = vnamei(iindex)(1: & trimlen(vnamei(iindex)))//statschars(j_loop) VTYPE3D(num_stats*(v_loop-1)+j_loop) = M3REAL UNITS3D(num_stats*(v_loop-1)+j_loop) = unitsi(iindex) VDESC3D(num_stats*(v_loop-1)+j_loop) = vdesci(iindex) enddo enddo do i_loop = 1, num_spe*num_stats vnameo(i_loop) = VNAME3D(i_loop) vtypeo(i_loop) = VTYPE3D(i_loop) unitso(i_loop) = UNITS3D(i_loop) vdesco(i_loop) = VDESC3D(i_loop) enddo ! !******** open output file ! TSTEP3D = 240000 if ( .not. open3( output, FSCREA3, this_program ) ) & then MESG = 'Could not open file "' // & output( 1: TRIMLEN( output) ) // & '" for output' CALL M3EXIT( this_program, 0, 0, MESG, 2 ) end if !*!*!*!*! Open SA output file if ( .not. open3( saoutput, FSCREA3, this_program ) ) & then MESG = 'Could not open file "' // & output( 1: TRIMLEN( saoutput) ) // & '" for output' CALL M3EXIT( this_program, 0, 0, MESG, 2 ) end if ! !******** allocate storage for data ! allocate (data(ncols,nrows,nlays)) allocate (sa_data(ncols,nrows,nlays)) allocate (eighthr_dat(ncols,nrows,nlays,8)) allocate (eighthr_max(ncols,nrows,nlays)) allocate (eighthr_sum(ncols,nrows,nlays)) allocate (sa_eighthr_dat(ncols,nrows,nlays,8)) allocate (sa_eighthr_max(ncols,nrows,nlays)) allocate (sa_eighthr_sum(ncols,nrows,nlays)) allocate (data_buf(ncols,nrows,nlays,4)) allocate (sa_data_buf(ncols,nrows,nlays,4)) allocate (zone_col(4,ncols*nrows)) allocate (zone_row(4,ncols*nrows)) c input zone grid c inpu zone definition read(10,*) buf1, buf2, buf3, buf4, buf5 write(*,*)' line ', buf1, buf2, buf3, buf4, buf5 do i = 1, 4 num_zn(i) = 0 enddo z1 = 1 z2 = 1 z3 = 1 z4 = 1 do i = 1, ncols*nrows read(10,*,end=112) i_loop, j_loop, k_loop, buf1, t_loop if ( t_loop .eq. -8 ) then zone_col(4, z4) = i_loop zone_row(4, z4) = j_loop z4 = z4+1 num_zn(4) = num_zn(4)+1 else if ( t_loop .eq. -7 ) then zone_col(3, z3) = i_loop zone_row(3, z3) = j_loop z3 = z3+1 num_zn(3) = num_zn(3)+1 else if ( t_loop .eq. -6 ) then zone_col(2, z2) = i_loop zone_row(2, z2) = j_loop z2 = z2+1 num_zn(2) = num_zn(2)+1 else if ( t_loop .eq. -5 ) then zone_col(1, z1) = i_loop zone_row(1, z1) = j_loop z1 = z1+1 num_zn(1) = num_zn(1)+1 end if enddo 112 continue write(*,*)'zone summary: Eastern:Centeral:Mountain:Pacific', * (num_zn(i),i=1,4) write(*,*)' z1 ', num_zn(1), 'col row', * zone_col(1,num_zn(1)), zone_row(1,num_zn(1)) write(*,*)' z2 ', num_zn(2), 'col row', * zone_col(2,num_zn(2)), zone_row(2,num_zn(2)) write(*,*)' z3 ', num_zn(3), 'col row', * zone_col(3,num_zn(3)), zone_row(3,num_zn(3)) write(*,*)' z4 ', num_zn(4), 'col row', * zone_col(4,num_zn(4)), zone_row(4,num_zn(4)) o_day = ENVINT( 'OUT_DAY', ' ', 0, STATUS ) osteps = o_day*24 write(*,*)'total output # of day:hour ', o_day, osteps c loop on species do v_loop = 1,num_spe c get jdate/jtime : output date/time jdate = sdate jtime = stime c get rdate/rtime : reading date/time rdate = jdate rtime = stime + 50000 ! 00 eastern time = GMT + 5 wdate = jdate wtime = 0 ! local time c----------------------------------------------------------- c fetch 7 hours data c first hour data: pre-fetch 4 hours data hr_id = 1 do t_loop = 1, 4 if ( .not. & read3(input,vnamei(spec_lookup(v_loop)), & ALLAYS3,rdate,rtime,data(1,1,1)) & ) then infil_id = infil_id+1 write(infil_ch, '(I2.2)') infil_id input(10:11) = infil_ch sainput(8:9) = infil_ch if ( .not. open3( input, FSREAD3, this_program ) ) & then MESG = 'Could not open file "' // & input( 1: TRIMLEN( input) ) // & '" for input' CALL M3EXIT( this_program, 0, 0, MESG, 2 ) end if if ( .not. & read3(input,vnamei(spec_lookup(v_loop)), & ALLAYS3,rdate,rtime,data(1,1,1)) & ) then mesg = 'Error reading '// & vnamei(spec_lookup(v_loop))// & 'from file '// & input( 1: TRIMLEN( input ) ) call m3exit( this_program, 0, 0, mesg, 2 ) end if end if do i_loop = 1, ncols do j_loop = 1, nrows do k_loop = 1, nlays data_buf(i_loop, j_loop, k_loop, t_loop) = & data(i_loop, j_loop, k_loop) enddo enddo enddo !*!*!*!*!*! SA if ( .not. & read3(sainput,vnamei(spec_lookup(v_loop)), & ALLAYS3,rdate,rtime,sa_data(1,1,1)) & ) then ! infil_id = infil_id+1 write(infil_ch, '(I2.2)') infil_id input(10:11) = infil_ch sainput(8:9) = infil_ch if ( .not. open3( sainput, FSREAD3, this_program ) ) & then MESG = 'Could not open file "' // & input( 1: TRIMLEN( sainput) ) // & '" for input' CALL M3EXIT( this_program, 0, 0, MESG, 2 ) end if if ( .not. & read3(sainput,vnamei(spec_lookup(v_loop)), & ALLAYS3,rdate,rtime,sa_data(1,1,1)) & ) then mesg = 'Error reading '// & vnamei(spec_lookup(v_loop))// & 'from file '// & input( 1: TRIMLEN( sainput ) ) call m3exit( this_program, 0, 0, mesg, 2 ) end if end if do i_loop = 1, ncols do j_loop = 1, nrows do k_loop = 1, nlays sa_data_buf(i_loop, j_loop, k_loop, t_loop) = & sa_data(i_loop, j_loop, k_loop) enddo enddo enddo !*!*!*! call nextime(rdate,rtime,tstep) end do ! end t_loop for first hour data write(*,*)' complete first four hour data fetch ...' zone_ptr = 1 call zone_merge(eighthr_dat, data_buf, & ncols, nrows, nlays, * zone_ptr, * zone_col, zone_row, num_zn, hr_id) call sa_zone_merge(sa_eighthr_dat, sa_data_buf, & ncols, nrows, nlays, * zone_ptr, * zone_col, zone_row, num_zn, hr_id) write(*,*)' rdate:rtime ------------- ',rdate, rtime c get next 6 hour data do l_loop = 2, 7 if ( .not. & read3(input,vnamei(spec_lookup(v_loop)), & ALLAYS3,rdate,rtime,data(1,1,1)) & ) then c ... end of this day, continue next day ... infil_id = infil_id+1 write(infil_ch, '(I2.2)') infil_id input(10:11) = infil_ch sainput(8:9) = infil_ch if ( .not. open3( input, FSREAD3, this_program ) ) & then MESG = 'Could not open file "' // & input( 1: TRIMLEN( input) ) // & '" for input' CALL M3EXIT( this_program, 0, 0, MESG, 2 ) end if if ( .not. & read3(input,vnamei(spec_lookup(v_loop)), & ALLAYS3,rdate,rtime,data(1,1,1)) & ) then mesg = 'Error reading '// & vnamei(spec_lookup(v_loop))// & 'from file '// & input( 1: TRIMLEN( input ) ) call m3exit( this_program, 0, 0, mesg, 2 ) end if end if do i_loop = 1, ncols do j_loop = 1, nrows do k_loop = 1, nlays data_buf(i_loop, j_loop, k_loop, zone_ptr) = & data(i_loop, j_loop, k_loop) enddo enddo enddo !*!*!*!*! SA part if ( .not. & read3(sainput,vnamei(spec_lookup(v_loop)), & ALLAYS3,rdate,rtime,sa_data(1,1,1)) & ) then c ... end of this day, continue next day ... ! infil_id = infil_id+1 write(infil_ch, '(I2.2)') infil_id input(10:11) = infil_ch sainput(8:9) = infil_ch if ( .not. open3( sainput, FSREAD3, this_program ) ) & then MESG = 'Could not open file "' // & input( 1: TRIMLEN( sainput) ) // & '" for input' CALL M3EXIT( this_program, 0, 0, MESG, 2 ) end if if ( .not. & read3(sainput,vnamei(spec_lookup(v_loop)), & ALLAYS3,rdate,rtime,sa_data(1,1,1)) & ) then mesg = 'Error reading '// & vnamei(spec_lookup(v_loop))// & 'from file '// & input( 1: TRIMLEN( sainput ) ) call m3exit( this_program, 0, 0, mesg, 2 ) end if end if do i_loop = 1, ncols do j_loop = 1, nrows do k_loop = 1, nlays sa_data_buf(i_loop, j_loop, k_loop, zone_ptr) = & sa_data(i_loop, j_loop, k_loop) enddo enddo enddo !*!*!*! call nextime(rdate,rtime,tstep) c update ptr zone_ptr = mod(zone_ptr+1, 4) if (zone_ptr .eq. 0) zone_ptr = 4 hr_id = hr_id + 1 if (hr_id .gt. 8) hr_id = mod(hr_id,8) call zone_merge(eighthr_dat, data_buf, & ncols, nrows, nlays, * zone_ptr, * zone_col, zone_row, num_zn, hr_id) call sa_zone_merge(sa_eighthr_dat, sa_data_buf, & ncols, nrows, nlays, * zone_ptr, * zone_col, zone_row, num_zn, hr_id) end do ! end 7 hour fetch write(*,*)'complete 7 hour fetch, begin 8hr processing ...' c----------------------------------------------------------- c begin 8 hr process do t_loop = 1,osteps write (logdev,*) t_loop, jdate,jtime if ( .not. & read3(input,vnamei(spec_lookup(v_loop)), & ALLAYS3,rdate,rtime,data(1,1,1)) & ) then c ... end of this day, continue next day ... infil_id = infil_id+1 write(infil_ch, '(I2.2)') infil_id input(10:11) = infil_ch if ( .not. open3( input, FSREAD3, this_program ) ) & then MESG = 'Could not open file "' // & input( 1: TRIMLEN( input) ) // & '" for input' CALL M3EXIT( this_program, 0, 0, MESG, 2 ) end if if ( .not. & read3(input,vnamei(spec_lookup(v_loop)), & ALLAYS3,rdate,rtime,data(1,1,1)) & ) then mesg = 'Error reading '// & vnamei(spec_lookup(v_loop))// & 'from file '// & input( 1: TRIMLEN( input ) ) call m3exit( this_program, 0, 0, mesg, 2 ) end if end if do i_loop = 1, ncols do j_loop = 1, nrows do k_loop = 1, nlays data_buf(i_loop, j_loop, k_loop, zone_ptr) = & data(i_loop, j_loop, k_loop) enddo enddo enddo !*!*!*!*! SA part if ( .not. & read3(sainput,vnamei(spec_lookup(v_loop)), & ALLAYS3,rdate,rtime,sa_data(1,1,1)) & ) then c ... end of this day, continue next day ... ! infil_id = infil_id+1 write(infil_ch, '(I2.2)') infil_id input(10:11) = infil_ch sainput(8:9) = infil_ch if ( .not. open3( sainput, FSREAD3, this_program ) ) & then MESG = 'Could not open file "' // & input( 1: TRIMLEN( sainput) ) // & '" for input' CALL M3EXIT( this_program, 0, 0, MESG, 2 ) end if if ( .not. & read3(sainput,vnamei(spec_lookup(v_loop)), & ALLAYS3,rdate,rtime,sa_data(1,1,1)) & ) then mesg = 'Error reading '// & vnamei(spec_lookup(v_loop))// & 'from file '// & input( 1: TRIMLEN( sainput ) ) call m3exit( this_program, 0, 0, mesg, 2 ) end if end if do i_loop = 1, ncols do j_loop = 1, nrows do k_loop = 1, nlays sa_data_buf(i_loop, j_loop, k_loop, zone_ptr) = & sa_data(i_loop, j_loop, k_loop) enddo enddo enddo !*!*!*!*! call nextime(rdate,rtime,tstep) zone_ptr = mod(zone_ptr+1, 4) if (zone_ptr .eq. 0) zone_ptr = 4 hr_id = hr_id + 1 if (hr_id .gt. 8) hr_id = mod(hr_id,8) call zone_merge(eighthr_dat, data_buf, & ncols, nrows, nlays, * zone_ptr, * zone_col, zone_row, num_zn, hr_id) call sa_zone_merge(sa_eighthr_dat, sa_data_buf, & ncols, nrows, nlays, * zone_ptr, * zone_col, zone_row, num_zn, hr_id) c 8 hr max eighthr_sum(1:ncols,1:nrows,1:nlays) = 0.0 sa_eighthr_sum(1:ncols,1:nrows,1:nlays) = 0.0 do l_loop = 1,8 do k_loop = 1,nlays do j_loop = 1,nrows do i_loop = 1,ncols eighthr_sum(i_loop, j_loop, k_loop) = & eighthr_sum(i_loop, j_loop, k_loop) + & eighthr_dat(i_loop, j_loop, k_loop,l_loop) sa_eighthr_sum(i_loop, j_loop, k_loop) = & sa_eighthr_sum(i_loop, j_loop, k_loop) + & sa_eighthr_dat(i_loop, j_loop, k_loop,l_loop) enddo enddo enddo enddo do k_loop = 1, nlays do j_loop = 1,nrows do i_loop = 1,ncols eighthr_sum(i_loop, j_loop, k_loop) = & eighthr_sum(i_loop, j_loop, k_loop)/8.0 sa_eighthr_sum(i_loop, j_loop, k_loop) = & sa_eighthr_sum(i_loop, j_loop, k_loop)/8.0 enddo enddo enddo c ... starting of the day if (jtime .eq. 0) then eighthr_max(1:ncols,1:nrows,1:nlays) & = eighthr_sum(1:ncols, 1:nrows, 1:nlays) sa_eighthr_max(1:ncols,1:nrows,1:nlays) & = sa_eighthr_sum(1:ncols, 1:nrows, 1:nlays) endif do k_loop = 1, nlays do j_loop = 1,nrows do i_loop = 1,ncols if (eighthr_sum(i_loop, j_loop, k_loop) .gt. & eighthr_max(i_loop, j_loop, k_loop)) then eighthr_max(i_loop, j_loop, k_loop) = & eighthr_sum(i_loop, j_loop, k_loop) !! keep SA 8hr max at same time as regular model output time 8hr max sa_eighthr_max(i_loop, j_loop, k_loop) = & sa_eighthr_sum(i_loop, j_loop, k_loop) endif enddo enddo enddo write(*,*)' rdate:rtime output ',rdate, rtime if (mod(t_loop,24) .eq. 0) then if ( .not. & write3(output, & vnameo(num_stats*(v_loop-1)+1), & wdate,wtime, & eighthr_max(1,1,1))) then mesg = 'Error writing 8 hourly_max to file '// & output( 1: TRIMLEN( output ) ) call m3exit( this_program, 0, 0, mesg, 2 ) end if if ( .not. & write3(saoutput, & vnameo(num_stats*(v_loop-1)+1), & wdate,wtime, & sa_eighthr_max(1,1,1))) then mesg = 'Error writing 8 hourly_max to file '// & output( 1: TRIMLEN( saoutput ) ) call m3exit( this_program, 0, 0, mesg, 2 ) end if call nextime(wdate, wtime, 24*tstep) ! next time step end if call nextime(jdate, jtime, tstep) ! next time step enddo enddo ! !******** must shut down IO/API before exiting (necessary when writing to IO/API files) ! if ( .not. shut3() ) then write (*,*) 'Error shutting down IOAPI' stop 'ERROR: Abnormal Termination' endif write (logdev,*) this_program,': Normal Completion' stop 990 end subroutine updatemax(ncols,nrows, nlays, & maxvals,data) implicit none integer ncols, nrows, nlays real maxvals(ncols,nrows,nlays) real data(ncols, nrows, nlays) integer i_loop, j_loop, k_loop do k_loop = 1, nlays do j_loop = 1, nrows do i_loop = 1, ncols if (data(i_loop, j_loop, k_loop) .gt. & maxvals(i_loop,j_loop, k_loop)) then maxvals(i_loop, j_loop, k_loop)= & data(i_loop, j_loop, k_loop) endif enddo enddo enddo return end subroutine parseline(maxlinelength,maxwords, nwords, & maxwordlength, & line,words,keyword) ! !assumed delimiter between words is one or more blanks ' ' ! integer maxwords, maxwordlength, nwords, maxlinelength character(len=maxwordlength) words(maxwords) character(len=maxlinelength) line character(len=maxwordlength) keyword ! !******** locals ! integer tlength,lbks, istart, iend, i_loop logical inword !true if in a word, false otherwise character*1 char,prevchar ! !******** external function calls ! external lblank, trimlen integer trimlen, lblank ! ! !******* take a string of text (line) and split it into words ! subject to the following limitations: ! the maximum number of words in the string is maxwords ! the maximum lenght of any word is maxwordlength ! the maximum length of line is maxlinelength tlength = trimlen(line(1:maxlinelength)) if (tlength .eq. maxlinelength) then write (*,*) 'ERROR: MAXLINELENGTH REACHED' write (*,*) 'increase the maximum length of the record' stop endif lbks = lblank(line(1:maxlinelength)) nwords = 0 inword = .true. istart = lbks+1 iend = lbks+1 prevchar = ' ' do i_loop = lbks+1,tlength+1 char = line(i_loop:i_loop) if ((char .eq. ' ') .and. (prevchar .ne. ' ')) then inword = .false. nwords = nwords + 1 words(nwords) = line(istart:iend) else if((char .eq. ' ') .and. (prevchar .eq. ' ')) then ! !************ do nothing ! else if (inword) then iend = iend + 1 else inword = .true. istart = i_loop iend = i_loop endif endif endif prevchar = char enddo keyword = words(1) keyword = keyword(2:trimlen(keyword)) return end c merge data subroutine zone_merge(eighthr_dat, data_buf, & ncols, nrows, nlays, * zone_ptr, * zone_col, zone_row, num_zn, hr_id) real eighthr_dat(ncols, nrows, nlays, 8) real data_buf(ncols, nrows, nlays,4) integer zone_col(4, ncols*nrows) integer zone_row(4, ncols*nrows) integer zone_ptr, hr_id integer num_zn(4) integer i_loop, j_loop, k_loop integer i, zn integer c, r c merge data do i = 0, 3 zn = mod(zone_ptr+i-1, 4) + 1 do k_loop = 1, nlays do j_loop = 1, num_zn(i+1) c = zone_col(i+1, j_loop) r = zone_row(i+1, j_loop) eighthr_dat(c, r, k_loop, hr_id) = & data_buf(c, r, k_loop, zn) enddo enddo enddo return end subroutine sa_zone_merge(sa_eighthr_dat, sa_data_buf, & ncols, nrows, nlays, * zone_ptr, * zone_col, zone_row, num_zn, hr_id) real sa_eighthr_dat(ncols, nrows, nlays, 8) real sa_data_buf(ncols, nrows, nlays,4) integer zone_col(4, ncols*nrows) integer zone_row(4, ncols*nrows) integer zone_ptr, hr_id integer num_zn(4) integer i_loop, j_loop, k_loop integer i, zn integer c, r c merge data do i = 0, 3 zn = mod(zone_ptr+i-1, 4) + 1 do k_loop = 1, nlays do j_loop = 1, num_zn(i+1) c = zone_col(i+1, j_loop) r = zone_row(i+1, j_loop) sa_eighthr_dat(c, r, k_loop, hr_id) = & sa_data_buf(c, r, k_loop, zn) enddo enddo enddo return end