1 2 SUBROUTINE MEREAD 3 !*********************************************************************** 4 ! MEREAD Module of the AMS/EPA Regulatory Model - AERMOD - EVENT 5 ! 6 ! PURPOSE: Controls Extraction and Quality Assurance of 7 ! One Day of Meteorological Data 8 ! 9 ! PROGRAMMER: ROGER BRODE, JEFF WANG 10 ! 11 ! DATE: March 2, 1992 12 ! 13 ! MODIFIED: To remove support for unformatted meteorological 14 ! data files. 15 ! R.W. Brode, PES, Inc., 4/10/2000 16 ! 17 ! MODIFIED: To incorporate modifications to date processing 18 ! for Y2K compliance, including use of date window 19 ! variables (ISTRT_WIND and ISTRT_CENT) and calculation 20 ! of 10-digit date variable (FULLDATE) with 4-digit 21 ! year for date comparisons. 22 ! Also modified calls to METDAT insteaad of EV_METDAT 23 ! to allow use of same routine for both normal and 24 ! EVENT processing. 25 ! R.W. Brode, PES, Inc., 5/12/99 26 ! 27 ! INPUTS: Meteorology File Specifications 28 ! 29 ! OUTPUTS: Arrays of Meteorological Variables for One Day 30 ! 31 ! CALLED FROM: EVLOOP 32 !*********************************************************************** 33 34 ! Variable Declarations 35 USE MAIN1 36 IMPLICIT NONE 37 CHARACTER MODNAM*12 38 39 SAVE 40 REAL :: DAY , AFVM1 , C1 , C2 , C3 , STEFB , RN , ES25 41 INTEGER :: I , J , IHR , IJDAY , IDATCHK , JUYI , JUSI , JSYI , & 42 & JSSI , JFLAG , LEVEL , METVER 43 CHARACTER(LEN=8) :: CUSI , CSSI , COSI 44 CHARACTER(LEN=6) :: SPEC1 , SPEC2 , SPEC3 45 CHARACTER(LEN=132) :: BUFFER 46 47 ! Variable Initializations 48 MODNAM = 'MEREAD' ! 0 49 PATH = 'MX' 50 51 !---- Constants used in the computation of QSW 52 C1 = 5.31E-13 53 C2 = 60.0 54 C3 = 1.12 55 STEFB = 5.67E-08 56 57 ! READ Meteorology Data Based on Format -- 58 ! When DRY deposition is modeled, U-star, L, and z0 (surface 59 ! roughness length) are read in addition to the standard RAMMET 60 ! data. These must be provided at the end of each hourly record 61 ! for the FORMATTED ASCII, CARD, and FREE options. 62 ! 63 ! When WET deposition is modeled, ipcode (precip. 64 ! code) and prate (precip. rate in mm/hr) must also be added to 65 ! each hourly record. 66 ! The format statement allows for all additional data: 67 68 ! Calculate the MMDDHH variable to check for end of the year 69 IDATCHK = KURDAT - INT(KURDAT/1000000)*1000000 70 IF ( (IMONTH.EQ.12 .AND. IDAY.EQ.31 .AND. IHOUR.EQ.24) .OR. & 71 & IDATCHK.EQ.123124 ) THEN 72 ! End of year has been reached - check for presence of header 73 ! record at beginning of next year for multi-year data files. 74 READ (MFUNIT,'(A132)',ERR=998,END=1000,IOSTAT=IOERRN) BUFFER ! 0 75 READ (BUFFER,1900,ERR=998,IOSTAT=IOERRN) ALAT , ALON , SPEC1 , & 76 & CUSI , SPEC2 , CSSI , SPEC3 , COSI , METVER 77 1900 FORMAT (2A10,T31,A6,T38,A8,T48,A6,T55,A8,T65,A6,T72,A8,T94,I5) 78 ! Convert character IDs to integers 79 CALL STONUM(CUSI,8,FNUM,IMIT) ! 0 80 IF ( IMIT.EQ.1 ) THEN 81 JUSI = NINT(FNUM) ! 0 82 ELSE 83 JUSI = 0 ! 0 84 ENDIF 85 CALL STONUM(CSSI,8,FNUM,IMIT) ! 0 86 IF ( IMIT.EQ.1 ) THEN 87 JSSI = NINT(FNUM) ! 0 88 ELSE 89 JSSI = 0 ! 0 90 ENDIF 91 92 IF ( JSSI.NE.IDSURF .OR. JUSI.NE.IDUAIR ) THEN ! 0 93 ! Station IDs don't match runstream input, assume that header 94 ! record is missing. Backspace met file and continue processing. 95 BACKSPACE MFUNIT ! 0 96 ELSEIF ( INDEX(BUFFER,':').EQ.0 ) THEN 97 ! Station IDs match, but record does not contain colon. 98 ! Assume it must be regular met data record, so backspace met file. 99 BACKSPACE MFUNIT ! 0 100 ENDIF 101 102 GOTO 1001 ! 0 103 104 ! Error reading 'header record' - assume that header record is 105 ! missing. Backspace met file and continue processing. 106 998 BACKSPACE MFUNIT ! 0 107 108 ENDIF 109 110 1001 CONTINUE ! 0 111 112 HOUR_LOOP:DO IHR = 1 , NHR 113 ! 114 !---- READ surface scaling meteorology data based on format 115 ! 116 IF ( LDPART .OR. LWPART .OR. LDGAS .OR. LWGAS ) THEN ! 0 117 ! Read record from ASCII scalar parameter file using FREE format 118 ! with deposition variables 119 ! 120 READ (MFUNIT,*,END=1000,ERR=99,IOSTAT=IOERRN) IYEAR , & 121 & IMONTH , IDAY , IJDAY , IHOUR , ASFCHF(IHR) , & 122 & AUSTAR(IHR) , AWSTAR(IHR) , AVPTGZI(IHR) , & 123 & AZICONV(IHR) , AZIMECH(IHR) , AOBULEN(IHR) , & 124 & ASFCZ0(IHR) , ABOWEN(IHR) , AALBEDO(IHR) , AUREF(IHR) & 125 & , AWDREF(IHR) , AUREFHT(IHR) , ATA(IHR) , ATREFHT(IHR)& 126 & , IAPCODE(IHR) , APRATE(IHR) , ARH(IHR) , ASFCP(IHR) ,& 127 & NACLOUD(IHR) 128 ! 129 ! Calculate solar irradiance, QSW, from Heat Flux, Bowen ratio, 130 ! albedo and cloud cover, for use in gas deposition algorithm. 131 IF ( ASFCHF(IHR).LT.0.0 .OR. ATA(IHR).LT.0.0 .OR. & 132 & AOBULEN(IHR).EQ.-99999.0 ) THEN 133 ! Hour is stable or missing 134 AQSW(IHR) = 0.0 ! 0 135 ELSE 136 RN = (1.+1./ABOWEN(IHR))*ASFCHF(IHR)/0.9 ! 0 137 AQSW(IHR) = (RN*(1.+C3)-C1*ATA(IHR)**6+STEFB*ATA(IHR) & 138 & **4-C2*0.1*NACLOUD(IHR))/(1.-AALBEDO(IHR)) 139 ENDIF 140 ! 141 ! Save precipitation rates for two previous hours 142 IF ( IHR.EQ.1 ) THEN ! 0 143 APREC2(IHR) = APRATE(NHR-1) ! 0 144 APREC1(IHR) = APRATE(NHR) 145 ELSEIF ( IHR.EQ.2 ) THEN 146 APREC2(IHR) = APRATE(NHR) ! 0 147 APREC1(IHR) = APRATE(IHR-1) 148 ELSE 149 APREC2(IHR) = APRATE(IHR-2) ! 0 150 APREC1(IHR) = APRATE(IHR-1) 151 ENDIF 152 153 ! Set variables for dry deposition 154 IF ( LDPART .OR. LDGAS ) THEN ! 0 155 IF ( ATA(IHR).LT.0.0 .OR. APRATE(IHR).LT.0.0 ) THEN ! 0 156 AWNEW(IHR) = AWOLD(IHR) ! 0 157 ELSE 158 ! ... Compute saturation vapor pressure based on CMAQ formula 159 AESTA(IHR) = 0.6112*EXP(19.83-5417.4/ATA(IHR)) ! 0 160 ES25 = 3.167 161 AWNEW(IHR) = WOLD + APREC1(IHR) - 0.5*F2*AESTA(IHR) & 162 & /ES25 163 WOLD = AWNEW(IHR) 164 AF2(IHR) = AWNEW(IHR)/200. 165 IF ( AF2(IHR).LE.0.01 ) AF2(IHR) = 0.01 166 IF ( AF2(IHR).GT.1.0 ) AF2(IHR) = 1.0 167 F2 = AF2(IHR) 168 ENDIF 169 ENDIF 170 171 ELSE 172 ! Read record from ASCII scalar parameter file without deposition 173 ! parameters, using FREE format 174 ! 175 READ (MFUNIT,*,END=1000,ERR=99,IOSTAT=IOERRN) IYEAR , & 176 & IMONTH , IDAY , IJDAY , IHOUR , ASFCHF(IHR) , & 177 & AUSTAR(IHR) , AWSTAR(IHR) , AVPTGZI(IHR) , & 178 & AZICONV(IHR) , AZIMECH(IHR) , AOBULEN(IHR) , & 179 & ASFCZ0(IHR) , ABOWEN(IHR) , AALBEDO(IHR) , AUREF(IHR) & 180 & , AWDREF(IHR) , AUREFHT(IHR) , ATA(IHR) , ATREFHT(IHR) 181 ! 182 ENDIF 183 184 ! Set the stability logical variables 185 IF ( AOBULEN(IHR).GT.0.0 ) THEN ! 0 186 UNSTAB = .FALSE. ! 0 187 STABLE = .TRUE. 188 ELSE 189 UNSTAB = .TRUE. ! 0 190 STABLE = .FALSE. 191 ENDIF 192 193 !---- Initialize the profile data to missing; 194 ! READ profile data based on format 195 ! 196 CALL PFLINI() ! 0 197 LEVEL = 1 198 JFLAG = 0 199 IF ( PROFRM.EQ.'FREE' ) THEN 200 ! Read record from ASCII profile file using FREE format; compute 201 ! sigma_V from sigma_A and wind speed 202 203 DO WHILE ( JFLAG.EQ.0 ) ! 0 204 READ (MPUNIT,*,END=1000,ERR=98,IOSTAT=IOERRN) KYEAR , & 205 & KMONTH , KDAY , KHOUR , PFLHT(LEVEL) , JFLAG , & 206 & PFLWD(LEVEL) , PFLWS(LEVEL) , PFLTA(LEVEL) , & 207 & PFLSA(LEVEL) , PFLSW(LEVEL) 208 209 ! Convert the data to the required units 210 CALL PFLCNV(LEVEL) ! 0 211 212 ! Set the number of profile levels to current index, store 213 ! the 'top of profile' flag, and increment level if not at top 214 ! Check that the level does not exceed the maximum allowable 215 NPLVLS = LEVEL 216 ANPLVLS(IHR) = LEVEL 217 AIFLAG(IHR,LEVEL) = JFLAG 218 APFLHT(IHR,LEVEL) = PFLHT(LEVEL) 219 APFLWD(IHR,LEVEL) = PFLWD(LEVEL) 220 APFLWS(IHR,LEVEL) = PFLWS(LEVEL) 221 APFLTA(IHR,LEVEL) = PFLTA(LEVEL) 222 APFLSA(IHR,LEVEL) = PFLSA(LEVEL) 223 APFLSV(IHR,LEVEL) = PFLSV(LEVEL) 224 APFLSW(IHR,LEVEL) = PFLSW(LEVEL) 225 IF ( JFLAG.EQ.0 ) THEN 226 LEVEL = LEVEL + 1 ! 0 227 228 IF ( LEVEL.GT.MXPLVL ) THEN 229 IF ( .NOT.PFLERR ) THEN ! 0 230 ! WRITE Error Message: Number of profile levels 231 ! exceeds maximum allowable 232 WRITE (DUMMY,'(I8)') MXPLVL ! 0 233 CALL ERRHDL(PATH,MODNAM,'E','465',DUMMY) 234 PFLERR = .TRUE. 235 RUNERR = .TRUE. 236 ENDIF 237 238 ! Limit the number of levels to the maximum allowable 239 LEVEL = MXPLVL ! 0 240 ENDIF 241 242 ENDIF 243 244 ENDDO 245 246 ! Compute the vertical potential temperature gradient profile 247 IF ( .NOT.RUNERR ) THEN ! 0 248 NTGLVL = 0 ! 0 249 CALL COMPTG() 250 ANTGLVL(IHR) = NTGLVL 251 DO I = 1 , NTGLVL 252 APFLTG(IHR,I) = PFLTG(I) ! 0 253 APFLTGZ(IHR,I) = PFLTGZ(I) 254 ENDDO 255 ENDIF 256 257 ! 258 ELSE 259 ! READ record from ASCII profile file using the default format OR 260 ! the format specified by the user; compute sigma_V from sigma_A 261 ! and wind speed 262 ! 263 DO WHILE ( JFLAG.EQ.0 ) ! 0 264 READ (MPUNIT,PROFRM,END=1000,ERR=98,IOSTAT=IOERRN) & 265 & KYEAR , KMONTH , KDAY , KHOUR , PFLHT(LEVEL) , & 266 & JFLAG , PFLWD(LEVEL) , PFLWS(LEVEL) , PFLTA(LEVEL) & 267 & , PFLSA(LEVEL) , PFLSW(LEVEL) 268 269 ! Convert the data to the required units 270 CALL PFLCNV(LEVEL) ! 0 271 272 ! Set the number of profile levels to current index, store 273 ! the 'top of profile' flag, and increment level if not at top 274 NPLVLS = LEVEL 275 ANPLVLS(IHR) = LEVEL 276 AIFLAG(IHR,LEVEL) = JFLAG 277 APFLHT(IHR,LEVEL) = PFLHT(LEVEL) 278 APFLWD(IHR,LEVEL) = PFLWD(LEVEL) 279 APFLWS(IHR,LEVEL) = PFLWS(LEVEL) 280 APFLTA(IHR,LEVEL) = PFLTA(LEVEL) 281 APFLSA(IHR,LEVEL) = PFLSA(LEVEL) 282 APFLSV(IHR,LEVEL) = PFLSV(LEVEL) 283 APFLSW(IHR,LEVEL) = PFLSW(LEVEL) 284 IF ( JFLAG.EQ.0 ) THEN 285 LEVEL = LEVEL + 1 ! 0 286 287 IF ( LEVEL.GT.MXPLVL ) THEN 288 IF ( .NOT.PFLERR ) THEN ! 0 289 ! WRITE Error Message: Number of profile levels 290 ! exceeds maximum allowable 291 WRITE (DUMMY,'(I8)') MXPLVL ! 0 292 CALL ERRHDL(PATH,MODNAM,'E','465',DUMMY) 293 PFLERR = .TRUE. 294 RUNERR = .TRUE. 295 ENDIF 296 297 ! Limit the number of levels to the maximum allowable 298 LEVEL = MXPLVL ! 0 299 ENDIF 300 301 ENDIF 302 303 ENDDO 304 305 ! Compute the vertical potential temperature gradient profile 306 IF ( .NOT.RUNERR ) THEN ! 0 307 NTGLVL = 0 ! 0 308 CALL COMPTG() 309 ANTGLVL(IHR) = NTGLVL 310 DO I = 1 , NTGLVL 311 APFLTG(IHR,I) = PFLTG(I) ! 0 312 APFLTGZ(IHR,I) = PFLTGZ(I) 313 ENDDO 314 ENDIF 315 316 317 ENDIF 318 319 ENDDO HOUR_LOOP 320 321 ! Set the date variables 322 CALL SET_DATES ! 0 323 324 GOTO 999 325 326 ! WRITE Error Messages: Error Reading Met Data File 327 328 98 CALL ERRHDL(PATH,MODNAM,'E','510','PROFFILE') ! 0 329 RUNERR = .TRUE. 330 GOTO 999 331 332 99 CALL ERRHDL(PATH,MODNAM,'E','510','SURFFILE') ! 0 333 RUNERR = .TRUE. 334 GOTO 999 335 336 1000 EOF = .TRUE. ! 0 337 ! Set the date variables 338 CALL SET_DATES 339 340 999 CONTINUE ! 0 341 END
HyperKWIC - Version 1.00DD executed at 20:00 on 1 Mar 2018 | Personal or Academic or Evaluation User | Free for Non-Commercial, Non-Government Use