/[MITgcm]/MITgcm/verification/fizhi-gridalt-hs/code/fizhi_init_veg.F
ViewVC logotype

Diff of /MITgcm/verification/fizhi-gridalt-hs/code/fizhi_init_veg.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph | View Patch Patch

revision 1.1 by molod, Tue Aug 24 19:33:15 2004 UTC revision 1.2 by jmc, Thu Mar 22 14:25:21 2012 UTC
# Line 3  C $Name$ Line 3  C $Name$
3    
4  #include "FIZHI_OPTIONS.h"  #include "FIZHI_OPTIONS.h"
5    
6        subroutine fizhi_init_veg(mythid,vegdata,im,jm,Nsx,Nsy,Nxg,Nyg,        SUBROUTINE FIZHI_INIT_VEG(myThid,vegdata,im,jm,nSx,nSy,Nxg,Nyg,
7       . maxtyp,nchp,nchptot,nchpland,lons,lats,surftype,tilefrac,igrd,       & maxtyp,nchp,nchptot,nchpland,lons,lats,surftype,tilefrac,
8       . ityp,chfr,chlt,chlon)       & igrd,ityp,chfr,chlt,chlon)
9  C***********************************************************************  C***********************************************************************
10  C Subroutine fizhi_init_veg - routine to read in the land surface types,  C Subroutine fizhi_init_veg - routine to read in the land surface types,
11  C      interpolate to the models grid, and set up tile space for use by  C      interpolate to the models grid, and set up tile space for use by
12  C      the land surface model, the albedo calculation and the surface  C      the land surface model, the albedo calculation and the surface
13  C      roughness calculation.  C      roughness calculation.
14  C  C
15  C INPUT:  C INPUT:
16  C  C
17  C mythid   - thread number (processor number)  C myThid   - thread number (processor number)
18  C vegdata  - Character*40 Vegetation Dataset name  C vegdata  - Character*40 Vegetation Dataset name
19  C im       - longitude dimension  C im       - longitude dimension
20  C jm       - latitude dimension (number of lat. points)  C jm       - latitude dimension (number of lat. points)
21  C Nsx      - Number of processors in x-direction  C nSx      - Number of processors in x-direction
22  C Nsy      - Number of processors in y-direction  C nSy      - Number of processors in y-direction
23  C maxtyp   - maximum allowable number of land surface types per grid box  C maxtyp   - maximum allowable number of land surface types per grid box
24  C nchp     - integer per-processor number of tiles in tile space  C nchp     - integer per-processor number of tiles in tile space
25  C lons     - longitude in degrees [im,jm,nSx,nSy]  C lons     - longitude in degrees [im,jm,nSx,nSy]
# Line 27  C lats     - latitude in degrees [im,jm, Line 27  C lats     - latitude in degrees [im,jm,
27  C  C
28  C OUTPUT:  C OUTPUT:
29  C  C
30  C surftype - integer array of land surface types [im,jm,maxtyp,Nsx,Nsy]  C surftype - integer array of land surface types [im,jm,maxtyp,nSx,nSy]
31  C tilefrac - real array of corresponding land surface type fractions  C tilefrac - real array of corresponding land surface type fractions
32  C            [im,jm,maxtyp,Nsx,Nsy]  C            [im,jm,maxtyp,nSx,nSy]
33  C igrd     - integer array in tile space of grid point number for each  C igrd     - integer array in tile space of grid point number for each
34  C            tile [nchp,Nsx,Nsy]  C            tile [nchp,nSx,nSy]
35  C ityp     - integer array in tile space of land surface type for each  C ityp     - integer array in tile space of land surface type for each
36  C            tile [nchp,Nsx,Nsy]  C            tile [nchp,nSx,nSy]
37  C chfr     - real array in tile space of land surface type fraction for  C chfr     - real array in tile space of land surface type fraction for
38  C            each tile [nchp,Nsx,Nsy]  C            each tile [nchp,nSx,nSy]
39  C  C
40  C NOTES:  C NOTES:
41  C       Vegetation type as follows:  C       Vegetation type as follows:
# Line 46  C                  4:  GROUND COVER Line 46  C                  4:  GROUND COVER
46  C                  5:  BROADLEAF SHRUBS  C                  5:  BROADLEAF SHRUBS
47  C                  6:  DWARF TREES (TUNDRA)  C                  6:  DWARF TREES (TUNDRA)
48  C                  7:  BARE SOIL  C                  7:  BARE SOIL
49  C                  8:  DESERT      C                  8:  DESERT
50  C                  9:  GLACIER  C                  9:  GLACIER
51  C                 10:  DARK DESERT  C                 10:  DARK DESERT
52  C                100:  OCEAN  C                100:  OCEAN
53  C***********************************************************************  C***********************************************************************
54        implicit none        IMPLICIT NONE
55  #include "EEPARAMS.h"  #include "EEPARAMS.h"
56    
57        integer mythid,im,jm,maxtyp,nchp,Nsx,Nsy,Nxg,Nyg        INTEGER myThid,im,jm,maxtyp,nchp,nSx,nSy,Nxg,Nyg
58        integer nchptot(Nsx,Nsy),nchpland(Nsx,Nsy)        INTEGER nchptot(nSx,nSy), nchpland(nSx,nSy)
59        integer surftype(im,jm,maxtyp,Nsx,Nsy)        INTEGER surftype(im,jm,maxtyp,nSx,nSy)
60        integer igrd(nchp,Nsx,Nsy),ityp(nchp,Nsx,Nsy)        INTEGER igrd(nchp,nSx,nSy),ityp(nchp,nSx,nSy)
61        _RL tilefrac(im,jm,maxtyp,Nsx,Nsy)        _RL tilefrac(im,jm,maxtyp,nSx,nSy)
62        _RL lats(im,jm,nSx,nSy), lons(im,jm,nSx,nSy)        _RL lats(im,jm,nSx,nSy), lons(im,jm,nSx,nSy)
63        _RL chfr(nchp,Nsx,Nsy),chlt(nchp,Nsx,Nsy),chlon(nchp,Nsx,Nsy)        _RL chfr(nchp,nSx,nSy),chlt(nchp,nSx,nSy),chlon(nchp,nSx,nSy)
       character*40 vegdata  
