C----extract a IOAPI file for a certain time period include 'PARMS3.EXT' ! i/o API include 'FDESC3.EXT' ! i/o API include 'IODECL3.EXT' ! i/o API parameter (maxspecies=100) character*300 afile1, afile2,aline character varaerosol(maxspecies)*16 integer begdate,begtime,enddate,endtime,nowtstep real, dimension (:,:,:) :: work, temp, press namelist /control/begdate,begtime,enddate,endtime,varaerosol real c303,c302 parameter(C303=19.83,C302=5417.4) c --- end WORK_AREAS declarations ESAT(TEMK)=.611*EXP(C303-C302/TEMK) ! for calculating saturated water vapor pressure QSAT(ESAT1,PCB)=ESAT1*.622/(PCB-ESAT1) ! TEMK is ambient temperature in K, PCB is the pressue in KPa ! QSAT is the saturated humidity in kg/kg c logical dele varaerosol(1:maxspecies)=' ' open(7,file='extract-1.ini') read(7,control) do L=1,maxspecies if(varaerosol(L).eq.' ') exit enddo nspecies=L-1 begtime=begtime*10000 ! convert to HHMMSS endtime=endtime*10000 if(.not.OPEN3('INPUT',FSREAD3,'ave-ioapi')) stop if(.not.DESC3('INPUT')) stop imax=ncols3d jmax=nrows3d kmax=nlays3d nowtstep=tstep3d c-----create the file nvars3d=nspecies+3 do L=1,nspecies vname3d(L)=varaerosol(L) units3d(L)='ug/m3' vdesc3d(L)='STEM concentration of '//trim(vname3d(L)) enddo vname3d(nspecies+1)='P' units3d(nspecies+1)='Pa' vdesc3d(nspecies+1)='Pressure' vname3d(nspecies+2)='T' units3d(nspecies+2)='K' vdesc3d(nspecies+2)='Temperature' vname3d(nspecies+3)='RH' units3d(nspecies+3)='%' vdesc3d(nspecies+3)='Relative Humidity' sdate3d=begdate stime3d=begtime if(.not.OPEN3('OUTPUT',FSUNKN3,'ave-ioapi')) stop if(.not.OPEN3('METEO3D',FSREAD3,'extract-1')) stop if(.not.DESC3('METEO3D')) stop if(imax.ne.ncols3d.or.jmax.ne.nrows3d.or.kmax.gt.nlays3d) then print*,'meteo3d file dimension does not match' stop endif kmet=nlays3d allocate(work(imax,jmax,kmet),STAT=ierr) allocate(press(imax,jmax,kmet),STAT=ierr) allocate(temp(imax,jmax,kmet),STAT=ierr) ****** begin read and output nowdate=begdate nowtime=begtime do while(.not.(nowdate.ge.enddate.and.nowtime.gt.endtime)) do L=1,nspecies if(.not.read3('INPUT',varaerosol(L),ALLAYS3,nowdate,nowtime, 1 work)) stop if(.not.write3('OUTPUT',varaerosol(L),nowdate,nowtime, 1 work)) stop enddo if(.not.INTERP3('METEO3D','P','extract-1',nowdate,nowtime, 1 imax*jmax*kmet,press)) stop if(.not.write3('OUTPUT','P',nowdate,nowtime, 1 press)) stop if(.not.INTERP3('METEO3D','T','extract-1',nowdate,nowtime, 1 imax*jmax*kmet,temp)) stop if(.not.write3('OUTPUT','T',nowdate,nowtime, 1 temp)) stop if(.not.INTERP3('METEO3D','VAPOR','extract-1',nowdate,nowtime, ! water vapor 1 imax*jmax*kmet,work)) stop work(1:imax,1:jmax,1:kmet)=100.*work(1:imax,1:jmax,1:kmet)/ ! computing relative humudity 1 QSAT(ESAT(temp(1:imax,1:jmax,1:kmet)), 2 press(1:imax,1:jmax,1:kmet)/1000) ! convert to KPa if(.not.write3('OUTPUT','RH',nowdate,nowtime, 1 work)) stop call nextime(nowdate,nowtime,nowtstep) enddo ! end of time loop *** construct output file header iflag=shut3() end