Changeset 3c1da52 in flexpart.git for src/readwind_nests.f90
- Timestamp:
- Jul 15, 2015, 9:45:31 AM (9 years ago)
- Branches:
- svn-petra
- Children:
- c07b09e
- Parents:
- 31674b5
- git-author:
- Ignacio Pisso <Ignacio.Pisso@…> (07/15/15 09:23:56)
- git-committer:
- Ignacio Pisso <Ignacio.Pisso@…> (07/15/15 09:45:31)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
src/readwind_nests.f90
r31674b5 r3c1da52 34 34 ! * 35 35 !***************************************************************************** 36 ! Changes, Bernd C. Krueger, Feb. 2001: * 36 ! CHANGES: 37 ! 38 ! 2/2001, Bernc C. Krueger: * 37 39 ! Variables tthn and qvhn (on eta coordinates) in common block * 38 ! CHANGE: 11/01/2008, Harald Sodemann, GRIB1/2 input with ECMWF grib_api * 39 ! CHANGE: 03/12/2008, Harald Sodemann, update to f90 with ECMWF grib_api * 40 ! 2/2015, PS: add missing paramter iret in call to grib subr 40 ! 11/01/2008, Harald Sodemann, GRIB1/2 input with ECMWF grib_api * 41 ! 03/12/2008, Harald Sodemann, update to f90 with ECMWF grib_api * 42 ! 2/2015, PS: add missing paramter iret in call to grib subr 43 ! 3/2015, PS: add verbosity output 44 ! 3/2015, PS: add GRIB2 coding for deta/dt dp/deta used by ZAMG 41 45 ! 42 46 !***************************************************************************** … … 84 88 character(len=24) :: gribErrorMsg = 'Error reading grib file' 85 89 character(len=20) :: gribFunction = 'readwind_nests' 90 91 call verboutput(verbosity,1,' entered readwind_nests') 86 92 87 93 do l=1,numbnests … … 100 106 ! 101 107 102 5 call grib_open_file(ifile,path(numpath+2*(l-1)+1) & 108 5 continue 109 call verboutput(verbosity,1,' call grib_open_file '//& 110 path(numpath+2*(l-1)+1)(1:length(numpath+2*(l-1)+1))//trim(wfnamen(l,indj))) 111 call grib_open_file(ifile,path(numpath+2*(l-1)+1) & 103 112 (1:length(numpath+2*(l-1)+1))//trim(wfnamen(l,indj)),'r') 104 113 if (iret.ne.GRIB_SUCCESS) then … … 127 136 if (gribVer.eq.1) then ! GRIB Edition 1 128 137 129 !print*,'GRiB Edition 1' 130 !read the grib2 identifiers 131 call grib_get_int(igrib,'indicatorOfParameter',isec1(6),iret) 132 call grib_check(iret,gribFunction,gribErrorMsg) 133 call grib_get_int(igrib,'level',isec1(8),iret) 134 call grib_check(iret,gribFunction,gribErrorMsg) 135 136 !change code for etadot to code for omega 137 if (isec1(6).eq.77) then 138 isec1(6)=135 139 endif 140 141 conversion_factor=1. 142 138 !print*,'GRiB Edition 1' 139 140 ! read the grib2 identifiers 141 142 call grib_get_int(igrib,'indicatorOfParameter',isec1(6),iret) 143 call grib_check(iret,gribFunction,gribErrorMsg) 144 call grib_get_int(igrib,'level',isec1(8),iret) 145 call grib_check(iret,gribFunction,gribErrorMsg) 146 147 !change code for etadot to code for omega 148 if (isec1(6).eq.77) then 149 isec1(6)=135 150 endif 151 152 conversion_factor=1. 143 153 144 154 else 145 155 146 !print*,'GRiB Edition 2' 147 !read the grib2 identifiers 148 call grib_get_int(igrib,'discipline',discipl,iret) 149 call grib_check(iret,gribFunction,gribErrorMsg) 150 call grib_get_int(igrib,'parameterCategory',parCat,iret) 151 call grib_check(iret,gribFunction,gribErrorMsg) 152 call grib_get_int(igrib,'parameterNumber',parNum,iret) 153 call grib_check(iret,gribFunction,gribErrorMsg) 154 call grib_get_int(igrib,'typeOfFirstFixedSurface',typSurf,iret) 155 call grib_check(iret,gribFunction,gribErrorMsg) 156 call grib_get_int(igrib,'level',valSurf,iret) 157 call grib_check(iret,gribFunction,gribErrorMsg) 158 call grib_get_int(igrib,'paramId',parId,iret) !added by mc to make it consisitent with new readwind.f90 159 call grib_check(iret,gribFunction,gribErrorMsg) !added by mc to make it consisitent with new readwind.f90 160 161 !print*,discipl,parCat,parNum,typSurf,valSurf 162 163 !convert to grib1 identifiers 164 isec1(6)=-1 165 isec1(7)=-1 166 isec1(8)=-1 167 isec1(8)=valSurf ! level 168 conversion_factor=1. 169 if ((parCat.eq.0).and.(parNum.eq.0).and.(typSurf.eq.105)) then ! T 170 isec1(6)=130 ! indicatorOfParameter 171 elseif ((parCat.eq.2).and.(parNum.eq.2).and.(typSurf.eq.105)) then ! U 172 isec1(6)=131 ! indicatorOfParameter 173 elseif ((parCat.eq.2).and.(parNum.eq.3).and.(typSurf.eq.105)) then ! V 174 isec1(6)=132 ! indicatorOfParameter 175 elseif ((parCat.eq.1).and.(parNum.eq.0).and.(typSurf.eq.105)) then ! Q 176 isec1(6)=133 ! indicatorOfParameter 177 elseif ((parCat.eq.3).and.(parNum.eq.0).and.(typSurf.eq.1)) then !SP 178 isec1(6)=134 ! indicatorOfParameter 179 elseif ((parCat.eq.2).and.(parNum.eq.32)) then ! W, actually eta dot ! 180 isec1(6)=135 ! indicatorOfParameter 181 elseif ((parCat.eq.128).and.(parNum.eq.77)) then ! W, actually eta dot !added by mc to make it consisitent with new readwind.f90 182 isec1(6)=135 ! indicatorOfParameter !added by mc to make it consisitent with new readwind.f90 183 elseif ((parCat.eq.3).and.(parNum.eq.0).and.(typSurf.eq.101)) then !SLP 184 isec1(6)=151 ! indicatorOfParameter 185 elseif ((parCat.eq.2).and.(parNum.eq.2).and.(typSurf.eq.103)) then ! 10U 186 isec1(6)=165 ! indicatorOfParameter 187 elseif ((parCat.eq.2).and.(parNum.eq.3).and.(typSurf.eq.103)) then ! 10V 188 isec1(6)=166 ! indicatorOfParameter 189 elseif ((parCat.eq.0).and.(parNum.eq.0).and.(typSurf.eq.103)) then ! 2T 190 isec1(6)=167 ! indicatorOfParameter 191 elseif ((parCat.eq.0).and.(parNum.eq.6).and.(typSurf.eq.103)) then ! 2D 192 isec1(6)=168 ! indicatorOfParameter 193 elseif ((parCat.eq.1).and.(parNum.eq.11).and.(typSurf.eq.1)) then ! SD 194 isec1(6)=141 ! indicatorOfParameter 195 conversion_factor=1000. !added by mc to make it consisitent with new readwind.f90 196 elseif ((parCat.eq.6).and.(parNum.eq.1) .or. parId .eq. 164) then ! CC !added by mc to make it consisitent with new readwind.f90 197 isec1(6)=164 ! indicatorOfParameter 198 elseif ((parCat.eq.1).and.(parNum.eq.9) .or. parId .eq. 142) then ! LSP !added by mc to make it consisitent with new readwind.f90 199 isec1(6)=142 ! indicatorOfParameter 200 elseif ((parCat.eq.1).and.(parNum.eq.10)) then ! CP 201 isec1(6)=143 ! indicatorOfParameter 202 conversion_factor=1000. !added by mc to make it consisitent with new readwind.f90 203 elseif ((parCat.eq.0).and.(parNum.eq.11).and.(typSurf.eq.1)) then ! SHF 204 isec1(6)=146 ! indicatorOfParameter 205 elseif ((parCat.eq.4).and.(parNum.eq.9).and.(typSurf.eq.1)) then ! SR 206 isec1(6)=176 ! indicatorOfParameter 207 elseif ((parCat.eq.2).and.(parNum.eq.17) .or. parId .eq. 180) then ! EWSS !added by mc to make it consisitent with new readwind.f90 208 isec1(6)=180 ! indicatorOfParameter 209 elseif ((parCat.eq.2).and.(parNum.eq.18) .or. parId .eq. 181) then ! NSSS !added by mc to make it consisitent with new readwind.f90 210 isec1(6)=181 ! indicatorOfParameter 211 elseif ((parCat.eq.3).and.(parNum.eq.4)) then ! ORO 212 isec1(6)=129 ! indicatorOfParameter 213 elseif ((parCat.eq.3).and.(parNum.eq.7) .or. parId .eq. 160) then ! SDO !added by mc to make it consisitent with new readwind.f90 214 isec1(6)=160 ! indicatorOfParameter 215 elseif ((discipl.eq.2).and.(parCat.eq.0).and.(parNum.eq.0).and. & 216 (typSurf.eq.1)) then ! LSM 217 isec1(6)=172 ! indicatorOfParameter 218 else 219 print*,'***WARNING: undefined GRiB2 message found!',discipl, & 220 parCat,parNum,typSurf 221 endif 222 if(parId .ne. isec1(6) .and. parId .ne. 77) then !added by mc to make it consisitent with new readwind.f90 223 write(*,*) 'parId',parId, 'isec1(6)',isec1(6) ! 224 ! stop 225 endif 156 !print*,'GRiB Edition 2' 157 158 ! read the grib2 identifiers 159 160 call grib_get_int(igrib,'discipline',discipl,iret) 161 call grib_check(iret,gribFunction,gribErrorMsg) 162 call grib_get_int(igrib,'parameterCategory',parCat,iret) 163 call grib_check(iret,gribFunction,gribErrorMsg) 164 call grib_get_int(igrib,'parameterNumber',parNum,iret) 165 call grib_check(iret,gribFunction,gribErrorMsg) 166 call grib_get_int(igrib,'typeOfFirstFixedSurface',typSurf,iret) 167 call grib_check(iret,gribFunction,gribErrorMsg) 168 call grib_get_int(igrib,'level',valSurf,iret) 169 call grib_check(iret,gribFunction,gribErrorMsg) 170 call grib_get_int(igrib,'paramId',parId,iret) !added by mc to make it consistent with new readwind.f90 171 call grib_check(iret,gribFunction,gribErrorMsg) !added by mc to make it consistent with new readwind.f90 172 173 !print*,discipl,parCat,parNum,typSurf,valSurf 174 175 !convert to grib1 identifiers 176 isec1(6)=-1 177 isec1(7)=-1 178 isec1(8)=-1 179 isec1(8)=valSurf ! level 180 conversion_factor=1. 181 if ((parCat.eq.0).and.(parNum.eq.0).and.(typSurf.eq.105)) then ! T 182 isec1(6)=130 ! indicatorOfParameter 183 elseif ((parCat.eq.2).and.(parNum.eq.2).and.(typSurf.eq.105)) then ! U 184 isec1(6)=131 ! indicatorOfParameter 185 elseif ((parCat.eq.2).and.(parNum.eq.3).and.(typSurf.eq.105)) then ! V 186 isec1(6)=132 ! indicatorOfParameter 187 elseif ((parCat.eq.1).and.(parNum.eq.0).and.(typSurf.eq.105)) then ! Q 188 isec1(6)=133 ! indicatorOfParameter 189 elseif ((parCat.eq.3).and.(parNum.eq.0).and.(typSurf.eq.1)) then !SP 190 isec1(6)=134 ! indicatorOfParameter 191 elseif ((parCat.eq.2).and.(parNum.eq.32)) then ! W, actually detadtdpdeta 192 isec1(6)=135 ! indicatorOfParameter 193 elseif ((parCat.eq.128).and.(parNum.eq.77)) then ! W, actually detadtdpdeta 194 isec1(6)=135 ! indicatorOfParameter !added by mc to make it consistent with new readwind.f90 195 elseif ((parCat.eq.2).and.(parNum.eq.8)) then ! Omega, actually detadtdpdeta 196 isec1(6)=135 ! indicatorOfParameter !added by mc to make it consistent with new readwind.f90 197 elseif ((parCat.eq.3).and.(parNum.eq.0).and.(typSurf.eq.101)) then !SLP 198 isec1(6)=151 ! indicatorOfParameter 199 elseif ((parCat.eq.2).and.(parNum.eq.2).and.(typSurf.eq.103)) then ! 10U 200 isec1(6)=165 ! indicatorOfParameter 201 elseif ((parCat.eq.2).and.(parNum.eq.3).and.(typSurf.eq.103)) then ! 10V 202 isec1(6)=166 ! indicatorOfParameter 203 elseif ((parCat.eq.0).and.(parNum.eq.0).and.(typSurf.eq.103)) then ! 2T 204 isec1(6)=167 ! indicatorOfParameter 205 elseif ((parCat.eq.0).and.(parNum.eq.6).and.(typSurf.eq.103)) then ! 2D 206 isec1(6)=168 ! indicatorOfParameter 207 elseif ((parCat.eq.1).and.(parNum.eq.11).and.(typSurf.eq.1)) then ! SD 208 isec1(6)=141 ! indicatorOfParameter 209 conversion_factor=1000. !added by mc to make it consistent with new readwind.f90 210 elseif ((parCat.eq.6).and.(parNum.eq.1) .or. parId .eq. 164) then ! CC !added by mc to make it consistent with new readwind.f90 211 isec1(6)=164 ! indicatorOfParameter 212 elseif ((parCat.eq.1).and.(parNum.eq.9) .or. parId .eq. 142) then ! LSP !added by mc to make it consistent with new readwind.f90 213 isec1(6)=142 ! indicatorOfParameter 214 elseif ((parCat.eq.1).and.(parNum.eq.10)) then ! CP 215 isec1(6)=143 ! indicatorOfParameter 216 conversion_factor=1000. !added by mc to make it consistent with new readwind.f90 217 elseif ((parCat.eq.0).and.(parNum.eq.11).and.(typSurf.eq.1)) then ! SHF 218 isec1(6)=146 ! indicatorOfParameter 219 elseif ((parCat.eq.4).and.(parNum.eq.9).and.(typSurf.eq.1)) then ! SR 220 isec1(6)=176 ! indicatorOfParameter 221 elseif ((parCat.eq.2).and.(parNum.eq.17) .or. parId .eq. 180) then ! EWSS !added by mc to make it consistent with new readwind.f90 222 isec1(6)=180 ! indicatorOfParameter 223 elseif ((parCat.eq.2).and.(parNum.eq.18) .or. parId .eq. 181) then ! NSSS !added by mc to make it consistent with new readwind.f90 224 isec1(6)=181 ! indicatorOfParameter 225 elseif ((parCat.eq.3).and.(parNum.eq.4)) then ! ORO 226 isec1(6)=129 ! indicatorOfParameter 227 elseif ((parCat.eq.3).and.(parNum.eq.7) .or. parId .eq. 160) then ! SDO !added by mc to make it consistent with new readwind.f90 228 isec1(6)=160 ! indicatorOfParameter 229 elseif ((discipl.eq.2).and.(parCat.eq.0).and.(parNum.eq.0).and. & 230 (typSurf.eq.1)) then ! LSM 231 isec1(6)=172 ! indicatorOfParameter 232 else 233 print*,'***WARNING: undefined GRiB2 message found!',discipl, & 234 parCat,parNum,typSurf 235 endif 236 if(parId .ne. isec1(6) .and. parId .ne. 77) then !added by mc to make it consistent with new readwind.f90 237 write(*,*) 'parId',parId, 'isec1(6)',isec1(6) ! 238 ! stop 239 endif 226 240 227 241 endif … … 246 260 ! CHECK GRID SPECIFICATIONS 247 261 if(isec2(2).ne.nxn(l)) stop & 248 'READWIND: NX NOT CONSISTENT FOR A NESTING LEVEL'262 'READWIND: NX NOT CONSISTENT FOR A NESTING LEVEL' 249 263 if(isec2(3).ne.nyn(l)) stop & 250 'READWIND: NY NOT CONSISTENT FOR A NESTING LEVEL'251 if(isec2(12)/2-1.ne.nlev_ec) stop 'READWIND: VERTICAL DISCRET&252 IZATION NOT CONSISTENT FOR A NESTING LEVEL'264 'READWIND: NY NOT CONSISTENT FOR A NESTING LEVEL' 265 if(isec2(12)/2-1.ne.nlev_ec) stop & 266 'READWIND: VERTICAL DISCRETIZATION NOT CONSISTENT FOR A NESTING LEVEL' 253 267 endif ! ifield 254 268 255 269 !HSO get the second part of the grid dimensions only from GRiB1 messages 256 if (isec1(6) .eq. 167 .and. (gotGrid.eq.0)) then ! !added by mc to make it consis itent with new readwind.f90270 if (isec1(6) .eq. 167 .and. (gotGrid.eq.0)) then ! !added by mc to make it consistent with new readwind.f90 257 271 call grib_get_real8(igrib,'longitudeOfFirstGridPointInDegrees', & 258 272 xauxin,iret) … … 294 308 zsec4(nxn(l)*(nyn(l)-j-1)+i+1) 295 309 if(isec1(6).eq.141) sdn(i,j,1,n,l)= &!! SNOW DEPTH 296 zsec4(nxn(l)*(nyn(l)-j-1)+i+1)/conversion_factor !added by mc to make it consis itent with new readwind.f90!310 zsec4(nxn(l)*(nyn(l)-j-1)+i+1)/conversion_factor !added by mc to make it consistent with new readwind.f90! 297 311 if(isec1(6).eq.151) msln(i,j,1,n,l)= &!! SEA LEVEL PRESS. 298 312 zsec4(nxn(l)*(nyn(l)-j-1)+i+1) … … 312 326 endif 313 327 if(isec1(6).eq.143) then !! CONVECTIVE PREC. 314 convprecn(i,j,1,n,l)=zsec4(nxn(l)*(nyn(l)-j-1)+i+1)/conversion_factor !added by mc to make it consis itent with new readwind.f90328 convprecn(i,j,1,n,l)=zsec4(nxn(l)*(nyn(l)-j-1)+i+1)/conversion_factor !added by mc to make it consistent with new readwind.f90 315 329 if (convprecn(i,j,1,n,l).lt.0.) convprecn(i,j,1,n,l)=0. 316 330 endif … … 421 435 end do 422 436 437 call verboutput(verbosity,1,' leaving readwind_nests') 438 423 439 return 440 424 441 888 write(*,*) ' #### FLEXPART MODEL ERROR! WINDFIELD #### ' 425 442 write(*,*) ' #### ',wfnamen(l,indj),' FOR NESTING LEVEL #### ' … … 431 448 write(*,*) ' #### ',wfnamen(l,indj),' #### ' 432 449 write(*,*) ' #### CANNOT BE OPENED FOR NESTING LEVEL ',l,'####' 450 call verboutput(verbosity,1,' leaving readwind') 433 451 434 452 end subroutine readwind_nests
Note: See TracChangeset
for help on using the changeset viewer.