Ignore:
Timestamp:
May 23, 2014, 11:48:41 AM (10 years ago)
Author:
igpis
Message:

version 9.2 beta. Changes from HH, AST, MC, NIK, IP. Changes in vert transform. New SPECIES input includes scavenging coefficients

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/src/readwind_nests.f90

    r4 r24  
    5151  integer :: igrib
    5252  integer :: gribVer,parCat,parNum,typSurf,valSurf,discipl
     53  integer :: parId !!added by mc for making it consistent with new readwind.f90
    5354  integer :: gotGrid
    5455  !HSO  end
     
    6970  integer :: isec1(56),isec2(22+nxmaxn+nymaxn)
    7071  real(kind=4) :: zsec4(jpunp)
     72  real(kind=4) :: xaux,yaux
    7173  real(kind=8) :: xauxin,yauxin
    72   real(kind=4) :: xaux,yaux,xaux0,yaux0
     74  real,parameter :: eps=1.e-4
    7375  real :: ewss(0:nxmaxn-1,0:nymaxn-1),nsss(0:nxmaxn-1,0:nymaxn-1)
    7476  real :: plev1,pmean,tv,fu,hlev1,ff10m,fflev1
     77  real :: conversion_factor !added by mc to make it consistent with new gridchek.f90
    7578
    7679  logical :: hflswitch,strswitch
     
    134137  endif
    135138
     139  conversion_factor=1.
     140
     141
    136142  else
    137143
     
    148154  call grib_get_int(igrib,'level',valSurf,iret)
    149155  call grib_check(iret,gribFunction,gribErrorMsg)
     156  call grib_get_int(igrib,'paramId',parId,iret) !added by mc to make it consisitent with new readwind.f90
     157  call grib_check(iret,gribFunction,gribErrorMsg) !added by mc to make it consisitent with new readwind.f90
    150158
    151159  !print*,discipl,parCat,parNum,typSurf,valSurf
     
    156164  isec1(8)=-1
    157165  isec1(8)=valSurf     ! level
     166   conversion_factor=1.
    158167  if ((parCat.eq.0).and.(parNum.eq.0).and.(typSurf.eq.105)) then ! T
    159168    isec1(6)=130         ! indicatorOfParameter
     
    166175  elseif ((parCat.eq.3).and.(parNum.eq.0).and.(typSurf.eq.1)) then !SP
    167176    isec1(6)=134         ! indicatorOfParameter
    168   elseif ((parCat.eq.2).and.(parNum.eq.32)) then ! W, actually eta dot
     177  elseif ((parCat.eq.2).and.(parNum.eq.32)) then ! W, actually eta dot !
    169178    isec1(6)=135         ! indicatorOfParameter
     179  elseif ((parCat.eq.128).and.(parNum.eq.77)) then ! W, actually eta dot !added by mc to make it consisitent with new readwind.f90
     180    isec1(6)=135         ! indicatorOfParameter    !added by mc to make it consisitent with new readwind.f90
    170181  elseif ((parCat.eq.3).and.(parNum.eq.0).and.(typSurf.eq.101)) then !SLP
    171182    isec1(6)=151         ! indicatorOfParameter
     
    180191  elseif ((parCat.eq.1).and.(parNum.eq.11).and.(typSurf.eq.1)) then ! SD
    181192    isec1(6)=141         ! indicatorOfParameter
    182   elseif ((parCat.eq.6).and.(parNum.eq.1)) then ! CC
     193    conversion_factor=1000. !added by mc to make it consisitent with new readwind.f90
     194  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
    183195    isec1(6)=164         ! indicatorOfParameter
    184   elseif ((parCat.eq.1).and.(parNum.eq.9)) then ! LSP
     196  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
    185197    isec1(6)=142         ! indicatorOfParameter
    186198  elseif ((parCat.eq.1).and.(parNum.eq.10)) then ! CP
    187199    isec1(6)=143         ! indicatorOfParameter
     200    conversion_factor=1000. !added by mc to make it consisitent with new readwind.f90
    188201  elseif ((parCat.eq.0).and.(parNum.eq.11).and.(typSurf.eq.1)) then ! SHF
    189202    isec1(6)=146         ! indicatorOfParameter
    190203  elseif ((parCat.eq.4).and.(parNum.eq.9).and.(typSurf.eq.1)) then ! SR
    191204    isec1(6)=176         ! indicatorOfParameter
    192   elseif ((parCat.eq.2).and.(parNum.eq.17)) then ! EWSS
     205  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
    193206    isec1(6)=180         ! indicatorOfParameter
    194   elseif ((parCat.eq.2).and.(parNum.eq.18)) then ! NSSS
     207  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
    195208    isec1(6)=181         ! indicatorOfParameter
    196209  elseif ((parCat.eq.3).and.(parNum.eq.4)) then ! ORO
    197210    isec1(6)=129         ! indicatorOfParameter
    198   elseif ((parCat.eq.3).and.(parNum.eq.7)) then ! SDO
     211   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
    199212    isec1(6)=160         ! indicatorOfParameter
    200213  elseif ((discipl.eq.2).and.(parCat.eq.0).and.(parNum.eq.0).and. &
     
    202215    isec1(6)=172         ! indicatorOfParameter
    203216  else
    204     print*,'***ERROR: undefined GRiB2 message found!',discipl, &
     217    print*,'***WARNING: undefined GRiB2 message found!',discipl, &
    205218         parCat,parNum,typSurf
     219  endif
     220  if(parId .ne. isec1(6) .and. parId .ne. 77) then !added by mc to make it consisitent with new readwind.f90
     221    write(*,*) 'parId',parId, 'isec1(6)',isec1(6)  !
     222!    stop
    206223  endif
    207224
     
    227244  ! CHECK GRID SPECIFICATIONS
    228245  if(isec2(2).ne.nxn(l)) stop &
    229        'READWIND: NX NOT CONSISTENT FOR A NESTING LEVEL'
     246  'READWIND: NX NOT CONSISTENT FOR A NESTING LEVEL'
    230247  if(isec2(3).ne.nyn(l)) stop &
    231        'READWIND: NY NOT CONSISTENT FOR A NESTING LEVEL'
     248  'READWIND: NY NOT CONSISTENT FOR A NESTING LEVEL'
    232249  if(isec2(12)/2-1.ne.nlev_ec) stop 'READWIND: VERTICAL DISCRET&
    233        &IZATION NOT CONSISTENT FOR A NESTING LEVEL'
     250  IZATION NOT CONSISTENT FOR A NESTING LEVEL'
    234251  endif ! ifield
    235252
    236253  !HSO  get the second part of the grid dimensions only from GRiB1 messages
    237   if ((gribVer.eq.1).and.(gotGrid.eq.0)) then
     254 if (isec1(6) .eq. 167 .and. (gotGrid.eq.0)) then ! !added by mc to make it consisitent with new readwind.f90
    238255    call grib_get_real8(igrib,'longitudeOfFirstGridPointInDegrees', &
    239256         xauxin,iret)
     
    242259         yauxin,iret)
    243260    call grib_check(iret,gribFunction,gribErrorMsg)
     261    if (xauxin.gt.180.) xauxin=xauxin-360.0
     262    if (xauxin.lt.-180.) xauxin=xauxin+360.0
     263
    244264    xaux=xauxin
    245265    yaux=yauxin
    246     xaux0=xlon0n(l)
    247     yaux0=ylat0n(l)
    248     if(xaux.lt.0.) xaux=xaux+360.
    249     if(yaux.lt.0.) yaux=yaux+360.
    250     if(xaux0.lt.0.) xaux0=xaux0+360.
    251     if(yaux0.lt.0.) yaux0=yaux0+360.
    252     if(xaux.ne.xaux0) &
    253          stop 'READWIND: LOWER LEFT LONGITUDE NOT CONSISTENT FOR A NES&
    254          &TING LEVEL'
    255     if(yaux.ne.yaux0) &
    256          stop 'READWIND: LOWER LEFT LATITUDE NOT CONSISTENT FOR A NEST&
    257          &ING LEVEL'
     266    if (abs(xaux-xlon0n(l)).gt.eps) &
     267    stop 'READWIND: LOWER LEFT LONGITUDE NOT CONSISTENT FOR A NESTING LEVEL'
     268    if (abs(yaux-ylat0n(l)).gt.eps) &
     269    stop 'READWIND: LOWER LEFT LATITUDE NOT CONSISTENT FOR A NESTING LEVEL'
    258270    gotGrid=1
    259271  endif
     
    280292             zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
    281293        if(isec1(6).eq.141) sdn(i,j,1,n,l)= &!! SNOW DEPTH
    282              zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
     294             zsec4(nxn(l)*(nyn(l)-j-1)+i+1)/conversion_factor !added by mc to make it consisitent with new readwind.f90!
    283295        if(isec1(6).eq.151) msln(i,j,1,n,l)= &!! SEA LEVEL PRESS.
    284296             zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
     
    298310        endif
    299311        if(isec1(6).eq.143) then                         !! CONVECTIVE PREC.
    300           convprecn(i,j,1,n,l)=zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
     312          convprecn(i,j,1,n,l)=zsec4(nxn(l)*(nyn(l)-j-1)+i+1)/conversion_factor !added by mc to make it consisitent with new readwind.f90
    301313          if (convprecn(i,j,1,n,l).lt.0.) convprecn(i,j,1,n,l)=0.
    302314        endif
Note: See TracChangeset for help on using the changeset viewer.
hosted by ZAMG