Skip to content

Commit

Permalink
Merge branch 'develop' into feature/msienkie/acftBC_cleanup_tails
Browse files Browse the repository at this point in the history
  • Loading branch information
gmao-msienkie committed Dec 6, 2021
2 parents d61fc58 + b3a880f commit 028dbde
Show file tree
Hide file tree
Showing 4 changed files with 175 additions and 35 deletions.
2 changes: 1 addition & 1 deletion src/Applications/NCEP_Etc/NCEP_bkgecov/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -25,7 +25,7 @@ esma_add_library (${this}
DEPENDENCIES NCEP_bacio_r4i4 NCEP_w3_r4i4 NCEP_sp_r4i4 GMAO_gfio_r8 GMAO_hermes GMAO_transf ${MKL_LIBRARIES}
)

target_compile_definitions(${this} PRIVATE _LAPACK_ gmao_intf)
add_definitions (-D_LAPACK_ -Dgmao_intf)

ecbuild_add_executable(TARGET calcstats.x SOURCES statsmain.F90 LIBS ${this} ${MKL_LIBRARIES})
ecbuild_add_executable(TARGET write_berror_global.x SOURCES write_berror_global.f90 LIBS ${this} ${MKL_LIBRARIES})
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,7 @@
! note: biasrm=.false. means do bias correction - go figure!
fnm0=>>>NMODES<<<,vectype=5,biasrm=.false.,laddoz=.false.,
LBAL=.true.,
nreaders=8,
nreaders=16,
readperts=.true.,
! Debug only:
! calchrzscl=.false.,
Expand Down
65 changes: 38 additions & 27 deletions src/Applications/NCEP_Etc/NCEP_bkgecov/ut_gen_berrcov.j
Original file line number Diff line number Diff line change
@@ -1,10 +1,10 @@
#!/bin/csh -x
#SBATCH --account=g0613
#SBATCH --constraint=hasw
#SBATCH --ntasks=1
#SBATCH --constraint=sky
#_SBATCH --ntasks=8
#_SBATCH --ntasks=24
#_SBATCH --ntasks=96
#_SBATCH --ntasks=384
#SBATCH --ntasks=672
#_SBATCH --ntasks-per-node=24
#SBATCH --time=12:00:00
#SBATCH --qos=dastest
Expand All @@ -22,7 +22,7 @@

setenv DRYRUN #echo
setenv SKIPSETTING 0
setenv MYNCPUS 4
setenv MYNCPUS 8
# Set initial time and number of samples required
setenv EXPID x0039_p5_REPLAY_L132
setenv EXPID x41Lrt
Expand All @@ -32,33 +32,38 @@
setenv EXPID f525_p7_fp
setenv EXPID f522_fp
setenv EXPID f525_fp
setenv EXPID f5271_fp
#setenv FVHOME /home/dao_ops/$EXPID
#setenv ARCROOT /home/dao_ops/$EXPID/run/.../archive/prog
#setenv FVHOME /discover/nobackup/projects/gmao/obsdev/rtodling/$EXPID
#setenv FVHOME /discover/nobackup/projects/gmao/obsdev/rtodling/x0041
setenv FVHOME /discover/nobackup/projects/gmao/dadev/rtodling/x0043rt
setenv FVHOME /discover/nobackup/projects/gmao/dadev/rtodling/prePP
setenv FVROOT `cat $FVHOME/.FVROOT`
setenv PLAINDIR 0
setenv ARCROOT /archive/u/$user/$EXPID/prog
if ( $EXPID == "x0041" ) then
setenv ARCROOT /archive/u/dao_it/$EXPID/prog
endif
if ( $EXPID == "f521_fp" || $EXPID == "f522_fp" || $EXPID == "f525_fp" || $EXPID == "f525_p5_fp" || $EXPID == "f525_p7_fp" ) then
if ( $EXPID == "f521_fp" || $EXPID == "f522_fp" || $EXPID == "f525_fp" || $EXPID == "f525_p5_fp" || $EXPID == "f525_p7_fp" || $EXPID == "f5271_fp" ) then
setenv ARCROOT /home/dao_ops/$EXPID/run/.../archive/prog
endif
if ($EXPID == "x0039_p5_REPLAY_L132") then
setenv ARCROOT /gpfsm/dnb78s1/projects/p18/ltakacs/REPLAY_Experiments/x0039_p5_REPLAY_L132/forecasts/Regular_RPLFRQ:7200_ANA:x0039_p5
setenv PLAINDIR 1
endif
setenv DMGET 0
setenv GET_SET 1
setenv DMGET 0
setenv GET_SET 0

