fizhi-gridalt-hs/code/fizhi_init_veg.F

bac0c7c499 Andr*0001 #include "FIZHI_OPTIONS.h"
                0002
ac632d43e3 Jean*0003       SUBROUTINE FIZHI_INIT_VEG(myThid,vegdata,im,jm,nSx,nSy,Nxg,Nyg,
                0004      & maxtyp,nchp,nchptot,nchpland,lons,lats,surftype,tilefrac,
                0005      & igrd,ityp,chfr,chlt,chlon)
bac0c7c499 Andr*0006 C***********************************************************************
                0007 C Subroutine fizhi_init_veg - routine to read in the land surface types,
ac632d43e3 Jean*0008 C      interpolate to the models grid, and set up tile space for use by
                0009 C      the land surface model, the albedo calculation and the surface
bac0c7c499 Andr*0010 C      roughness calculation.
                0011 C
                0012 C INPUT:
ac632d43e3 Jean*0013 C
                0014 C myThid   - thread number (processor number)
bac0c7c499 Andr*0015 C vegdata  - Character*40 Vegetation Dataset name
                0016 C im       - longitude dimension
                0017 C jm       - latitude dimension (number of lat. points)
ac632d43e3 Jean*0018 C nSx      - Number of processors in x-direction
                0019 C nSy      - Number of processors in y-direction
bac0c7c499 Andr*0020 C maxtyp   - maximum allowable number of land surface types per grid box
                0021 C nchp     - integer per-processor number of tiles in tile space
                0022 C lons     - longitude in degrees [im,jm,nSx,nSy]
                0023 C lats     - latitude in degrees [im,jm,nSx,nSy]
                0024 C
                0025 C OUTPUT:
                0026 C
ac632d43e3 Jean*0027 C surftype - integer array of land surface types [im,jm,maxtyp,nSx,nSy]
                0028 C tilefrac - real array of corresponding land surface type fractions
                0029 C            [im,jm,maxtyp,nSx,nSy]
                0030 C igrd     - integer array in tile space of grid point number for each
                0031 C            tile [nchp,nSx,nSy]
                0032 C ityp     - integer array in tile space of land surface type for each
                0033 C            tile [nchp,nSx,nSy]
                0034 C chfr     - real array in tile space of land surface type fraction for
                0035 C            each tile [nchp,nSx,nSy]
bac0c7c499 Andr*0036 C
                0037 C NOTES:
                0038 C       Vegetation type as follows:
                0039 C                  1:  BROADLEAF EVERGREEN TREES
                0040 C                  2:  BROADLEAF DECIDUOUS TREES
                0041 C                  3:  NEEDLELEAF TREES
                0042 C                  4:  GROUND COVER
                0043 C                  5:  BROADLEAF SHRUBS
                0044 C                  6:  DWARF TREES (TUNDRA)
                0045 C                  7:  BARE SOIL
ac632d43e3 Jean*0046 C                  8:  DESERT
bac0c7c499 Andr*0047 C                  9:  GLACIER
                0048 C                 10:  DARK DESERT
                0049 C                100:  OCEAN
                0050 C***********************************************************************
ac632d43e3 Jean*0051       IMPLICIT NONE
bac0c7c499 Andr*0052 #include "EEPARAMS.h"
                0053
ac632d43e3 Jean*0054       INTEGER myThid,im,jm,maxtyp,nchp,nSx,nSy,Nxg,Nyg
                0055       INTEGER nchptot(nSx,nSy), nchpland(nSx,nSy)
                0056       INTEGER surftype(im,jm,maxtyp,nSx,nSy)
                0057       INTEGER igrd(nchp,nSx,nSy),ityp(nchp,nSx,nSy)
                0058       _RL tilefrac(im,jm,maxtyp,nSx,nSy)
bac0c7c499 Andr*0059       _RL lats(im,jm,nSx,nSy), lons(im,jm,nSx,nSy)
ac632d43e3 Jean*0060       _RL chfr(nchp,nSx,nSy),chlt(nchp,nSx,nSy),chlon(nchp,nSx,nSy)
bac0c7c499 Andr*0061
ac632d43e3 Jean*0062 C-    local variables:
                0063       CHARACTER*40 vegdata
                0064       INTEGER i,j,k,bi,bj
bac0c7c499 Andr*0065
ac632d43e3 Jean*0066       CHARACTER*15 aim_landfile
bac0c7c499 Andr*0067       _RS  aim_landFr(-1:34,-1:34,6,1)
ac632d43e3 Jean*0068       DATA aim_landfile /'landFrc.2f2.bin'/
                0069
                0070       WRITE(standardMessageUnit,'(2A)') ' FIZHI_INIT_VEG: ',
                0071      &   'defining surface type and fraction: ----------------------'
                0072
bac0c7c499 Andr*0073       CALL READ_REC_XY_RS(aim_LandFile,aim_landFr,1,0,myThid)
                0074