64    
65        integer i,j,k,bi,bj  C-    local variables:
66          CHARACTER*40 vegdata
67          INTEGER i,j,k,bi,bj
68    
69        character *15 aim_landfile        CHARACTER*15 aim_landfile
70        _RS  aim_landFr(-1:34,-1:34,6,1)        _RS  aim_landFr(-1:34,-1:34,6,1)
71        data aim_landfile /'landFrc.2f2.bin'/        DATA aim_landfile /'landFrc.2f2.bin'/
72    
73          WRITE(standardMessageUnit,'(2A)') ' FIZHI_INIT_VEG: ',
74         &   'defining surface type and fraction: ----------------------'
75    
76        CALL READ_REC_XY_RS(aim_LandFile,aim_landFr,1,0,myThid)        CALL READ_REC_XY_RS(aim_LandFile,aim_landFr,1,0,myThid)
77    
78        DO BJ = myByLo(myThid), myByHi(myThid)        DO bj = myByLo(myThid), myByHi(myThid)
79        DO BI = myBxLo(myThid), myBxHi(myThid)         DO bi = myBxLo(myThid), myBxHi(myThid)
80    
81        do j = 1,jm          do j = 1,jm
82        do i = 1,im           do i = 1,im
83         if(aim_landfr(i,j,bi,bj).gt.0.1) then            if(aim_landfr(i,j,bi,bj).gt.0.1) then
84          surftype(i,j,1,bi,bj) = 1             surftype(i,j,1,bi,bj) = 1
85          tilefrac(i,j,1,bi,bj) = 0.5             tilefrac(i,j,1,bi,bj) = 0.5 _d 0
86          surftype(i,j,2,bi,bj) = 2             surftype(i,j,2,bi,bj) = 2
87          tilefrac(i,j,2,bi,bj) = 0.5             tilefrac(i,j,2,bi,bj) = 0.5 _d 0
88         else            else
89          surftype(i,j,1,bi,bj) = 100             surftype(i,j,1,bi,bj) = 100
90          tilefrac(i,j,1,bi,bj) = 0.99             tilefrac(i,j,1,bi,bj) = 0.99 _d 0
91          surftype(i,j,2,bi,bj) = 100             surftype(i,j,2,bi,bj) = 100
92          tilefrac(i,j,2,bi,bj) = 0.01             tilefrac(i,j,2,bi,bj) = 0.01 _d 0
93         endif            endif
94        enddo           enddo
95        enddo          enddo
96        do k = 3,maxtyp          do k = 3,maxtyp
97        do j = 1,jm           do j = 1,jm
98        do i = 1,im            do i = 1,im
99         surftype(i,j,k,bi,bj) = 0             surftype(i,j,k,bi,bj) = 0
100         tilefrac(i,j,k,bi,bj) = 0.             tilefrac(i,j,k,bi,bj) = 0.
101        enddo            enddo
102        enddo           enddo
103        enddo          enddo
104    
105        ENDDO         ENDDO
106        ENDDO        ENDDO
107    
108  c     create chip arrays for :  C     create chip arrays for :
109  c      igrd :  grid index  C      igrd :  grid index
110  c      ityp :  veg. type  C      ityp :  veg. type
111  c      chfr :  vegetation fraction  C      chfr :  vegetation fraction
112  c      chlon:  chip longitude  C      chlon:  chip longitude
113  c      chlt :  chip latitude  C      chlt :  chip latitude
114    
115  c     nchpland<=nchptot is the actual number of land chips  C     nchpland<=nchptot is the actual number of land chips
116          WRITE(standardMessageUnit,'(2A)') ' FIZHI_INIT_VEG: ',
117        DO BJ = myByLo(myThid), myByHi(myThid)       &   'setting surface Tiles:'
118          DO BI = myBxLo(myThid), myBxHi(myThid)  
119          DO bj = myByLo(myThid), myByHi(myThid)
120  c         land points         DO bi = myBxLo(myThid), myBxHi(myThid)
121  c         -----------  
122            nchpland(bi,bj) = 0  C-       initialise grid index array:
123            do k=1,maxtyp           do i=1,nchp
124              do j=1,jm             igrd(i,bi,bj) = 1
125                do i=1,im           enddo
126                  if(surftype(i,j,k,bi,bj).lt.100 .and.  
127       .               tilefrac(i,j,k,bi,bj).gt.0.) then  C-       land points:
128                    nchpland(bi,bj)  = nchpland(bi,bj) + 1           nchpland(bi,bj) = 0
129                    igrd (nchpland(bi,bj),bi,bj) = i + (j-1)*im           do k=1,maxtyp
130                    ityp (nchpland(bi,bj),bi,bj) = surftype(i,j,k,bi,bj)            do j=1,jm
131                    chfr (nchpland(bi,bj),bi,bj) = tilefrac(i,j,k,bi,bj)             do i=1,im
132                    chlon(nchpland(bi,bj),bi,bj) = lons(i,j,bi,bj)               if(surftype(i,j,k,bi,bj).lt.100 .and.
133                    chlt (nchpland(bi,bj),bi,bj) = lats(i,j,bi,bj)       &            tilefrac(i,j,k,bi,bj).gt.0.) then
134                  endif                 nchpland(bi,bj)  = nchpland(bi,bj) + 1
135                enddo                 igrd (nchpland(bi,bj),bi,bj) = i + (j-1)*im
136              enddo                 ityp (nchpland(bi,bj),bi,bj) = surftype(i,j,k,bi,bj)
137                   chfr (nchpland(bi,bj),bi,bj) = tilefrac(i,j,k,bi,bj)
138                   chlon(nchpland(bi,bj),bi,bj) = lons(i,j,bi,bj)
139                   chlt (nchpland(bi,bj),bi,bj) = lats(i,j,bi,bj)
140                 endif
141               enddo
142            enddo            enddo
143                       enddo
144  c         ocean points  
145  c         ------------  C-       ocean points:
146            nchptot(bi,bj) = nchpland(bi,bj)           nchptot(bi,bj) = nchpland(bi,bj)
147                       do k=1,maxtyp
148            do k=1,maxtyp            do j=1,jm
149              do j=1,jm             do i=1,im
150                do i=1,im               if(surftype(i,j,k,bi,bj).ge.100 .and.
151                  if(surftype(i,j,k,bi,bj).ge.100 .and.       &            tilefrac(i,j,k,bi,bj).gt.0.) then
152       .               tilefrac(i,j,k,bi,bj).gt.0.) then                 nchptot(bi,bj)  = nchptot(bi,bj) + 1
153                    nchptot(bi,bj)  = nchptot(bi,bj) + 1                 igrd (nchptot(bi,bj),bi,bj) = i + (j-1)*im
154                    igrd (nchptot(bi,bj),bi,bj) = i + (j-1)*im                 ityp (nchptot(bi,bj),bi,bj) = surftype(i,j,k,bi,bj)
155                    ityp (nchptot(bi,bj),bi,bj) = surftype(i,j,k,bi,bj)                 chfr (nchptot(bi,bj),bi,bj) = tilefrac(i,j,k,bi,bj)
156                    chfr (nchptot(bi,bj),bi,bj) = tilefrac(i,j,k,bi,bj)                 chlon(nchptot(bi,bj),bi,bj) = lons(i,j,bi,bj)
157                    chlon(nchptot(bi,bj),bi,bj) = lons(i,j,bi,bj)                 chlt (nchptot(bi,bj),bi,bj) = lats(i,j,bi,bj)
158                    chlt (nchptot(bi,bj),bi,bj) = lats(i,j,bi,bj)               endif
159                  endif             enddo
               enddo  
             enddo  
160            enddo            enddo
161                       enddo
          print *,'Number of Total Tiles for bi=',bi,': ',nchptot(bi,bj)  
          print *,'Number of Land  Tiles for bi=',bi,': ',nchpland(bi,bj)  
162    
163          ENDDO           WRITE(standardMessageUnit,'(2(A,I4),2(A,I10))') '  bi=', bi,
164         &    ', bj=', bj, ', # of Land Tiles=', nchpland(bi,bj),
165         &                 ', Total # of Tiles=', nchptot(bi,bj)
166    
167           ENDDO
168        ENDDO        ENDDO
169    
170                    WRITE(standardMessageUnit,'(2A)') ' FIZHI_INIT_VEG: done'
171    
172        RETURN        RETURN
173        END        END
Legend:
Removed from v.1.1  
changed lines
  Added in v.1.2
  ViewVC Help
Powered by ViewVC 1.1.22