setenv DO_ACQUIRE 1
setenv GEN_NMCDIFFS 1
setenv DO_ACQUIRE 0
setenv GEN_NMCDIFFS 0

setenv GET_BERROR 0
setenv GET_BERROR 1
setenv BERROR_NMODES 25

set JCAP = 512
set NLON = 576
set NLAT = 361

# Basic settings (weak dependency on version of DAS)
# --------------------------------------------------
set path = ( . $FVROOT/bin $path )
Expand All @@ -72,9 +77,9 @@ setenv FCSTWORK /discover/nobackup/projects/gmao/obsdev/$user/fcst4berrcov.$EXPI
setenv FCSTWORK /discover/nobackup/projects/gmao/dadev/$user/fcst4berrcov.$EXPID

if ($?I_MPI_ROOT ) then
setenv MPIRUN_CALCSTATS "mpirun -np 32 calcstats.x"
setenv MPIRUN_CALCSTATS "mpirun -np 96 calcstats.x"
else
setenv MPIRUN_CALCSTATS "mpirun_exec -np 32 calcstats.x"
setenv MPIRUN_CALCSTATS "mpirun_exec -np 96 calcstats.x"
endif
#
# initial verification date and number of samples (will get samples ahead of initial date)
Expand All @@ -85,6 +90,8 @@ set vnymd0 = 20181101
set vnymd0 = 20200203
set vnymd0 = 20200401
set vnymd0 = 20200301
set vnymd0 = 20210203 # start
set vnymd0 = 20210305
#if ( $EXPID == "f522_fp") then
# set vnymd0 = 20190304
# @ nsamples = 350
Expand All @@ -94,8 +101,8 @@ set vnymd0 = 20200301
#endif
set vnhms0 = 000000
@ nsamples = 31
setenv FCSTWRK $FCSTWORK.$vnymd0.$vnhms0
mkdir -p $FCSTWRK/BERROR.WORK
#setenv FCSTWRK $FCSTWORK.$vnymd0.$vnhms0
setenv FCSTWRK $FCSTWORK.all
set diren = `dirname $FCSTWRK`
set spool = "-s $diren/spool "

Expand All @@ -106,8 +113,8 @@ if ( $GET_SET ) then
set foffset_hr = 3 # offset in initial time from synoptic hour
set nmc_hrv = 24 # begin time of NMC method / verification time
set nmc_hrf = 48 # end time of NMC method / fcst time
set vmn = 00 # either black or 00
set vmn = # either black or 00
set vmn = 00 # either black or 00

# ----------------------------------
# No user change from this part down
Expand Down Expand Up @@ -181,8 +188,13 @@ if ( $GET_SET ) then
set inidate = ( `tick $inidate $gap_sc` )
end

set lst = `cat $acqfile`
if (-e missing.txt ) /bin/rm missing.txt
touch missing.txt
foreach fn ($lst)
if (! -e $fn ) echo $fn >> missing.txt
end
if ( $DMGET ) then
set lst = `cat $acqfile`
echo $lst
dmget $lst
exit
Expand Down Expand Up @@ -295,32 +307,31 @@ endif # GEN_NMCDIFFS
if ( $GET_BERROR ) then

if ( ! -d $FCSTWRK/BERROR.WORK ) mkdir -p $FCSTWRK/BERROR.WORK
cd $FCSTWRK/BERROR.WORK

# Get positioned and set namelist parameters
# ------------------------------------------
cd $FCSTWRK/BERROR.WORK/