ac632d43e3 Jean*0075       DO bj = myByLo(myThid), myByHi(myThid)
                0076        DO bi = myBxLo(myThid), myBxHi(myThid)
bac0c7c499 Andr*0077
ac632d43e3 Jean*0078         do j = 1,jm
                0079          do i = 1,im
                0080           if(aim_landfr(i,j,bi,bj).gt.0.1) then
                0081            surftype(i,j,1,bi,bj) = 1
                0082            tilefrac(i,j,1,bi,bj) = 0.5 _d 0
                0083            surftype(i,j,2,bi,bj) = 2
                0084            tilefrac(i,j,2,bi,bj) = 0.5 _d 0
                0085           else
                0086            surftype(i,j,1,bi,bj) = 100
                0087            tilefrac(i,j,1,bi,bj) = 0.99 _d 0
                0088            surftype(i,j,2,bi,bj) = 100
                0089            tilefrac(i,j,2,bi,bj) = 0.01 _d 0
                0090           endif
                0091          enddo
                0092         enddo
                0093         do k = 3,maxtyp
                0094          do j = 1,jm
                0095           do i = 1,im
                0096            surftype(i,j,k,bi,bj) = 0
                0097            tilefrac(i,j,k,bi,bj) = 0.
                0098           enddo
                0099          enddo
                0100         enddo
                0101
                0102        ENDDO
bac0c7c499 Andr*0103       ENDDO
                0104
ac632d43e3 Jean*0105 C     create chip arrays for :
                0106 C      igrd :  grid index
                0107 C      ityp :  veg. type
                0108 C      chfr :  vegetation fraction
                0109 C      chlon:  chip longitude
                0110 C      chlt :  chip latitude
                0111
                0112 C     nchpland<=nchptot is the actual number of land chips
                0113       WRITE(standardMessageUnit,'(2A)') ' FIZHI_INIT_VEG: ',
                0114      &   'setting surface Tiles:'
                0115
                0116       DO bj = myByLo(myThid), myByHi(myThid)
                0117        DO bi = myBxLo(myThid), myBxHi(myThid)
                0118
                0119 C-       initialise grid index array:
                0120          do i=1,nchp
                0121            igrd(i,bi,bj) = 1
                0122          enddo
                0123
                0124 C-       land points:
                0125          nchpland(bi,bj) = 0
                0126          do k=1,maxtyp
                0127           do j=1,jm
                0128            do i=1,im
                0129              if(surftype(i,j,k,bi,bj).lt.100 .and.
                0130      &            tilefrac(i,j,k,bi,bj).gt.0.) then
                0131                nchpland(bi,bj)  = nchpland(bi,bj) + 1
                0132                igrd (nchpland(bi,bj),bi,bj) = i + (j-1)*im
                0133                ityp (nchpland(bi,bj),bi,bj) = surftype(i,j,k,bi,bj)
                0134                chfr (nchpland(bi,bj),bi,bj) = tilefrac(i,j,k,bi,bj)
                0135                chlon(nchpland(bi,bj),bi,bj) = lons(i,j,bi,bj)
                0136                chlt (nchpland(bi,bj),bi,bj) = lats(i,j,bi,bj)
                0137              endif
                0138            enddo
bac0c7c499 Andr*0139           enddo
ac632d43e3 Jean*0140          enddo
                0141
                0142 C-       ocean points:
                0143          nchptot(bi,bj) = nchpland(bi,bj)
                0144          do k=1,maxtyp
                0145           do j=1,jm
                0146            do i=1,im
                0147              if(surftype(i,j,k,bi,bj).ge.100 .and.
                0148      &            tilefrac(i,j,k,bi,bj).gt.0.) then
                0149                nchptot(bi,bj)  = nchptot(bi,bj) + 1
                0150                igrd (nchptot(bi,bj),bi,bj) = i + (j-1)*im
                0151                ityp (nchptot(bi,bj),bi,bj) = surftype(i,j,k,bi,bj)
                0152                chfr (nchptot(bi,bj),bi,bj) = tilefrac(i,j,k,bi,bj)
                0153                chlon(nchptot(bi,bj),bi,bj) = lons(i,j,bi,bj)
                0154                chlt (nchptot(bi,bj),bi,bj) = lats(i,j,bi,bj)
                0155              endif
                0156            enddo
bac0c7c499 Andr*0157           enddo
ac632d43e3 Jean*0158          enddo
                0159
                0160          WRITE(standardMessageUnit,'(2(A,I4),2(A,I10))') '  bi=', bi,
                0161      &    ', bj=', bj, ', # of Land Tiles=', nchpland(bi,bj),
                0162      &                 ', Total # of Tiles=', nchptot(bi,bj)
bac0c7c499 Andr*0163
ac632d43e3 Jean*0164        ENDDO
bac0c7c499 Andr*0165       ENDDO
                0166
ac632d43e3 Jean*0167       WRITE(standardMessageUnit,'(2A)') ' FIZHI_INIT_VEG: done'
                0168
bac0c7c499 Andr*0169       RETURN
                0170       END