SUBROUTINE VFRAC(glat1,glon1,SM, & SICE,VEGFRC,im,jm,nm) !******************************************************************************************* ! !!!!!!!!!==========================!!!!!!!!!! ! !!!!!!!!! DRAGAN, August, 2011 !!!!!!!!!! ! !!!!!!!!!==========================!!!!!!!!!! ! Subroutine modified! There were some problems in calling subroutine 'GETGB'. ! If you try to call it twice, it doesn't read well the data from grib file, ! so the values for FPAR2 were not good, and after that values for vegfrc, as well. ! So, each month was read independently, the values for FPAR were calculated ! and saved in the files '../data_out/FPAR_12mon'//cn//'.bin'. ! For initial conditions of veg. greeness for specific date, monthly data of ! two consecutive months are read, interpolated and transfered to E-grid. ! There is also a problem with calling subroutine 'putveg'. It also can't be called ! more then once, so it was impossible to use it in a 'do' loop, and read montlhy ! values of veg. greenes for a whole year. Monthly values were also created manually, month by month, ! and saved in the files '../data_out/VGREEN_12mon'//cn3//'.bin'. Later in the forecast part, ! these data are just read from these files, when it is necessary to update veg. greeness. ! The problem with the subroutines 'GETGB' and 'putveg' can probably be solved in a better ! way, but it would take too much time, that's why all the data were created manualy ! and they can be used on a regular basis. !******************************************************************************************** include 'const.h' ! include masks.comm integer:: im,jm,nm real:: glat1(0:im+1,0:jm+1,nm), glon1(0:im+1,0:jm+1,nm) real:: glat(0:im+1,0:jm+1,nm), glon(0:im+1,0:jm+1,nm) real:: sm(0:im+1,0:jm+1,nm),sice(0:im+1,0:jm+1,nm),vegfrc(0:im+1,0:jm+1,nm) !real*4 dum2(0:im+1,0:jm+1,nm),dum2b(0:im+1,0:jm+1,nm) real:: dum2(0:im+1,0:jm+1,nm),dum2b(0:im+1,0:jm+1,nm) !integer*4 JULD,MON2,MON1,day1,day2,JULM(13) integer:: JULD,MON2,MON1,day1,day2,JULM(13) real*4:: wght1,wght2,rday,dtr !real*4 TLM0D,TPH0D,DLMD,DPHD real:: vegfrm(0:im+1,0:jm+1,nm) ! monthly veg greeness integer:: imo,imo2 character(len=112)::green, green2, green_mon character(len=2)::cn, cn2, cn3 real:: vegfrc1(0:im+1,0:jm+1,nm),vegfrc2(0:im+1,0:jm+1,nm) ! real,dimension(:,:,:),allocatable:: vegfrm data JULM/0,31,59,90,120,151,181,212,243,273,304,334,365/ INTEGER KPDS(200),KGDS(200),JPDS(200),JGDS(200),KF,KNUM REAL*4 FPAR1(2500,1250),FPAR2(2500,1250),FSVE LOGICAL BITMAP(2500,1250) ! lveg=.true. ! variable veg. greeness lveg=.false. ! FOURTH, MONTHLY GREEN VEG FRACTION, INTERPOLATE TO DAY OF YEAR. ! THE MONTHLY FIELDS ASSUMED VALID AT 15TH OF MONTH. ! VALUES RANGE FROM 0.001 TO 1.0 OVER LAND, 0.0 OVER WATER. ! INTERPOLATE IN SPACE. ! ! READ FROM A FILE THAT HAS 13 RECORDS, ONE PER MONTH WITH ! JANUARY OCCURRING BOTH AS RECORD 1 AND AGAIN AS RECORD 13, ! THE LATTER TO SIMPLIFY TIME INTERPOLATION FOR DAYS ! BETWEEN DEC 16 AND JAN 15. WE TREAT JAN 1 TO JAN 15 ! AS JULIAN DAYS 366 TO 380 BELOW, I.E WRAP AROUND YEAR. ! *** THE FOLLOWING PART WAS REVISED BY F. CHEN 7/96 TO REFLECT ! A NEW NESDIS VEGETATION FRACTION PRODUCT (FIVE-YEAR ! CLIMATOLOGY WITH 0.144 DEGREE RESOLUTION ! FROM 89.928S, 180W TO 89.928N, 180E) ! glon=-glon1 glat=glat1 print *,glat(1,1,1),glon(1,1,1) REWIND 38 ! ! **** DO TIME INTERPOLATION **** JULD=JULM(IDAT(1))+IDAT(2) IF(JULD.LE.15) JULD=JULD+365 MON2=IDAT(1) IF(IDAT(2).GT.15) MON2=MON2+1 IF(MON2.EQ.1) MON2=13 MON1=MON2-1 ! **** ASSUME DATA VALID AT 15TH OF MONTH DAY2=JULM(MON2)+15 DAY1=JULM(MON1)+15 RDAY=JULD WGHT1=(DAY2-RDAY)/(DAY2-DAY1) WGHT2=(RDAY-DAY1)/(DAY2-DAY1) ! write(0,*)day2,day1,wght1,wght2 !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ !!!!!!! CALL BAOPEN(38,'data_in/static2/veg.eta.grb',IERR) !!!!!!! JPDS = -1 !!!!!!! CALL GETGB(38,0,2500*1250,MON1-1,JPDS,JGDS,KF,KNUM, & !!!!!!! KPDS,KGDS,BITMAP,FPAR1,IERR) !!!!!!! CALL GETGB(38,0,2500*1250,imo-1,JPDS,JGDS,KF,KNUM, & !!!!!!! KPDS,KGDS,BITMAP,FPAR1,IERR) ! WRITE(6,*) 'AFTER GETGB FOR MONTH ', MON1,' IRET=', IERR !!!!!!! write(66)FPAR1 !original part of the code, commented !!!!!!! close(66) !DRAGAN, August 2011 !!!!!!!!! CALL PUTVEG(GLAT,GLON, & !!!!!!!!! FPAR1,SM,SICE,vegfrc,im,jm,nm) !!!!!!!! JPDS = -1 !!!!!!!!! CALL GETGB(38,0,2500*1250,MON2-1,JPDS,JGDS,KF,KNUM, & !!!!!!!!! KPDS,KGDS,BITMAP,FPAR2,IERR) !!!!!!!!!! WRITE(6,*) 'AFTER GETGB FOR MONTH ', MON2,' IRET=', IERR !!!!!!!!!! CALL PUTVEG(GLAT,GLON, & !!!!!!!!!! FPAR1,SM,SICE,vegfrc,im,jm,nm) write(0,*)'MON1=',MON1,'MON2=',MON2 !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ write(cn,'(i2.2)') MON1 ! green='../data_out/VGREEN_12mon'//cn//'.bin' green='../data_in/init/FPAR_12mon'//cn//'.bin' open(65,file=green,form='unformatted',status='unknown') read(65)FPAR1 close(65) write(cn2,'(i2.2)') MON2 ! green2='../data_out/VGREEN_12mon'//cn2//'.bin' green2='../data_in/init/FPAR_12mon'//cn//'.bin' ! declaring the initial conditions of veg. greenes for specific date, ! reading the data of two consecutive months, created manualy, open(77,file=green2,form='unformatted',status='unknown') ! and then interpolated and transfered to E-grid !DRAGAN, August 2011 read(77)FPAR2 close(77) ! write(0,*)'FPAR1=',FPAR1 DO JJ=1,1250 DO I=1,2500 ! print *,wght1,fpar1(i,jj),wght2,fpar2(i,jj) FSVE=WGHT1*FPAR1(I,JJ)+WGHT2*FPAR2(I,JJ) FPAR1(I,JJ)=FSVE END DO END DO ! ** SPACE INTERPOLATION OF FPAR1 TO E GRID !============================================ CALL PUTVEG(GLAT,GLON, & FPAR1,SM,SICE,vegfrc,im,jm,nm) print *,im/2,jm/2 write(0,*) 'veg fracs ', vegfrc(im/2,jm/2,1),vegfrc(im/2+5,jm/2+5,1) !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ !============================DRAGAN, August 2011=============================== ! READING THE FILES OF MONTHLY VEG. GREENESS if(lveg) then do imo=1,12 write(cn3,'(i2.2)') imo green_mon='VGREEN_12mon'//cn3//'.bin' open(55,file=green_mon,form='unformatted',status='unknown') !DRAGAN, August 2011 read(55)vegfrm close(55) enddo endif !=============================================================================== RETURN END