mkdir samples
cd samples
ln -sf $FCSTWRK/NMC48m24/*f48m24*nc4 .
cd -

# Get set to run berror stats code
# --------------------------------
if (-e infiles ) /bin/rm infiles
cat $FCSTWRK/shrt_fcst.txt >> infiles
cat $FCSTWRK/long_fcst.txt >> infiles
# cat $FCSTWRK/shrt_fcst.txt >> infiles
# cat $FCSTWRK/long_fcst.txt >> infiles
ls samples/*f48m24*nc4 > infiles
ln -sf infiles fort.10

# just for the record ...
ls -lrt samples

# Figure out dimensions of forecast files (assumes that checking one suffices)
# ----------------------------------------------------------------------------
set fcst_files = (`cat samples`)
set fcst_files = (`cat infiles`)
set fcst_res = (`getgfiodim.x $fcst_files[1]`)

# wire for now
set JCAP = 512
set NSIG = $fcst_res[3]
set NLON = 576
set NLAT = 361

# Prepare resource files
set this_param = $FVROOT/etc/berror_stats.nml.tmpl
Expand Down
141 changes: 135 additions & 6 deletions src/Applications/NCEP_Etc/NCEP_bkgecov/write_berror_global.f90
Original file line number Diff line number Diff line change
Expand Up @@ -66,19 +66,29 @@ program write_berror_global
call berror_read_(ivars)
endif

if (mlat/=ivars%nlat.or.mlon/=ivars%nlon) then
print *, 'cannot interpolate yet ...'
if (mlat/=ivars%nlat.and.mlon/=ivars%nlon.and.msig/=ivars%nsig) then
print *, 'cannot interpolate all three dims at one, try horz than vert ...'
call exit(1)
endif
if (mlat/=ivars%nlat.or.mlon/=ivars%nlon) then
write(6,'(a)') ' Horizontally interpolating error covariance fields ...'
call nc_berror_vars_init(xvars,ilon,ilat,isig)
call nc_berror_vars_copy(ivars,xvars)
call nc_berror_vars_final(ivars)
call nc_berror_vars_init(ivars,mlon,mlat,isig)
call hinterp_berror_vars_(xvars,ivars)
write(6,'(a)') ' Finish horizontal interpolation.'
call nc_berror_vars_final(xvars)
endif
if (msig/=ivars%nsig) then
write(6,'(a)') ' Interpolating error covariance fields ...'
write(6,'(a)') ' Vertically interpolating error covariance fields ...'
call nc_berror_vars_init(xvars,ilon,ilat,isig)
call nc_berror_vars_copy(ivars,xvars)
call nc_berror_vars_final(ivars)
call nc_berror_vars_init(ivars,ilon,ilat,msig)
call nc_berror_vars_copy(xvars,ivars)
call nc_berror_vars_copy(xvars,ivars) ! copy lat/lon fields
call vinterp_berror_vars_(xvars,ivars)
write(6,'(a)') ' Finish interpolation.'
write(6,'(a)') ' Finish vertical interpolation.'
call nc_berror_vars_final(xvars)
endif
call berror_write_(ivars,merra2current)
Expand Down Expand Up @@ -602,10 +612,129 @@ subroutine vinterp_berror_vars_(ivars,ovars)
call spline( plevi, plevo, ivars%cvln (j,:), ovars%cvln (j,:) )
enddo
deallocate(aux)

deallocate(plevi)
deallocate(plevo)

end subroutine vinterp_berror_vars_

subroutine hinterp_berror_vars_(ivars,ovars)

use m_spline, only: spline
use m_set_eta, only: set_eta
use m_set_eta, only: get_ref_plevs
use m_const, only: pstd
implicit none

type(nc_berror_vars) ivars
type(nc_berror_vars) ovars

real(4),allocatable,dimension(:,:) :: aux
real(4),allocatable,dimension(:) :: lati,lato
real(4),allocatable,dimension(:) :: loni,lono
integer i,j,k,k2
real dlon,dlat

if( ivars%nsig/=ovars%nsig ) then
print *, 'hinterp_berror_vars_: error, nsig must equal'
call exit(1)
endif

! Input levels
! ------------
allocate(lati(ivars%nlat))
allocate(loni(ivars%nlon))
dlat = 180./(ivars%nlat-1)
do j=1,ivars%nlat
lati(j) = -90.0 + (j-1.0)*dlat
enddo
dlon = 360./ivars%nlon
do i=1,ivars%nlon
loni(i) = i*dlon ! GSI def
enddo

! Output levels
! -------------
allocate(lato(ovars%nlat))
allocate(lono(ovars%nlon))
dlat = 180./(ovars%nlat-1)
do j=1,ovars%nlat
lato(j) = -90.0 + (j-1.0)*dlat
enddo
dlon = 360./ovars%nlon
do i=1,ovars%nlon
lono(i) = i*dlon ! GSI def
enddo


do k=1,ivars%nsig ! very, very parallelizable

do k2=1,ivars%nsig
call spline( lati, lato, ivars%tcon(:,k2,k), ovars%tcon(:,k2,k) )
enddo

call spline( lati, lato, ivars%vpcon (:,k), ovars%vpcon (:,k) )
call spline( lati, lato, ivars%pscon (:,k), ovars%pscon (:,k) )
call spline( lati, lato, ivars%sfvar (:,k), ovars%sfvar (:,k) )
call spline( lati, lato, ivars%sfhln (:,k), ovars%sfhln (:,k) )
call spline( lati, lato, ivars%sfvln (:,k), ovars%sfvln (:,k) )
call spline( lati, lato, ivars%vpvar (:,k), ovars%vpvar (:,k) )
call spline( lati, lato, ivars%vphln (:,k), ovars%vphln (:,k) )
call spline( lati, lato, ivars%vpvln (:,k), ovars%vpvln (:,k) )
call spline( lati, lato, ivars%tvar (:,k), ovars%tvar (:,k) )
call spline( lati, lato, ivars%thln (:,k), ovars%thln (:,k) )
call spline( lati, lato, ivars%tvln (:,k), ovars%tvln (:,k) )
call spline( lati, lato, ivars%qvar (:,k), ovars%qvar (:,k) )
call spline( lati, lato, ivars%nrhvar(:,k), ovars%nrhvar(:,k) )
call spline( lati, lato, ivars%qhln (:,k), ovars%qhln (:,k) )
call spline( lati, lato, ivars%qvln (:,k), ovars%qvln (:,k) )
call spline( lati, lato, ivars%qivar (:,k), ovars%qivar (:,k) )
call spline( lati, lato, ivars%qihln (:,k), ovars%qihln (:,k) )
call spline( lati, lato, ivars%qivln (:,k), ovars%qivln (:,k) )
call spline( lati, lato, ivars%qlvar (:,k), ovars%qlvar (:,k) )
call spline( lati, lato, ivars%qlhln (:,k), ovars%qlhln (:,k) )
call spline( lati, lato, ivars%qlvln (:,k), ovars%qlvln (:,k) )
call spline( lati, lato, ivars%qrvar (:,k), ovars%qrvar (:,k) )
call spline( lati, lato, ivars%qrhln (:,k), ovars%qrhln (:,k) )
call spline( lati, lato, ivars%qrvln (:,k), ovars%qrvln (:,k) )
call spline( lati, lato, ivars%qsvar (:,k), ovars%qsvar (:,k) )
call spline( lati, lato, ivars%qshln (:,k), ovars%qshln (:,k) )
call spline( lati, lato, ivars%qsvln (:,k), ovars%qsvln (:,k) )
call spline( lati, lato, ivars%ozvar (:,k), ovars%ozvar (:,k) )
call spline( lati, lato, ivars%ozhln (:,k), ovars%ozhln (:,k) )
call spline( lati, lato, ivars%ozvln (:,k), ovars%ozvln (:,k) )
call spline( lati, lato, ivars%cvar (:,k), ovars%cvar (:,k) )
call spline( lati, lato, ivars%chln (:,k), ovars%chln (:,k) )
call spline( lati, lato, ivars%cvln (:,k), ovars%cvln (:,k) )
enddo
call spline( lati, lato, ivars%psvar, ovars%psvar )
call spline( lati, lato, ivars%pshln, ovars%pshln )

! Now handle horizontal 2d fields
allocate(aux(ovars%nlat,ivars%nlon))

! varsst ...
do i=1,ivars%nlon
call spline( lati, lato, ivars%varsst(:,i), aux(:,i) )
enddo
do j=1,ovars%nlat
call spline( loni, lono, aux(j,:), ovars%varsst(j,:) )
enddo
! corlsst
do i=1,ivars%nlon
call spline( lati, lato, ivars%corlsst(:,i), aux(:,i) )
enddo
do j=1,ovars%nlat
call spline( loni, lono, aux(j,:), ovars%corlsst(j,:) )
enddo

deallocate(aux)

deallocate(lati,loni)
deallocate(lato,lono)


end subroutine hinterp_berror_vars_

subroutine be_write_nc_(fname,ivars)

use m_set_eta, only: set_eta
Expand Down

0 comments on commit 028dbde

Please sign in to comment.