/[MITgcm]/MITgcm/pkg/fizhi/fizhi_init_veg.F
ViewVC logotype

Diff of /MITgcm/pkg/fizhi/fizhi_init_veg.F

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

revision 1.9 by molod, Wed Jun 9 20:33:37 2004 UTC revision 1.15 by molod, Fri Jul 23 22:32:28 2004 UTC
# Line 1  Line 1 
1  C $Header$  C $Header$
2  C $Name$  C $Name$
3    
4    #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,lons,lats,surftype,tilefrac,igrd,ityp,chfr,chlt,chlon)       . maxtyp,nchp,nchptot,nchpland,lons,lats,surftype,tilefrac,
8         . 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
# Line 51  C*************************************** Line 54  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,nchptot,nchpland,Nsx,Nsy,Nxg,Nyg
58        integer surftype(im,jm,maxtyp,Nsx,Nsy)        integer surftype(im,jm,maxtyp,Nsx,Nsy)
59        integer igrd(nchp,Nsx,Nsy),ityp(nchp,Nsx,Nsy)        integer igrd(nchp,Nsx,Nsy),ityp(nchp,Nsx,Nsy)
60        real tilefrac(im,jm,maxtyp,Nsx,Nsy)        _RL tilefrac(im,jm,maxtyp,Nsx,Nsy)
61        real lats(im,jm,nSx,nSy), lons(im,jm,nSx,nSy)        _RL lats(im,jm,nSx,nSy), lons(im,jm,nSx,nSy)
62        real 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)
63        character*40 vegdata        character*40 vegdata
64        integer imdata,jmdata,Nxgdata,Nygdata        integer imdata,jmdata,Nxgdata,Nygdata
65        integer nchplocal,nchpland,biglobal,bjglobal        integer biglobal,bjglobal
66    
67        integer*4 im_32, jm_32, Nxg_32, Nyg_32        integer*4 im_32, jm_32, Nxg_32, Nyg_32
68        integer*4 iveg_32(im,jm,maxtyp,Nxg,Nyg)        integer*4 iveg_32(im,jm,maxtyp,Nxg,Nyg)
# Line 67  C*************************************** Line 70  C***************************************
70    
71        integer i,j,k,bi,bj,ierr1,kveg        integer i,j,k,bi,bj,ierr1,kveg
72    
73    #ifdef ALLOW_MNC
74          character*(MAX_LEN_FNAM) fizhi_veg_bn
75    #endif
76    
77    C     Allow for MDSIO format if someday needed
78    #ifdef ALLOW_MDSIO
79          IF ( .FALSE. ) THEN
80    
81        call mdsfindunit( kveg, myThid )        call mdsfindunit( kveg, myThid )
82        close(kveg)        close(kveg)
83        open(kveg,file=vegdata,form='unformatted',access='sequential',        open(kveg,file=vegdata,form='unformatted',access='sequential',
# Line 80  C*************************************** Line 91  C***************************************
91         call exit(101)         call exit(101)
92        endif        endif
93        close(kveg)        close(kveg)
94    
95        IF (myThid.eq.1) THEN        IF (myThid.eq.1) THEN
96        imdata = im_32        imdata = im_32
97        jmdata = jm_32        jmdata = jm_32
# Line 95  C*************************************** Line 107  C***************************************
107         print *         print *
108         call exit(102)         call exit(102)
109        ENDIF        ENDIF
110                                                                                          ENDIF
111    
112        DO BJ = myByLo(myThid), myByHi(myThid)        DO BJ = myByLo(myThid), myByHi(myThid)
113        DO BI = myBxLo(myThid), myBxHi(myThid)        DO BI = myBxLo(myThid), myBxHi(myThid)
114    
115        biglobal=bi+(myXGlobalLo-1)/im        biglobal=bi+(myXGlobalLo-1)/im
116        bjglobal=bj+(myYGlobalLo-1)/jm        bjglobal=bj+(myYGlobalLo-1)/jm
117                                                                                    #if defined( _BYTESWAPIO ) && defined( ALLOW_MDSIO )
118          call MDS_BYTESWAPR4(im*jm*maxtyp,veg_32(1,1,1,biglobal,bjglobal))
119    #endif
120        do k = 1,maxtyp        do k = 1,maxtyp
121        do j = 1,jm        do j = 1,jm
122        do i = 1,im        do i = 1,im
# Line 111  C*************************************** Line 126  C***************************************
126        enddo        enddo
127        enddo        enddo
128    
129  c   create chip arrays for :        ENDDO
130          ENDDO
131    
132          ENDIF
133    #endif
134    
135    #ifdef ALLOW_MNC
136          _BEGIN_MASTER( myThid )
137          
138          do i = 1,MAX_LEN_FNAM
139            fizhi_veg_bn(i:i) = ' '
140          enddo
141          
142    C     The following base name should be handled by some sort of input
143    C     name parameter in FIZHI_READPARMS() plus a possible size.
144          
145    C     Set the base name     1234567890
146          fizhi_veg_bn(1:10) = 'fizhi_veg '
147          
148          CALL MNC_CW_I_R('I', fizhi_veg_bn, 0,0,
149         &     'surftype',surftype,myThid)
150          CALL MNC_CW_RL_R('R', fizhi_veg_bn, 0,0,
151         &     'tilefrac', tilefrac, myThid)
152          
153          _END_MASTER( myThid )
154    #endif
155    
156    c     create chip arrays for :
157  c      igrd :  grid index  c      igrd :  grid index
158  c      ityp :  veg. type  c      ityp :  veg. type
159  c      chfr :  vegetation fraction  c      chfr :  vegetation fraction
160  c      chlon:  chip longitude  c      chlon:  chip longitude
161  c      chlt :  chip latitude  c      chlt :  chip latitude
162    
163  c  nchplnd<=nchplocal is the actual number of land chips  c     nchpland<=nchptot is the actual number of land chips
   
 c land points  
 c -----------  
       nchplnd = 0  
       do k=1,maxtyp  
       do j=1,jm  
       do i=1,im  
       if(surftype(i,j,k,bi,bj).lt.100.and.  
      .                                  tilefrac(i,j,k,bi,bj).gt.0.)then  
              nchplnd  = nchplnd + 1  
        igrd (nchplnd,bi,bj) = i + (j-1)*im  
        ityp (nchplnd,bi,bj) = surftype(i,j,k,bi,bj)  
        chfr (nchplnd,bi,bj) = tilefrac(i,j,k,bi,bj)  
        chlon(nchplnd,bi,bj) = lons(i,j,bi,bj)  
        chlt (nchplnd,bi,bj) = lats(i,j,bi,bj)  
       endif  
       enddo  
       enddo  
       enddo  
   
 c ocean points  
 c ------------  
       nchplocal = nchplnd  
   
       do k=1,maxtyp  
       do j=1,jm  
       do i=1,im  
       if(surftype(i,j,k,bi,bj).ge.100 .and.  
      .                                  tilefrac(i,j,k,bi,bj).gt.0.)then  
              nchplocal  = nchplocal + 1  
        igrd (nchplocal,bi,bj) = i + (j-1)*im  
        ityp (nchplocal,bi,bj) = surftype(i,j,k,bi,bj)  
        chfr (nchplocal,bi,bj) = tilefrac(i,j,k,bi,bj)  
        chlon(nchplocal,bi,bj) = lons(i,j,bi,bj)  
        chlt (nchplocal,bi,bj) = lats(i,j,bi,bj)  
       endif  
       enddo  
       enddo  
       enddo  
164    
165        print *, 'bi ',bi,' bj ',bj        DO BJ = myByLo(myThid), myByHi(myThid)
166        print *, 'Number of Total Tiles: ',nchplocal          DO BI = myBxLo(myThid), myBxHi(myThid)
       print *, 'Number of Land  Tiles: ',nchplnd  
       print *  
167    
168    c         land points
169    c         -----------
170              nchpland = 0
171              do k=1,maxtyp
172                do j=1,jm
173                  do i=1,im
174                    if(surftype(i,j,k,bi,bj).lt.100 .and.
175         .               tilefrac(i,j,k,bi,bj).gt.0.) then
176                      nchpland  = nchpland + 1
177                      igrd (nchpland,bi,bj) = i + (j-1)*im
178                      ityp (nchpland,bi,bj) = surftype(i,j,k,bi,bj)
179                      chfr (nchpland,bi,bj) = tilefrac(i,j,k,bi,bj)
180                      chlon(nchpland,bi,bj) = lons(i,j,bi,bj)
181                      chlt (nchpland,bi,bj) = lats(i,j,bi,bj)
182                    endif
183                  enddo
184                enddo
185              enddo
186              
187    c         ocean points
188    c         ------------
189              nchptot = nchpland
190              
191              do k=1,maxtyp
192                do j=1,jm
193                  do i=1,im
194                    if(surftype(i,j,k,bi,bj).ge.100 .and.
195         .               tilefrac(i,j,k,bi,bj).gt.0.) then
196                      nchptot  = nchptot + 1
197                      igrd (nchptot,bi,bj) = i + (j-1)*im
198                      ityp (nchptot,bi,bj) = surftype(i,j,k,bi,bj)
199                      chfr (nchptot,bi,bj) = tilefrac(i,j,k,bi,bj)
200                      chlon(nchptot,bi,bj) = lons(i,j,bi,bj)
201                      chlt (nchptot,bi,bj) = lats(i,j,bi,bj)
202                    endif
203                  enddo
204                enddo
205              enddo
206              
207              if(bi.eq.1.and.bj.eq.1)then
208              print *, 'Number of Total Tiles: ',nchptot
209              print *, 'Number of Land  Tiles: ',nchpland
210              print *
211              endif
212              
213            ENDDO
214        ENDDO        ENDDO
215        ENDDO        
   
216        RETURN        RETURN
217        END        END
Legend:
Removed from v.1.9  
changed lines
  Added in v.1.15
  ViewVC Help
Powered by ViewVC 1.1.22