117 CHARACTER(LEN=3) :: donst
118 INTEGER :: idim, jdim, lsoil, lugb, iy, im, id, ih, ialb
119 INTEGER :: isot, ivegsrc, lensfc, zsea1_mm, zsea2_mm, ierr
120 INTEGER :: nprocs, myrank, num_threads,
num_parthds, max_tasks
121 REAL :: fh, deltsfc, zsea1, zsea2
122 LOGICAL :: use_ufo, do_nsst, do_lndinc, do_sfccycle, frac_grid
124 NAMELIST/namcyc/ idim,jdim,lsoil,lugb,iy,im,id,ih,fh,&
125 deltsfc,ialb,use_ufo,donst, &
126 do_sfccycle,isot,ivegsrc,zsea1_mm, &
127 zsea2_mm, max_tasks, do_lndinc, frac_grid
129 DATA idim,jdim,lsoil/96,96,4/
130 DATA iy,im,id,ih,fh/1997,8,2,0,0./
131 DATA lugb/51/, deltsfc/0.0/, ialb/1/, max_tasks/99999/
132 DATA isot/1/, ivegsrc/2/, zsea1_mm/0/, zsea2_mm/0/
135 CALL mpi_comm_size(mpi_comm_world, nprocs, ierr)
136 CALL mpi_comm_rank(mpi_comm_world, myrank, ierr)
138 if (myrank==0)
call w3tagb(
'GLOBAL_CYCLE',2018,0179,0055,
'NP20')
143 print*,
"STARTING CYCLE PROGRAM ON RANK ", myrank
144 print*,
"RUNNING WITH ", nprocs,
"TASKS" 145 print*,
"AND WITH ", num_threads,
" THREADS." 154 print*,
"READ NAMCYC NAMELIST." 156 CALL baopenr(36,
"fort.36", ierr)
158 print*,
'FATAL ERROR READING FORT.36 NAMELIST. IERR: ', ierr
159 CALL mpi_abort(mpi_comm_world, 32, ierr)
161 READ(36, nml=namcyc, iostat=ierr)
163 print*,
'FATAL ERROR READING FORT.36 NAMELIST. IERR: ', ierr
164 CALL mpi_abort(mpi_comm_world, 33, ierr)
168 IF (max_tasks < 99999 .AND. myrank > (max_tasks - 1))
THEN 169 print*,
"USER SPECIFIED MAX NUMBER OF TASKS: ", max_tasks
170 print*,
"WILL NOT RUN CYCLE PROGRAM ON RANK: ", myrank
176 zsea1 = float(zsea1_mm) / 1000.0
177 zsea2 = float(zsea2_mm) / 1000.0
179 IF (donst ==
"YES")
THEN 186 IF (myrank==0) print*,
"LUGB,IDIM,JDIM,ISOT,IVEGSRC,LSOIL,DELTSFC,IY,IM,ID,IH,FH: ", &
187 lugb,idim,jdim,isot,ivegsrc,lsoil,deltsfc,iy,im,id,ih,fh
189 CALL sfcdrv(lugb,idim,jdim,lensfc,lsoil,deltsfc, &
190 iy,im,id,ih,fh,ialb, &
191 use_ufo,do_nsst,do_sfccycle,do_lndinc, &
192 frac_grid,zsea1,zsea2,isot,ivegsrc,myrank)
195 print*,
'CYCLE PROGRAM COMPLETED NORMALLY ON RANK: ', myrank
199 CALL mpi_barrier(mpi_comm_world, ierr)
201 if (myrank==0)
call w3tage(
'GLOBAL_CYCLE')
203 CALL mpi_finalize(ierr)
317 SUBROUTINE sfcdrv(LUGB, IDIM,JDIM,LENSFC,LSOIL,DELTSFC, &
318 IY,IM,ID,IH,FH,IALB, &
319 USE_UFO,DO_NSST,DO_SFCCYCLE,DO_LNDINC,&
320 FRAC_GRID,ZSEA1,ZSEA2,ISOT,IVEGSRC,MYRANK)
325 USE land_increments,
ONLY: gaussian_to_fv3_interp, &
326 add_increment_soil, &
327 add_jedi_increment_snow, &
328 calculate_landinc_mask, &
329 apply_land_da_adjustments_soil, &
330 apply_land_da_adjustments_snd, &
335 INTEGER,
INTENT(IN) :: IDIM, JDIM, LENSFC, LSOIL, IALB
336 INTEGER,
INTENT(IN) :: LUGB, IY, IM, ID, IH
337 INTEGER,
INTENT(IN) :: ISOT, IVEGSRC, MYRANK
339 LOGICAL,
INTENT(IN) :: USE_UFO, DO_NSST,DO_SFCCYCLE
340 LOGICAL,
INTENT(IN) :: DO_LNDINC, FRAC_GRID
342 REAL,
INTENT(IN) :: FH, DELTSFC, ZSEA1, ZSEA2
344 INTEGER,
PARAMETER :: NLUNIT=35
345 INTEGER,
PARAMETER :: SZ_NML=1
347 CHARACTER(LEN=5) :: TILE_NUM
348 CHARACTER(LEN=500) :: NST_FILE
349 CHARACTER(LEN=500) :: GSI_SOI_FILE,JEDI_SOI_FILE,JEDI_SNO_FILE
350 CHARACTER(LEN=4) :: INPUT_NML_FILE(SZ_NML)
353 INTEGER :: I_INDEX(LENSFC), J_INDEX(LENSFC)
354 INTEGER :: IDUM(IDIM,JDIM)
355 integer :: num_parthds, num_threads
360 real(kind=kind_io8) :: min_ice(lensfc)
362 REAL :: SLMASK(LENSFC), OROG(LENSFC)
363 REAL :: SIHFCS(LENSFC), SICFCS(LENSFC)
364 REAL :: SITFCS(LENSFC), TSFFCS(LENSFC)
365 REAL :: SWEFCS(LENSFC), ZORFCS(LENSFC)
366 REAL :: ALBFCS(LENSFC,4), TG3FCS(LENSFC)
367 REAL :: CNPFCS(LENSFC), SMCFCS(LENSFC,LSOIL)
368 REAL :: STCFCS(LENSFC,LSOIL), SLIFCS(LENSFC)
369 REAL :: AISFCS(LENSFC), F10M(LENSFC)
370 REAL :: VEGFCS(LENSFC), VETFCS(LENSFC)
371 REAL :: SOTFCS(LENSFC), ALFFCS(LENSFC,2)
372 REAL :: CVFCS(LENSFC), CVTFCS(LENSFC)
373 REAL :: CVBFCS(LENSFC), TPRCP(LENSFC)
374 REAL :: SRFLAG(LENSFC), SNDFCS(LENSFC)
375 REAL :: SLCFCS(LENSFC,LSOIL), VMXFCS(LENSFC)
376 REAL :: VMNFCS(LENSFC), T2M(LENSFC)
377 REAL :: Q2M(LENSFC), SLPFCS(LENSFC)
378 REAL :: ABSFCS(LENSFC), OROG_UF(LENSFC)
379 REAL :: USTAR(LENSFC), SOCFCS(LENSFC)
380 REAL :: FMM(LENSFC), FHH(LENSFC)
381 REAL :: RLA(LENSFC), RLO(LENSFC)
382 REAL(KIND=4) :: ZSOIL(LSOIL)
383 REAL :: SIG1T(LENSFC)
386 REAL,
ALLOCATABLE :: STC_BCK(:,:), SMC_BCK(:,:), SLC_BCK(:,:)
387 REAL,
ALLOCATABLE :: SLIFCS_FG(:), SICFCS_FG(:)
388 INTEGER,
ALLOCATABLE :: LANDINC_MASK_FG(:), LANDINC_MASK(:)
389 REAL,
ALLOCATABLE :: SND_BCK(:), SND_INC(:), SWE_BCK(:)
390 REAL(KIND=KIND_IO8),
ALLOCATABLE :: SLMASKL(:), SLMASKW(:), LANDFRAC(:)
392 TYPE(NSST_DATA) :: NSST
393 real,
dimension(idim,jdim) :: tf_clm,tf_trd,sal_clm
394 real,
dimension(lensfc) :: tf_clm_tile,tf_trd_tile,sal_clm_tile
395 INTEGER :: veg_type_landice
396 INTEGER,
DIMENSION(LENSFC) :: STC_UPDATED, SLC_UPDATED
397 REAL,
DIMENSION(LENSFC,LSOIL) :: STCINC, SLCINC
399 LOGICAL :: FILE_EXISTS, DO_SOI_INC_GSI, DO_SOI_INC_JEDI, DO_SNO_INC_JEDI
400 CHARACTER(LEN=3) :: RANKCH
401 INTEGER :: lsoil_incr
408 NAMELIST/namsfcd/ nst_file, lsoil_incr, do_sno_inc_jedi, do_soi_inc_jedi, do_soi_inc_gsi
410 DATA nst_file/
'NULL'/
412 do_sno_inc_jedi = .false.
413 do_soi_inc_gsi = .false.
414 do_soi_inc_jedi = .false.
420 input_nml_file =
"NULL" 422 CALL baopenr(37,
"fort.37", ierr)
424 print*,
'FATAL ERROR OPENING FORT.37 NAMELIST. IERR: ', ierr
425 CALL mpi_abort(mpi_comm_world, 30, ierr)
427 READ (37, nml=namsfcd, iostat=ierr)
429 print*,
'FATAL ERROR READING FORT.37 NAMELIST. IERR: ', ierr
430 CALL mpi_abort(mpi_comm_world, 31, ierr)
434 print*,
'IN ROUTINE SFCDRV,IDIM=',idim,
'JDIM=',jdim,
'FH=',fh
440 ALLOCATE(landfrac(lensfc))
442 print*,
'- RUNNING WITH FRACTIONAL GRID.' 443 CALL read_lat_lon_orog(rla,rlo,orog,orog_uf,tile_num,idim,jdim,lensfc,landfrac=landfrac)
453 i_index = reshape(idum, (/lensfc/))
459 j_index = reshape(idum, (/lensfc/) )
463 print*,
"WILL PROCESS NSST RECORDS." 464 ALLOCATE(nsst%C_0(lensfc))
465 ALLOCATE(nsst%C_D(lensfc))
466 ALLOCATE(nsst%D_CONV(lensfc))
467 ALLOCATE(nsst%DT_COOL(lensfc))
468 ALLOCATE(nsst%IFD(lensfc))
469 ALLOCATE(nsst%QRAIN(lensfc))
470 ALLOCATE(nsst%TREF(lensfc))
471 ALLOCATE(nsst%TFINC(lensfc))
472 ALLOCATE(nsst%W_0(lensfc))
473 ALLOCATE(nsst%W_D(lensfc))
474 ALLOCATE(nsst%XS(lensfc))
475 ALLOCATE(nsst%XT(lensfc))
476 ALLOCATE(nsst%XTTS(lensfc))
477 ALLOCATE(nsst%XU(lensfc))
478 ALLOCATE(nsst%XV(lensfc))
479 ALLOCATE(nsst%XZ(lensfc))
480 ALLOCATE(nsst%XZTS(lensfc))
481 ALLOCATE(nsst%Z_C(lensfc))
482 ALLOCATE(nsst%ZM(lensfc))
483 ALLOCATE(slifcs_fg(lensfc))
484 ALLOCATE(sicfcs_fg(lensfc))
489 IF (do_soi_inc_gsi .and. do_soi_inc_jedi)
THEN 491 print*,
'FATAL ERROR: Can not do gsi and jedi soil updates at the same time, choose one!' 492 CALL mpi_abort(mpi_comm_world, 15, ierr)
494 IF (do_soi_inc_gsi .or. do_soi_inc_jedi)
THEN 496 print*,
" APPLYING SOIL INCREMENTS FROM GSI OR JEDI" 497 ALLOCATE(stc_bck(lensfc, lsoil), smc_bck(lensfc, lsoil), slc_bck(lensfc,lsoil))
498 ALLOCATE(landinc_mask_fg(lensfc))
501 IF (do_sno_inc_jedi)
THEN 506 print*,
" APPLYING SNOW INCREMENTS FROM JEDI" 507 ALLOCATE(snd_bck(lensfc), snd_inc(lensfc), swe_bck(lensfc))
510 ALLOCATE(landinc_mask(lensfc))
511 if (ivegsrc == 2)
then 522 CALL read_data(lsoil,lensfc,do_nsst,do_sno_inc_jedi,do_soi_inc_jedi,.false.,is_noahmp=is_noahmp, &
523 tsffcs=tsffcs,smcfcs=smcfcs, &
524 swefcs=swefcs,stcfcs=stcfcs,tg3fcs=tg3fcs,zorfcs=zorfcs, &
525 cvfcs=cvfcs, cvbfcs=cvbfcs,cvtfcs=cvtfcs,albfcs=albfcs, &
526 vegfcs=vegfcs,slifcs=slifcs,cnpfcs=cnpfcs,f10m=f10m , &
527 vetfcs=vetfcs,sotfcs=sotfcs,alffcs=alffcs,ustar=ustar , &
528 fmm=fmm ,fhh=fhh ,sihfcs=sihfcs,sicfcs=sicfcs, &
529 sitfcs=sitfcs,tprcp=tprcp ,srflag=srflag,sndfcs=sndfcs, &
530 vmnfcs=vmnfcs,vmxfcs=vmxfcs,slcfcs=slcfcs,slpfcs=slpfcs, &
531 absfcs=absfcs,t2m=t2m ,q2m=q2m ,slmask=slmask, &
532 zsoil=zsoil, nsst=nsst)
534 IF (frac_grid .AND. .NOT. is_noahmp)
THEN 535 print *,
'FATAL ERROR: NOAH lsm update does not work with fractional grids.' 536 call mpi_abort(mpi_comm_world, 18, ierr)
539 IF (is_noahmp .AND. do_sno_inc_jedi)
THEN 540 print *,
'FATAL ERROR: Snow increment update does not work with NOAH_MP.' 541 call mpi_abort(mpi_comm_world, 29, ierr)
552 print*,
'USE UNFILTERED OROGRAPHY.' 559 IF(nint(slifcs(i)).EQ.2) aisfcs(i) = 1.
564 IF (.NOT. do_sfccycle )
THEN 566 print*,
"FIRST GUESS MASK ADJUSTED BY IFD RECORD" 568 WHERE(nint(nsst%IFD) == 3) slifcs_fg = 2.0
571 print*,
"SAVE FIRST GUESS MASK" 578 CALL calculate_landinc_mask(swefcs, vetfcs, sotfcs, &
579 lensfc,veg_type_landice, landinc_mask)
589 IF (do_sfccycle)
THEN 590 ALLOCATE(slmaskl(lensfc), slmaskw(lensfc))
592 set_mask :
IF (frac_grid)
THEN 595 IF(landfrac(i) > 0.0_kind_io8)
THEN 596 slmaskl(i) = ceiling(landfrac(i)-1.0e-6_kind_io8)
597 slmaskw(i) = floor(landfrac(i)+1.0e-6_kind_io8)
599 IF(nint(slmask(i)) == 1)
THEN 601 print*,
'FATAL ERROR: LAND FRAC AND SLMASK MISMATCH.' 602 CALL mpi_abort(mpi_comm_world, 27, ierr)
604 slmaskl(i) = 0.0_kind_io8
605 slmaskw(i) = 0.0_kind_io8
616 IF(nint(slmask(i)) == 1)
THEN 617 slmaskl(i) = 1.0_kind_io8
618 slmaskw(i) = 1.0_kind_io8
620 slmaskl(i) = 0.0_kind_io8
621 slmaskw(i) = 0.0_kind_io8
628 if(nint(slmask(i)) == 0)
then 629 min_ice(i) = 0.15_kind_io8
631 min_ice(i) = 0.0_kind_io8
637 num_threads = num_parthds()
639 print*,
"CALL SFCCYCLE TO UPDATE SURFACE FIELDS." 640 CALL sfccycle(lugb,lensfc,lsoil,sig1t,deltsfc, &
641 iy,im,id,ih,fh,rla,rlo, &
642 slmaskl,slmaskw, orog, orog_uf, use_ufo, do_nsst, &
643 sihfcs,sicfcs,sitfcs,sndfcs,slcfcs, &
644 vmnfcs,vmxfcs,slpfcs,absfcs, &
645 tsffcs,swefcs,zorfcs,albfcs,tg3fcs, &
646 cnpfcs,smcfcs,stcfcs,slifcs,aisfcs, &
647 vegfcs,vetfcs,sotfcs,socfcs,alffcs, &
648 cvfcs,cvbfcs,cvtfcs,myrank,num_threads, nlunit, &
649 sz_nml, input_nml_file, &
651 ialb,isot,ivegsrc,tile_num,i_index,j_index)
653 DEALLOCATE(slmaskl, slmaskw)
663 IF (nst_file ==
"NULL")
THEN 665 print*,
"NO GSI FILE. ADJUST IFD FOR FORMER ICE POINTS." 667 IF (sicfcs_fg(i) > 0.0 .AND. sicfcs(i) == 0)
THEN 674 print*,
"ADJUST TREF FROM GSI INCREMENT" 678 call get_tf_clm(rla,rlo,jdim,idim,iy,im,id,ih,tf_clm,tf_trd)
679 tf_clm_tile(:) = reshape(tf_clm, (/lensfc/) )
680 tf_trd_tile(:) = reshape(tf_trd, (/lensfc/) )
684 call get_sal_clm(rla,rlo,jdim,idim,iy,im,id,ih,sal_clm)
685 sal_clm_tile(:) = reshape(sal_clm, (/lensfc/) )
693 CALL adjust_nsst(rla,rlo,slifcs,slifcs_fg,tsffcs,sitfcs,sicfcs,sicfcs_fg,&
694 stcfcs,nsst,lensfc,lsoil,idim,jdim,zsea1,zsea2, &
695 tf_clm_tile,tf_trd_tile,sal_clm_tile,landfrac,frac_grid)
707 IF (do_sno_inc_jedi)
THEN 713 WRITE(rankch,
'(I3.3)') (myrank+1)
714 jedi_sno_file =
"snow_xainc." // rankch
716 INQUIRE(file=trim(jedi_sno_file), exist=file_exists)
717 IF (.not. file_exists)
then 718 print *,
'FATAL ERROR: snow increment (fv3 grid) update requested, & 719 but file does not exist : ', trim(jedi_sno_file)
720 call mpi_abort(mpi_comm_world, 10, ierr)
728 CALL read_data(lsoil,lensfc,.false.,.true.,.false.,.true.,sndfcs=snd_inc)
738 CALL add_jedi_increment_snow(snd_inc,landinc_mask,lensfc,sndfcs)
744 CALL apply_land_da_adjustments_snd(lsm, lensfc, landinc_mask, swe_bck, snd_bck, &
750 landinc_mask_fg = landinc_mask
752 IF (do_sfccycle .OR. do_sno_inc_jedi)
THEN 753 CALL calculate_landinc_mask(swefcs, vetfcs, sotfcs, lensfc, &
754 veg_type_landice, landinc_mask)
763 IF (do_soi_inc_gsi)
THEN 769 WRITE(rankch,
'(I3.3)') (myrank+1)
770 gsi_soi_file =
"sfcincr_gsi." // rankch
772 INQUIRE(file=trim(gsi_soi_file), exist=file_exists)
773 IF (.not. file_exists)
then 774 print *,
'FATAL ERROR: gsi soil increment (gaussian grid) update requested, & 775 but file does not exist : ', trim(gsi_soi_file)
776 call mpi_abort(mpi_comm_world, 10, ierr)
785 CALL gaussian_to_fv3_interp(lsoil_incr,rla,rlo,&
786 stcinc,slcinc,landinc_mask,lensfc,lsoil,idim,jdim,lsm,myrank)
791 CALL write_data(lensfc,idim,jdim,lsoil,do_nsst,.true.,nsst, &
792 stcinc=stcinc,slcinc=slcinc)
796 IF (do_soi_inc_jedi)
THEN 802 WRITE(rankch,
'(I3.3)') (myrank+1)
803 jedi_soi_file =
"soil_xainc." // rankch
805 INQUIRE(file=trim(jedi_soi_file), exist=file_exists)
806 IF (.not. file_exists)
then 807 print *,
'FATAL ERROR: soil increment (fv3 grid) update requested, but file & 808 does not exist: ', trim(jedi_soi_file)
809 call mpi_abort(mpi_comm_world, 10, ierr)
812 CALL read_data(lsoil,lensfc,.false.,.false.,.true., &
813 .true.,is_noahmp=is_noahmp, &
814 stcinc=stcinc,slcinc=slcinc)
818 IF (do_soi_inc_gsi .or. do_soi_inc_jedi)
THEN 825 CALL add_increment_soil(lsoil_incr,stcinc,slcinc,stcfcs,smcfcs,slcfcs,stc_updated, &
826 slc_updated,landinc_mask,landinc_mask_fg,lensfc,lsoil,lsm,myrank)
832 CALL apply_land_da_adjustments_soil(lsoil_incr, lsm, isot, ivegsrc,lensfc, lsoil, &
833 sotfcs, landinc_mask_fg, stc_bck, stcfcs, smcfcs, slcfcs, stc_updated, &
843 IF(
ALLOCATED(landinc_mask_fg))
DEALLOCATE(landinc_mask_fg)
844 IF(
ALLOCATED(landinc_mask))
DEALLOCATE(landinc_mask)
845 IF(
ALLOCATED(stc_bck))
DEALLOCATE(stc_bck)
846 IF(
ALLOCATED(smc_bck))
DEALLOCATE(smc_bck)
847 IF(
ALLOCATED(slc_bck))
DEALLOCATE(slc_bck)
848 IF(
ALLOCATED(snd_bck))
DEALLOCATE(snd_bck)
849 IF(
ALLOCATED(swe_bck))
DEALLOCATE(swe_bck)
850 IF(
ALLOCATED(snd_inc))
DEALLOCATE(snd_inc)
857 IF (lsm==lsm_noahmp)
THEN 859 CALL write_data(lensfc,idim,jdim,lsoil,do_nsst,.false.,nsst,vegfcs=vegfcs, &
860 slcfcs=slcfcs,smcfcs=smcfcs,stcfcs=stcfcs,&
861 sicfcs=sicfcs,sihfcs=sihfcs)
863 ELSEIF (lsm==lsm_noah)
THEN 866 do_nsst,.false.,nsst,slifcs=slifcs,tsffcs=tsffcs,vegfcs=vegfcs, &
867 swefcs=swefcs,tg3fcs=tg3fcs,zorfcs=zorfcs, &
868 albfcs=albfcs,alffcs=alffcs,cnpfcs=cnpfcs, &
869 f10m=f10m,t2m=t2m,q2m=q2m,vetfcs=vetfcs, &
870 sotfcs=sotfcs,ustar=ustar,fmm=fmm,fhh=fhh, &
871 sicfcs=sicfcs,sihfcs=sihfcs,sitfcs=sitfcs,tprcp=tprcp, &
872 srflag=srflag,swdfcs=sndfcs,vmnfcs=vmnfcs, &
873 vmxfcs=vmxfcs,slpfcs=slpfcs,absfcs=absfcs, &
874 slcfcs=slcfcs,smcfcs=smcfcs,stcfcs=stcfcs)
881 DEALLOCATE(nsst%D_CONV)
882 DEALLOCATE(nsst%DT_COOL)
884 DEALLOCATE(nsst%QRAIN)
885 DEALLOCATE(nsst%TREF)
886 DEALLOCATE(nsst%TFINC)
891 DEALLOCATE(nsst%XTTS)
895 DEALLOCATE(nsst%XZTS)
898 DEALLOCATE(slifcs_fg)
899 DEALLOCATE(sicfcs_fg)
937 SUBROUTINE adjust_nsst(RLA,RLO,SLMSK_TILE,SLMSK_FG_TILE,SKINT_TILE,&
938 SICET_TILE,sice_tile,sice_fg_tile,SOILT_TILE,NSST, &
939 LENSFC,LSOIL,IDIM,JDIM,ZSEA1,ZSEA2, &
940 tf_clm_tile,tf_trd_tile,sal_clm_tile,LANDFRAC, &
953 INTEGER,
INTENT(IN) :: LENSFC, LSOIL, IDIM, JDIM
955 LOGICAL,
INTENT(IN) :: FRAC_GRID
957 REAL,
INTENT(IN) :: SLMSK_TILE(LENSFC), SLMSK_FG_TILE(LENSFC), LANDFRAC(LENSFC)
958 real,
intent(in) :: tf_clm_tile(lensfc),tf_trd_tile(lensfc),sal_clm_tile(lensfc)
959 REAL,
INTENT(IN) :: ZSEA1, ZSEA2,sice_tile(lensfc),sice_fg_tile(lensfc)
960 REAL,
INTENT(IN) :: RLA(LENSFC), RLO(LENSFC)
961 REAL,
INTENT(INOUT) :: SKINT_TILE(LENSFC)
962 REAL,
INTENT(INOUT) :: SICET_TILE(LENSFC),SOILT_TILE(LENSFC,LSOIL)
964 TYPE(NSST_DATA) :: NSST
966 REAL,
PARAMETER :: TMAX=313.0,tzero=273.16
968 INTEGER :: IOPT, NRET, KGDS_GAUS(200)
969 INTEGER :: IGAUS, JGAUS, IJ, II, JJ, III, JJJ, KRAD
970 INTEGER :: ISTART, IEND, JSTART, JEND
972 INTEGER,
allocatable :: MASK_TILE(:),MASK_FG_TILE(:)
973 INTEGER :: ITILE, JTILE
974 INTEGER :: MAX_SEARCH, J, IERR
975 INTEGER :: IGAUSP1, JGAUSP1
976 integer :: nintp,nsearched,nice,nland
977 integer :: nfill,nfill_tice,nfill_clm
978 integer :: nset_thaw,nset_thaw_s,nset_thaw_i,nset_thaw_c
980 INTEGER,
ALLOCATABLE :: ID1(:,:), ID2(:,:), JDC(:,:)
985 REAL :: TREF_SAVE,WSUM,tf_ice,tf_thaw
986 REAL :: FILL, DTZM, GAUS_RES_KM, DTREF
987 REAL,
ALLOCATABLE :: XPTS(:), YPTS(:), LATS(:), LONS(:)
988 REAL,
ALLOCATABLE :: DUM2D(:,:), LATS_RAD(:), LONS_RAD(:)
989 REAL,
ALLOCATABLE :: AGRID(:,:,:), S2C(:,:,:)
998 kgds_gaus(7) = -90000
999 kgds_gaus(8) = nint(-360000./float(
idim_gaus))
1000 kgds_gaus(9) = nint((360.0 / float(
idim_gaus))*1000.0)
1007 print*,
'ADJUST NSST USING GSI INCREMENTS ON GAUSSIAN GRID' 1027 print*,
'FATAL ERROR: PROBLEM IN GDSWZD. STOP.' 1028 CALL mpi_abort(mpi_comm_world, 12, ierr)
1031 DEALLOCATE (xpts, ypts)
1039 lats_rad(j) = dum2d(1,
jdim_gaus-j+1) * 3.1415926 / 180.0
1045 lons_rad = dum2d(:,1) * 3.1415926 / 180.0
1049 ALLOCATE(agrid(idim,jdim,2))
1050 agrid(:,:,1) = reshape(rlo, (/idim,jdim/) )
1051 agrid(:,:,2) = reshape(rla, (/idim,jdim/) )
1052 agrid = agrid * 3.1415926 / 180.0
1054 ALLOCATE(id1(idim,jdim))
1055 ALLOCATE(id2(idim,jdim))
1056 ALLOCATE(jdc(idim,jdim))
1057 ALLOCATE(s2c(idim,jdim,4))
1066 lons_rad, lats_rad, id1, id2, jdc, s2c, agrid )
1068 DEALLOCATE(lons_rad, lats_rad, agrid)
1076 max_search = ceiling(500.0/gaus_res_km)
1079 print*,
'MAXIMUM SEARCH IS ',max_search,
' GAUSSIAN POINTS.' 1102 allocate(mask_tile(lensfc))
1103 allocate(mask_fg_tile(lensfc))
1105 IF(.NOT. frac_grid)
THEN 1106 mask_tile = nint(slmsk_tile)
1107 mask_fg_tile = nint(slmsk_fg_tile)
1110 WHERE(sice_tile > 0.0) mask_tile=2
1111 WHERE(landfrac == 1.0) mask_tile=1
1113 WHERE(sice_fg_tile > 0.0) mask_fg_tile=2
1114 WHERE(landfrac == 1.0) mask_fg_tile=1
1117 ij_loop :
DO ij = 1, lensfc
1122 tf_ice = tfreez(sal_clm_tile(ij)) + tzero
1127 IF (mask_tile(ij) == 1)
THEN 1135 if (mask_tile(ij) == 2)
then 1136 nsst%tref(ij)=tf_ice
1137 skint_tile(ij)=(1.0-sice_tile(ij))*nsst%tref(ij)+sice_tile(ij)*sicet_tile(ij)
1145 jtile = (ij-1) / idim + 1
1146 itile = mod(ij,idim)
1147 IF (itile==0) itile = idim
1155 IF (mask_fg_tile(ij) == 2 .AND. mask_tile(ij) == 0)
THEN 1159 call tf_thaw_set(nsst%tref,mask_fg_tile,itile,jtile,tf_ice,tf_clm_tile(ij),tf_thaw,idim,jdim, &
1160 nset_thaw_s,nset_thaw_i,nset_thaw_c)
1162 nset_thaw = nset_thaw + 1
1178 igaus = id1(itile,jtile)
1179 jgaus = jdc(itile,jtile)
1180 igausp1 = id2(itile,jtile)
1181 jgausp1 = jdc(itile,jtile)+1
1192 dtref = dtref + (s2c(itile,jtile,1) *
dtref_gaus(igaus,jgaus))
1193 wsum = wsum + s2c(itile,jtile,1)
1197 dtref = dtref + (s2c(itile,jtile,2) *
dtref_gaus(igausp1,jgaus))
1198 wsum = wsum + s2c(itile,jtile,2)
1202 dtref = dtref + (s2c(itile,jtile,3) *
dtref_gaus(igausp1,jgausp1))
1203 wsum = wsum + s2c(itile,jtile,3)
1207 dtref = dtref + (s2c(itile,jtile,4) *
dtref_gaus(igaus,jgausp1))
1208 wsum = wsum + s2c(itile,jtile,4)
1212 dtref = dtref / wsum
1214 tref_save = nsst%TREF(ij)
1215 nsst%TREF(ij) = nsst%TREF(ij) + dtref
1216 nsst%TREF(ij) = max(nsst%TREF(ij), tf_ice)
1217 nsst%TREF(ij) = min(nsst%TREF(ij), tmax)
1218 nsst%TFINC(ij) = nsst%TREF(ij) - tref_save
1220 CALL dtzm_point(nsst%XT(ij),nsst%XZ(ij),nsst%DT_COOL(ij), &
1221 nsst%Z_C(ij),zsea1,zsea2,dtzm)
1223 skint_tile(ij) = nsst%TREF(ij) + dtzm
1224 skint_tile(ij) = max(skint_tile(ij), tf_ice)
1225 skint_tile(ij) = min(skint_tile(ij), tmax)
1227 sicet_tile(ij) = skint_tile(ij)
1230 IF(.NOT. frac_grid) soilt_tile(ij,:) = skint_tile(ij)
1241 DO krad = 1, max_search
1243 istart = igaus - krad
1245 jstart = jgaus - krad
1248 DO jj = jstart, jend
1249 DO ii = istart, iend
1251 IF((jj == jstart) .OR. (jj == jend) .OR. &
1252 (ii == istart) .OR. (ii == iend))
THEN 1254 IF ((jj >= 1) .AND. (jj <=
jdim_gaus))
THEN 1271 IF (krad <= 2 .AND.
slmsk_gaus(iii,jjj) == 2) is_ice = .true.
1277 nsearched = nsearched + 1
1279 tref_save = nsst%TREF(ij)
1280 nsst%TREF(ij ) = nsst%TREF(ij) +
dtref_gaus(iii,jjj)
1281 nsst%TREF(ij) = max(nsst%TREF(ij), tf_ice)
1282 nsst%TREF(ij) = min(nsst%TREF(ij), tmax)
1283 nsst%TFINC(ij) = nsst%TREF(ij) - tref_save
1285 CALL dtzm_point(nsst%XT(ij),nsst%XZ(ij),nsst%DT_COOL(ij), &
1286 nsst%Z_C(ij),zsea1,zsea2,dtzm)
1288 skint_tile(ij) = nsst%TREF(ij) + dtzm
1289 skint_tile(ij) = max(skint_tile(ij), tf_ice)
1290 skint_tile(ij) = min(skint_tile(ij), tmax)
1292 sicet_tile(ij) = skint_tile(ij)
1293 IF(.NOT. frac_grid) soilt_tile(ij,:) = skint_tile(ij)
1316 nsst%TREF(ij) = tf_ice
1318 nfill_tice = nfill_tice + 1
1320 tref_save = nsst%TREF(ij)
1321 nsst%TREF(ij) = nsst%TREF(ij) + tf_trd_tile(ij)
1322 nsst%TREF(ij) = max(nsst%TREF(ij), tf_ice)
1323 nsst%TREF(ij) = min(nsst%TREF(ij), tmax)
1324 nsst%TFINC(ij) = nsst%TREF(ij) - tref_save
1326 nfill_clm = nfill_clm + 1
1329 CALL dtzm_point(nsst%XT(ij),nsst%XZ(ij),nsst%DT_COOL(ij), &
1330 nsst%Z_C(ij),zsea1,zsea2,dtzm)
1332 skint_tile(ij) = nsst%TREF(ij) + dtzm
1333 skint_tile(ij) = max(skint_tile(ij), tf_ice)
1334 skint_tile(ij) = min(skint_tile(ij), tmax)
1336 sicet_tile(ij) = skint_tile(ij)
1337 IF (.NOT. frac_grid) soilt_tile(ij,:) = skint_tile(ij)
1343 write(*,
'(a)')
'statistics of grids number processed for tile : ' 1344 write(*,
'(a,I8)')
' nintp = ',nintp
1345 write(*,
'(a,4I8)')
'nset_thaw,nset_thaw_s,nset_thaw_i,nset_thaw_c =',nset_thaw,nset_thaw_s,nset_thaw_i,nset_thaw_c
1346 write(*,
'(a,I8)')
' nsearched = ',nsearched
1347 write(*,
'(a,3I6)')
' nfill,nfill_tice,nfill_clm = ',nfill,nfill_tice,nfill_clm
1348 write(*,
'(a,I8)')
' nice = ',nice
1349 write(*,
'(a,I8)')
' nland = ',nland
1351 DEALLOCATE(id1, id2, jdc, s2c, mask_tile, mask_fg_tile)
1365 SUBROUTINE climo_trend(LATITUDE, MON, DAY, DELTSFC, DTREF)
1368 INTEGER,
INTENT(IN) :: MON, DAY
1370 REAL,
INTENT(IN) :: LATITUDE, DELTSFC
1371 REAL,
INTENT(OUT) :: DTREF
1373 INTEGER :: NUM_DAYS(12), MON2, MON1
1375 REAL,
TARGET :: SST_80_90(12)
1376 REAL,
TARGET :: SST_70_80(12)
1377 REAL,
TARGET :: SST_60_70(12)
1378 REAL,
TARGET :: SST_50_60(12)
1379 REAL,
TARGET :: SST_40_50(12)
1380 REAL,
TARGET :: SST_30_40(12)
1381 REAL,
TARGET :: SST_20_30(12)
1382 REAL,
TARGET :: SST_10_20(12)
1383 REAL,
TARGET :: SST_00_10(12)
1384 REAL,
TARGET :: SST_M10_00(12)
1385 REAL,
TARGET :: SST_M20_M10(12)
1386 REAL,
TARGET :: SST_M30_M20(12)
1387 REAL,
TARGET :: SST_M40_M30(12)
1388 REAL,
TARGET :: SST_M50_M40(12)
1389 REAL,
TARGET :: SST_M60_M50(12)
1390 REAL,
TARGET :: SST_M70_M60(12)
1391 REAL,
TARGET :: SST_M80_M70(12)
1392 REAL,
TARGET :: SST_M90_M80(12)
1394 REAL,
POINTER :: SST(:)
1396 DATA num_days /31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31/
1398 DATA sst_80_90 /271.466, 271.458, 271.448, 271.445, 271.519, 271.636, &
1399 272.023, 272.066, 272.001, 271.698, 271.510, 271.472/
1401 DATA sst_70_80 /272.149, 272.103, 272.095, 272.126, 272.360, 272.988, &
1402 274.061, 274.868, 274.415, 273.201, 272.468, 272.268/
1404 DATA sst_60_70 /274.240, 274.019, 273.988, 274.185, 275.104, 276.875, &
1405 279.005, 280.172, 279.396, 277.586, 275.818, 274.803/
1407 DATA sst_50_60 /277.277, 276.935, 277.021, 277.531, 279.100, 281.357, &
1408 283.735, 285.171, 284.399, 282.328, 279.918, 278.199/
1410 DATA sst_40_50 /281.321, 280.721, 280.850, 281.820, 283.958, 286.588, &
1411 289.195, 290.873, 290.014, 287.652, 284.898, 282.735/
1413 DATA sst_30_40 /289.189, 288.519, 288.687, 289.648, 291.547, 293.904, &
1414 296.110, 297.319, 296.816, 295.225, 292.908, 290.743/
1416 DATA sst_20_30 /294.807, 294.348, 294.710, 295.714, 297.224, 298.703, &
1417 299.682, 300.127, 300.099, 299.455, 297.953, 296.177/
1419 DATA sst_10_20 /298.878, 298.720, 299.033, 299.707, 300.431, 300.709, &
1420 300.814, 300.976, 301.174, 301.145, 300.587, 299.694/
1422 DATA sst_00_10 /300.415, 300.548, 300.939, 301.365, 301.505, 301.141, &
1423 300.779, 300.660, 300.818, 300.994, 300.941, 300.675/
1425 DATA sst_m10_00 /300.226, 300.558, 300.914, 301.047, 300.645, 299.870, &
1426 299.114, 298.751, 298.875, 299.294, 299.721, 299.989/
1428 DATA sst_m20_m10 /299.547, 299.985, 300.056, 299.676, 298.841, 297.788, &
1429 296.893, 296.491, 296.687, 297.355, 298.220, 298.964/
1431 DATA sst_m30_m20 /297.524, 298.073, 297.897, 297.088, 295.846, 294.520, &
1432 293.525, 293.087, 293.217, 293.951, 295.047, 296.363/
1434 DATA sst_m40_m30 /293.054, 293.765, 293.468, 292.447, 291.128, 289.781, &
1435 288.773, 288.239, 288.203, 288.794, 289.947, 291.553/
1437 DATA sst_m50_m40 /285.052, 285.599, 285.426, 284.681, 283.761, 282.826, &
1438 282.138, 281.730, 281.659, 281.965, 282.768, 283.961/
1440 DATA sst_m60_m50 /277.818, 278.174, 277.991, 277.455, 276.824, 276.229, &
1441 275.817, 275.585, 275.560, 275.687, 276.142, 276.968/
1443 DATA sst_m70_m60 /273.436, 273.793, 273.451, 272.813, 272.349, 272.048, &
1444 271.901, 271.838, 271.845, 271.889, 272.080, 272.607/
1446 DATA sst_m80_m70 /271.579, 271.578, 271.471, 271.407, 271.392, 271.391, &
1447 271.390, 271.391, 271.394, 271.401, 271.422, 271.486/
1449 DATA sst_m90_m80 /271.350, 271.350, 271.350, 271.350, 271.350, 271.350, &
1450 271.350, 271.350, 271.350, 271.350, 271.350, 271.350/
1453 IF (latitude > 80.0)
THEN 1455 ELSEIF (latitude > 70.0)
THEN 1457 ELSEIF (latitude > 60.0)
THEN 1459 ELSEIF (latitude > 50.0)
THEN 1461 ELSEIF (latitude > 40.0)
THEN 1463 ELSEIF (latitude > 30.0)
THEN 1465 ELSEIF (latitude > 20.0)
THEN 1467 ELSEIF (latitude > 10.0)
THEN 1469 ELSEIF (latitude > 0.0)
THEN 1471 ELSEIF (latitude > -10.0)
THEN 1473 ELSEIF (latitude > -20.0)
THEN 1475 ELSEIF (latitude > -30.0)
THEN 1477 ELSEIF (latitude > -40.0)
THEN 1479 ELSEIF (latitude > -50.0)
THEN 1481 ELSEIF (latitude > -60.0)
THEN 1483 ELSEIF (latitude > -70.0)
THEN 1485 ELSEIF (latitude > -80.0)
THEN 1493 IF(mon2 == 13) mon2 = 1
1495 dtref = (sst(mon2) - sst(mon1)) / num_days(mon1)
1498 IF (mon1 == 0) mon1=12
1500 dtref = (sst(mon2) - sst(mon1)) / num_days(mon1)
1503 dtref = dtref * (deltsfc / 24.0)
1518 SUBROUTINE dtzm_point(XT,XZ,DT_COOL,ZC,Z1,Z2,DTZM)
1521 real,
intent(in) :: xt,xz,dt_cool,zc,z1,z2
1522 real,
intent(out) :: dtzm
1524 real,
parameter :: zero = 0.0
1525 real,
parameter :: one = 1.0
1526 real,
parameter :: half = 0.5
1527 real :: dt_warm,dtw,dtc
1532 if ( xt > zero )
then 1533 dt_warm = (xt+xt)/xz
1536 dtw = dt_warm*(one-(z1+z2)/(xz+xz))
1537 elseif ( z1 < xz .and. z2 >= xz )
then 1538 dtw = half*(one-z1/xz)*dt_warm*(xz-z1)/(z2-z1)
1540 elseif ( z1 == z2 )
then 1542 dtw = dt_warm*(one-z1/xz)
1550 if ( zc > zero )
then 1553 dtc = dt_cool*(one-(z1+z2)/(zc+zc))
1554 elseif ( z1 < zc .and. z2 >= zc )
then 1555 dtc = half*(one-z1/zc)*dt_cool*(zc-z1)/(z2-z1)
1557 elseif ( z1 == z2 )
then 1559 dtc = dt_cool*(one-z1/zc)
1590 subroutine tf_thaw_set(tf_ij,mask_ij,itile,jtile,tice,tclm,tf_thaw,nx,ny, &
1591 nset_thaw_s,nset_thaw_i,nset_thaw_c)
1593 real,
dimension(nx*ny),
intent(in) :: tf_ij
1594 integer,
dimension(nx*ny),
intent(in) :: mask_ij
1595 real,
intent(in) :: tice,tclm
1596 integer,
intent(in) :: itile,jtile,nx,ny
1597 real,
intent(out) :: tf_thaw
1598 integer,
intent(inout) :: nset_thaw_s,nset_thaw_i,nset_thaw_c
1600 real,
parameter :: bmiss = -999.0
1601 real,
dimension(nx,ny) :: tf
1602 integer,
dimension(nx,ny) :: mask
1603 integer :: krad,max_search,istart,iend,jstart,jend
1604 integer :: ii,jj,iii,jjj
1609 mask(:,:) = reshape(mask_ij,(/nx,ny/) )
1610 tf(:,:) = reshape(tf_ij,(/nx,ny/) )
1614 do krad = 1, max_search
1616 istart = itile - krad
1618 jstart = jtile - krad
1621 do jj = jstart, jend
1622 do ii = istart, iend
1625 if ((jj == jstart) .or. (jj == jend) .or. &
1626 (ii == istart) .or. (ii == iend))
then 1628 if ((jj >= 1) .and. (jj <= ny))
then 1632 else if (ii >= (nx+1))
then 1643 if (krad <= 2 .and. mask(iii,jjj) == 2) is_ice = .true.
1645 if (mask(iii,jjj) == 0)
then 1646 tf_thaw = tf(iii,jjj)
1647 nset_thaw_s = nset_thaw_s + 1
1648 write(*,
'(a,I4,2F9.3)')
'nset_thaw_s,tf(iii,jjj),tclm : ',nset_thaw_s,tf(iii,jjj),tclm
1660 if ( tf_thaw == bmiss )
then 1663 nset_thaw_i = nset_thaw_i + 1
1664 write(*,
'(a,I4,F9.3)')
'nset_thaw_i,tf_ice : ',nset_thaw_i,tice
1666 tf_thaw = 0.8*tice+0.2*tclm
1667 nset_thaw_c = nset_thaw_c + 1
1668 write(*,
'(a,I4,2F9.3)')
'nset_thaw_c,tf_ice,tclm : ',nset_thaw_c,tice,tclm
1685 integer,
intent(in) :: ij
1687 real,
intent(in) :: tf_thaw
1689 type(nsst_data),
intent(inout) :: nsst
1693 nsst%d_conv(ij) = 0.0
1694 nsst%dt_cool(ij) = 0.0
1696 nsst%qrain(ij) = 0.0
1697 nsst%tref(ij) = tf_thaw
1727 subroutine get_tf_clm(xlats_ij,xlons_ij,ny,nx,iy,im,id,ih,tf_clm,tf_trd)
1732 real,
dimension(nx*ny),
intent(in) :: xlats_ij
1733 real,
dimension(nx*ny),
intent(in) :: xlons_ij
1734 real,
dimension(nx,ny),
intent(out) :: tf_clm
1735 real,
dimension(nx,ny),
intent(out) :: tf_trd
1736 integer,
intent(in) :: iy,im,id,ih,nx,ny
1738 real,
allocatable,
dimension(:,:) :: tf_clm0
1739 real,
allocatable,
dimension(:,:) :: tf_trd0
1740 real,
allocatable,
dimension(:) :: cxlats
1741 real,
allocatable,
dimension(:) :: cxlons
1743 real,
dimension(nx*ny) :: tf_clm_ij
1744 real,
dimension(nx*ny) :: tf_trd_ij
1746 integer :: nxc,nyc,mon1,mon2
1747 character (len=6),
parameter :: fin_tf_clm=
'sstclm' 1756 allocate( tf_clm0(nxc,nyc),tf_trd0(nxc,nyc),cxlats(nyc),cxlons(nxc) )
1760 call get_tf_clm_ta(tf_clm0,tf_trd0,cxlats,cxlons,nyc,nxc,mon1,mon2,wei1,wei2)
1764 if ( nx == nxc .and. ny == nyc )
then 1765 tf_clm(:,:) = tf_clm0(:,:)
1766 tf_trd(:,:) = tf_trd0(:,:)
1770 call intp_tile(tf_clm0, cxlats, cxlons, nyc, nxc, &
1771 tf_clm_ij,xlats_ij,xlons_ij,ny, nx)
1772 call intp_tile(tf_trd0, cxlats, cxlons, nyc, nxc, &
1773 tf_trd_ij,xlats_ij,xlons_ij,ny, nx)
1776 tf_clm(:,:) = reshape(tf_clm_ij, (/nx,ny/) )
1777 tf_trd(:,:) = reshape(tf_trd_ij, (/nx,ny/) )
1796 subroutine get_tf_clm_ta(tf_clm_ta,tf_clm_trend,xlats,xlons,nlat,nlon,mon1,mon2,wei1,wei2)
1801 integer,
intent(in) :: nlat,nlon,mon1,mon2
1802 real,
intent(in) :: wei1,wei2
1804 real,
dimension(nlon,nlat),
intent(out) :: tf_clm_ta,tf_clm_trend
1805 real,
dimension(nlat),
intent(out) :: xlats
1806 real,
dimension(nlon),
intent(out) :: xlons
1809 character (len=6),
parameter :: fin_tf_clm=
'sstclm' 1812 real,
dimension(nlon,nlat) :: tf_clm1,tf_clm2
1822 tf_clm_ta(:,:) = wei1*tf_clm1(:,:)+wei2*tf_clm2(:,:)
1826 tf_clm_trend(:,:) = (tf_clm2(:,:)-tf_clm1(:,:))/120.0
1828 write(*,
'(a,2f9.3)')
'tf_clm_ta, min, max : ',minval(tf_clm_ta),maxval(tf_clm_ta)
1829 write(*,
'(a,2f9.3)')
'tf_clm_trend, min, max : ',minval(tf_clm_trend),maxval(tf_clm_trend)
1844 subroutine get_sal_clm(xlats_ij,xlons_ij,ny,nx,iy,im,id,ih,sal_clm)
1848 real,
dimension(nx*ny),
intent(in) :: xlats_ij
1849 real,
dimension(nx*ny),
intent(in) :: xlons_ij
1850 real,
dimension(nx,ny),
intent(out) :: sal_clm
1851 integer,
intent(in) :: iy,im,id,ih,nx,ny
1853 real,
allocatable,
dimension(:,:) :: sal_clm0
1854 real,
allocatable,
dimension(:) :: cxlats
1855 real,
allocatable,
dimension(:) :: cxlons
1857 real,
dimension(nx*ny) :: sal_clm_ij
1859 integer :: nxc,nyc,mon1,mon2
1860 character (len=6),
parameter :: fin_sal_clm=
'salclm' 1869 allocate( sal_clm0(nxc,nyc),cxlats(nyc),cxlons(nxc) )
1873 call get_sal_clm_ta(sal_clm0,cxlats,cxlons,nyc,nxc,mon1,mon2,wei1,wei2)
1877 if ( nx == nxc .and. ny == nyc )
then 1878 sal_clm(:,:) = sal_clm0(:,:)
1882 call intp_tile(sal_clm0, cxlats, cxlons, nyc, nxc, &
1883 sal_clm_ij,xlats_ij,xlons_ij,ny, nx)
1887 sal_clm(:,:) = reshape(sal_clm_ij, (/nx,ny/) )
1904 subroutine get_sal_clm_ta(sal_clm_ta,xlats,xlons,nlat,nlon,mon1,mon2,wei1,wei2)
1910 integer,
intent(in) :: nlat,nlon,mon1,mon2
1911 real,
intent(in) :: wei1,wei2
1913 real,
dimension(nlon,nlat),
intent(out) :: sal_clm_ta
1914 real,
dimension(nlat),
intent(out) :: xlats
1915 real,
dimension(nlon),
intent(out) :: xlons
1918 character (len=6),
parameter :: fin_sal_clm=
'salclm' 1921 real,
dimension(nlon,nlat) :: sal_clm1,sal_clm2
1931 sal_clm_ta(:,:) = wei1*sal_clm1(:,:)+wei2*sal_clm2(:,:)
1932 write(*,
'(a,2f9.3)')
'sal_clm_ta, min, max : ',minval(sal_clm_ta),maxval(sal_clm_ta)
1949 subroutine intp_tile(tf_lalo,dlats_lalo,dlons_lalo,jdim_lalo,idim_lalo, &
1950 tf_tile,xlats_tile,xlons_tile,jdim_tile,idim_tile)
1957 real,
dimension(idim_lalo,jdim_lalo),
intent(in) :: tf_lalo
1958 real,
dimension(jdim_lalo),
intent(in) :: dlats_lalo
1959 real,
dimension(idim_lalo),
intent(in) :: dlons_lalo
1960 real,
dimension(jdim_tile*idim_tile),
intent(in) :: xlats_tile
1961 real,
dimension(jdim_tile*idim_tile),
intent(in) :: xlons_tile
1962 integer,
intent(in) :: jdim_lalo,idim_lalo,jdim_tile,idim_tile
1963 real,
dimension(jdim_tile*idim_tile),
intent(out) :: tf_tile
1966 real,
parameter :: deg2rad=3.1415926/180.0
1967 real,
dimension(jdim_lalo) :: xlats_lalo
1968 real,
dimension(idim_lalo) :: xlons_lalo
1970 integer :: itile,jtile
1972 integer :: ilalo,jlalo,ilalop1,jlalop1
1974 integer,
allocatable,
dimension(:,:) :: id1,id2,jdc
1975 real,
allocatable,
dimension(:,:,:) :: agrid,s2c
1978 print*,
'interpolate from lat/lon grids to any one grid with known lat/lon' 1980 xlats_lalo = dlats_lalo*deg2rad
1981 xlons_lalo = dlons_lalo*deg2rad
1983 allocate(agrid(idim_tile,jdim_tile,2))
1984 agrid(:,:,1) = reshape(xlons_tile, (/idim_tile,jdim_tile/) )
1985 agrid(:,:,2) = reshape(xlats_tile, (/idim_tile,jdim_tile/) )
1986 agrid = agrid*deg2rad
1988 allocate(id1(idim_tile,jdim_tile))
1989 allocate(id2(idim_tile,jdim_tile))
1990 allocate(jdc(idim_tile,jdim_tile))
1991 allocate(s2c(idim_tile,jdim_tile,4))
1999 call remap_coef( 1, idim_tile, 1, jdim_tile, idim_lalo, jdim_lalo, &
2000 xlons_lalo, xlats_lalo, id1, id2, jdc, s2c, agrid )
2002 do ij = 1, jdim_tile*idim_tile
2004 jtile = (ij-1)/idim_tile + 1
2005 itile = mod(ij,idim_tile)
2006 if (itile==0) itile = idim_tile
2008 ilalo = id1(itile,jtile)
2009 ilalop1 = id2(itile,jtile)
2010 jlalo = jdc(itile,jtile)
2011 jlalop1 = jdc(itile,jtile) + 1
2013 wsum = s2c(itile,jtile,1) + s2c(itile,jtile,2) + &
2014 s2c(itile,jtile,3) + s2c(itile,jtile,4)
2016 tf_tile(ij) = ( s2c(itile,jtile,1)*tf_lalo(ilalo,jlalo) + &
2017 s2c(itile,jtile,2)*tf_lalo(ilalop1,jlalo) + &
2018 s2c(itile,jtile,3)*tf_lalo(ilalop1,jlalop1) + &
2019 s2c(itile,jtile,4)*tf_lalo(ilalo,jlalop1) )/wsum
2022 deallocate(id1, id2, jdc, s2c)
2038 subroutine get_tim_wei(iy,im,id,ih,mon1,mon2,wei1,wei2)
2042 integer,
intent(in) :: iy,im,id,ih
2044 integer,
intent(out) :: mon1,mon2
2045 real,
intent(out) :: wei1,wei2
2049 integer :: mon,monend,monm,monp,jdow,jdoy,jday
2054 real,
dimension(13) :: dayhf
2055 data dayhf/15.5,45.0,74.5,105.0,135.5,166.0,196.5,227.5,258.0,288.5,319.0,349.5,380.5/
2068 call w3doxdat(jda,jdow,jdoy,jday)
2069 rjday=jdoy+jda(5)/24.
2070 if(rjday.lt.dayhf(1)) rjday=rjday+365.
2076 if( rjday >= dayhf(monm) .and. rjday < dayhf(monp) )
then 2083 print *,
'FATAL ERROR in get_tim_wei, wrong rjday',rjday
2087 wei1 = (dayhf(mon2)-rjday)/(dayhf(mon2)-dayhf(mon1))
2088 wei2 = (rjday-dayhf(mon1))/(dayhf(mon2)-dayhf(mon1))
2090 if( mon2 == 13 ) mon2=1
2092 write(*,
'(a,2i4,3f9.3)')
'mon1,mon2,rjday,wei1,wei2=',mon1,mon2,rjday,wei1,wei2
2105 real function tfreez(salinity)
2111 parameter(a1 = -0.0575)
2112 parameter(a2 = 1.710523e-3)
2113 parameter(a3 = -2.154996e-4)
2115 IF (salinity .LT. 0.)
THEN 2120 tfreez = sal*(a1+a2*sqrt(sal)+a3*sal)
integer function num_parthds()
Return the number of omp threads.
subroutine, public read_data(LSOIL, LENSFC, DO_NSST, DO_SNO_INC_JEDI, DO_SOI_INC_JEDI, INC_FILE, IS_NOAHMP, TSFFCS, SMCFCS, SWEFCS, STCFCS, TG3FCS, ZORFCS, CVFCS, CVBFCS, CVTFCS, ALBFCS, VEGFCS, SLIFCS, CNPFCS, F10M, VETFCS, SOTFCS, ALFFCS, USTAR, FMM, FHH, SIHFCS, SICFCS, SITFCS, TPRCP, SRFLAG, SNDFCS, VMNFCS, VMXFCS, SLCFCS, STCINC, SLCINC, SLPFCS, ABSFCS, T2M, Q2M, SLMASK, ZSOIL, NSST)
Read the first guess surface records and nsst records (if selected) for a single cubed-sphere tile...
subroutine, public get_tf_clm_dim(file_sst, mlat_sst, mlon_sst)
Get the i/j dimensions of RTG SST climatology file.
subroutine, public read_tf_clim_grb(file_sst, sst, rlats_sst, rlons_sst, mlat_sst, mlon_sst, mon)
Read a GRIB1 sst climatological analysis file.
integer, public idim_gaus
'i' dimension of GSI gaussian grid.
subroutine get_tf_clm(xlats_ij, xlons_ij, ny, nx, iy, im, id, ih, tf_clm, tf_trd)
Get the sst climatology at the valid time and on the target grid.
subroutine, public read_salclm_gfs_nc(filename, sal, xlats, xlons, nlat, nlon, itime)
Read the woa05 salinity monthly climatology file.
subroutine get_tf_clm_ta(tf_clm_ta, tf_clm_trend, xlats, xlons, nlat, nlon, mon1, mon2, wei1, wei2)
Get the reference temperature/sst climatology and its trend at analysis time.
subroutine, public remap_coef(is, ie, js, je, im, jm, lon, lat, id1, id2, jdc, s2c, agrid)
Generate the weights and index of the grids used in the bilinear interpolation.
integer, dimension(:,:), allocatable, public slmsk_gaus
GSI land mask on the gaussian grid.
subroutine dtzm_point(XT, XZ, DT_COOL, ZC, Z1, Z2, DTZM)
Compute the vertical mean of the NSST t-profile.
subroutine get_sal_clm_ta(sal_clm_ta, xlats, xlons, nlat, nlon, mon1, mon2, wei1, wei2)
Get climatological salinity at the analysis time.
Holds machine dependent constants for global_cycle.
subroutine sfcdrv(LUGB, IDIM, JDIM, LENSFC, LSOIL, DELTSFC, IY, IM, ID, IH, FH, IALB, USE_UFO, DO_NSST, DO_SFCCYCLE, DO_LNDINC, FRAC_GRID, ZSEA1, ZSEA2, ISOT, IVEGSRC, MYRANK)
Driver routine for updating the surface fields.
subroutine tf_thaw_set(tf_ij, mask_ij, itile, jtile, tice, tclm, tf_thaw, nx, ny, nset_thaw_s, nset_thaw_i, nset_thaw_c)
Set the background reference temperature (tf) for points where the ice has just melted.
subroutine, public get_dim_nc(filename, nlat, nlon)
Get the i/j dimensions of the data from a NetCDF file.
subroutine, public read_gsi_data(GSI_FILE, FILE_TYPE, LSOIL)
Read increment file from the GSI containing either the foundation temperature increments and mask...
subroutine intp_tile(tf_lalo, dlats_lalo, dlons_lalo, jdim_lalo, idim_lalo, tf_tile, xlats_tile, xlons_tile, jdim_tile, idim_tile)
Interpolate lon/lat grid data to the fv3 native grid (tf_lalo => tf_tile).
subroutine climo_trend(LATITUDE, MON, DAY, DELTSFC, DTREF)
If the tile point is an isolated water point that has no corresponding gsi water point, then tref is updated using the rtg sst climo trend.
real, dimension(:,:), allocatable, public dtref_gaus
GSI foundation temperature increment on the gaussian grid.
subroutine, public write_data(lensfc, idim, jdim, lsoil, do_nsst, inc_file, nsst, slifcs, tsffcs, vegfcs, swefcs, tg3fcs, zorfcs, albfcs, alffcs, cnpfcs, f10m, t2m, q2m, vetfcs, sotfcs, ustar, fmm, fhh, sicfcs, sihfcs, sitfcs, tprcp, srflag, swdfcs, vmnfcs, vmxfcs, slpfcs, absfcs, slcfcs, smcfcs, stcfcs, stcinc, slcinc)
Update surface records - and nsst records if selected - on a single cubed-sphere tile to a pre-existi...
Module containing utility routines.
This module contains routines that read and write data.
program sfc_drv
Stand alone surface/NSST cycle driver for the cubed-sphere grid.
subroutine adjust_nsst(RLA, RLO, SLMSK_TILE, SLMSK_FG_TILE, SKINT_TILE, SICET_TILE, sice_tile, sice_fg_tile, SOILT_TILE, NSST, LENSFC, LSOIL, IDIM, JDIM, ZSEA1, ZSEA2, tf_clm_tile, tf_trd_tile, sal_clm_tile, LANDFRAC, FRAC_GRID)
Read in gsi file with the updated reference temperature increments (on the gaussian grid)...
real function tfreez(salinity)
Compute the freezing point of water as a function of salinity.
integer, public jdim_gaus
'j' dimension of GSI gaussian grid.
subroutine get_tim_wei(iy, im, id, ih, mon1, mon2, wei1, wei2)
For a given date, determine the bounding months and the linear time interpolation weights...
subroutine get_sal_clm(xlats_ij, xlons_ij, ny, nx, iy, im, id, ih, sal_clm)
Get salinity climatology at the valid time on the target grid.
subroutine nsst_water_reset(nsst, ij, tf_thaw)
If the first guess was sea ice, but the analysis is open water, reset all nsst variables.
subroutine, public read_lat_lon_orog(RLA, RLO, OROG, OROG_UF, TILE_NUM, IDIM, JDIM, IJDIM, LANDFRAC)
Read latitude and longitude for the cubed-sphere tile from the 'grid' file.