SUBROUTINE read_AssPar (model, inp, out, Lwrite) ! !git $Id$ !svn $Id: read_asspar.F 1151 2023-02-09 03:08:53Z arango $ !================================================== Hernan G. Arango === ! Copyright (c) 2002-2023 The ROMS/TOMS Group ! ! Licensed under a MIT/X style license ! ! See License_ROMS.md ! !======================================================================= ! ! ! This routine reads and reports input data assimilation parameters. ! ! ! !======================================================================= ! USE mod_param USE mod_parallel USE mod_fourdvar USE mod_iounits USE mod_ncparam USE mod_scalars ! USE inp_decode_mod ! USE strings_mod, ONLY : FoundError USE strings_mod, ONLY : uppercase ! implicit none ! ! Imported variable declarations ! logical, intent(in) :: Lwrite ! integer, intent(in) :: model, inp, out ! ! Local variable declarations. ! logical :: Lvalue(4) integer :: Mval, Npts, Nval integer :: i, ib, igrid, itrc, k, ng, status integer :: Cdim, Clen, Rdim integer :: Ivalue(1) integer :: Nfiles(Ngrids) logical, dimension(MT) :: Ltracer real(r8), dimension(MT,Ngrids) :: Rtracer real(dp), dimension(nRval) :: Rval real(r8) :: Rvalue(1) character (len=1 ), parameter :: blank = ' ' character (len=40 ) :: KeyWord character (len=50 ) :: label character (len=256) :: fname, line character (len=256), dimension(nCval) :: Cval character (len=*), parameter :: MyFile = & & "ROMS/Utility/read_asspar.F" ! !----------------------------------------------------------------------- ! Initialize. !----------------------------------------------------------------------- ! igrid=1 Nfiles(1:Ngrids)=0 DO i=1,LEN(label) label(i:i)=blank END DO Cdim=SIZE(Cval,1) Clen=LEN(Cval(1)) Rdim=SIZE(Rval,1) ! !----------------------------------------------------------------------- ! Read in 4DVAR assimilation parameters. Then, load input data into ! module. !----------------------------------------------------------------------- ! DO WHILE (.TRUE.) READ (inp,'(a)',ERR=10,END=20) line status=decode_line(line, KeyWord, Nval, Cval, Rval) IF (status.gt.0) THEN SELECT CASE (TRIM(KeyWord)) CASE ('dTdz_min') Npts=load_r(Nval, Rval, Ngrids, dTdz_min) CASE ('ml_depth') Npts=load_r(Nval, Rval, Ngrids, ml_depth) DO ng=1,Ngrids ml_depth(ng)=ABS(ml_depth(ng)) END DO CASE ('LNM_depth') Npts=load_r(Nval, Rval, Ngrids, LNM_depth) DO ng=1,Ngrids LNM_depth(ng)=ABS(LNM_depth(ng)) END DO CASE ('LNM_flag') Npts=load_i(Nval, Rval, 1, Ivalue) LNM_flag=Ivalue(1) CASE ('GradErr') Npts=load_r(Nval, Rval, 1, Rvalue) GradErr=Rvalue(1) CASE ('HevecErr') Npts=load_r(Nval, Rval, 1, Rvalue) HevecErr=Rvalue(1) CASE ('LhessianEV') Npts=load_l(Nval, Cval, 1, Lvalue) LhessianEV=Lvalue(1) CASE ('LhotStart') Npts=load_l(Nval, Cval, 1, Lvalue) LhotStart=Lvalue(1) CASE ('Lprecond') Npts=load_l(Nval, Cval, 1, Lvalue) Lprecond=Lvalue(1) IF ( LhessianEV.and.Lprecond ) THEN LhessianEV=.FALSE. END IF CASE ('Lritz') Npts=load_l(Nval, Cval, 1, Lvalue) Lritz=Lvalue(1) CASE ('NritzEV') Npts=load_i(Nval, Rval, 1, Ivalue) NritzEV=Ivalue(1) CASE ('NpostI') Npts=load_i(Nval, Rval, 1, Ivalue) NpostI=Ivalue(1) CASE ('Nimpact') Npts=load_i(Nval, Rval, 1, Ivalue) Nimpact=Ivalue(1) CASE ('NextraObs') Npts=load_i(Nval, Rval, 1, Ivalue) NextraObs=Ivalue(1) IF (NextraObs.gt.0) THEN IF (.not.allocated(ExtraIndex)) THEN allocate ( ExtraIndex(NextraObs) ) END IF IF (.not.allocated(ExtraName)) THEN allocate ( ExtraName(NextraObs) ) END IF END IF CASE ('ExtraIndex') IF (NextraObs.gt.0) THEN Npts=load_i(Nval, Rval, NextraObs, ExtraIndex) END IF CASE ('ExtraName') IF (NextraObs.gt.0) THEN ExtraName(Nval)=TRIM(Cval(Nval)) END IF CASE ('tl_M2diff') Npts=load_r(Nval, Rval, Ngrids, tl_M2diff) CASE ('tl_M3diff') Npts=load_r(Nval, Rval, Ngrids, tl_M3diff) CASE ('tl_Tdiff') Npts=load_r(Nval, Rval, MT, Ngrids, Rtracer) DO ng=1,Ngrids DO itrc=1,NT(ng) tl_Tdiff(itrc,ng)=Rtracer(itrc,ng) END DO END DO CASE ('LdefNRM') Npts=load_l(Nval, Cval, 4, Ngrids, LdefNRM) CASE ('LwrtNRM') Npts=load_l(Nval, Cval, 4, Ngrids, LwrtNRM) CASE ('CnormM(isFsur)') IF (isFsur.eq.0) THEN IF (Master) WRITE (out,210) 'isFsur' exit_flag=5 RETURN END IF Npts=load_l(Nval, Cval, 1, Lvalue) Cnorm(2,isFsur)=Lvalue(1) CASE ('CnormM(isUbar)') IF (isUbar.eq.0) THEN IF (Master) WRITE (out,210) 'isUbar' exit_flag=5 RETURN END IF Npts=load_l(Nval, Cval, 1, Lvalue) Cnorm(2,isUbar)=Lvalue(1) CASE ('CnormM(isVbar)') IF (isVbar.eq.0) THEN IF (Master) WRITE (out,210) 'isVbar' exit_flag=5 RETURN END IF Npts=load_l(Nval, Cval, 1, Lvalue) Cnorm(2,isVbar)=Lvalue(1) CASE ('CnormM(isUvel)') IF (isUvel.eq.0) THEN IF (Master) WRITE (out,210) 'isUvel' exit_flag=5 RETURN END IF Npts=load_l(Nval, Cval, 1, Lvalue) Cnorm(2,isUvel)=Lvalue(1) CASE ('CnormM(isVvel)') IF (isVvel.eq.0) THEN IF (Master) WRITE (out,210) 'isVvel' exit_flag=5 RETURN END IF Npts=load_l(Nval, Cval, 1, Lvalue) Cnorm(2,isVvel)=Lvalue(1) CASE ('CnormM(isTvar)') IF (MAXVAL(isTvar).eq.0) THEN IF (Master) WRITE (out,210) 'isTvar' exit_flag=5 RETURN END IF Npts=load_l(Nval, Cval, MT, Ltracer) DO itrc=1,MT i=isTvar(itrc) Cnorm(2,i)=Ltracer(itrc) END DO CASE ('CnormI(isFsur)') Npts=load_l(Nval, Cval, 1, Lvalue) Cnorm(1,isFsur)=Lvalue(1) CASE ('CnormI(isUbar)') Npts=load_l(Nval, Cval, 1, Lvalue) Cnorm(1,isUbar)=Lvalue(1) CASE ('CnormI(isVbar)') Npts=load_l(Nval, Cval, 1, Lvalue) Cnorm(1,isVbar)=Lvalue(1) CASE ('CnormI(isUvel)') Npts=load_l(Nval, Cval, 1, Lvalue) Cnorm(1,isUvel)=Lvalue(1) CASE ('CnormI(isVvel)') Npts=load_l(Nval, Cval, 1, Lvalue) Cnorm(1,isVvel)=Lvalue(1) CASE ('CnormI(isTvar)') Npts=load_l(Nval, Cval, MT, Ltracer) DO itrc=1,MT i=isTvar(itrc) Cnorm(1,i)=Ltracer(itrc) END DO CASE ('balance(isSalt)') Npts=load_l(Nval, Cval, 1, Lvalue) balance(isTvar(isalt))=Lvalue(1) CASE ('balance(isFsur)') Npts=load_l(Nval, Cval, 1, Lvalue) balance(isFsur)=Lvalue(1) CASE ('balance(isVbar)') Npts=load_l(Nval, Cval, 1, Lvalue) balance(isVbar)=Lvalue(1) CASE ('balance(isVvel)') Npts=load_l(Nval, Cval, 1, Lvalue) balance(isVvel)=Lvalue(1) CASE ('Nmethod') Npts=load_i(Nval, Rval, Ngrids, Nmethod) CASE ('Rscheme') Npts=load_i(Nval, Rval, Ngrids, Rscheme) CASE ('Nrandom') Npts=load_i(Nval, Rval, 1, Ivalue) Nrandom=Ivalue(1) CASE ('Hgamma') Npts=load_r(Nval, Rval, 4, Hgamma) CASE ('Vgamma') Npts=load_r(Nval, Rval, 4, Vgamma) CASE ('HdecayM(isFsur)') Npts=load_r(Nval, Rval, Ngrids, Hdecay(2,isFsur,:)) CASE ('HdecayM(isUbar)') Npts=load_r(Nval, Rval, Ngrids, Hdecay(2,isUbar,:)) CASE ('HdecayM(isVbar)') Npts=load_r(Nval, Rval, Ngrids, Hdecay(2,isVbar,:)) CASE ('HdecayM(isUvel)') Npts=load_r(Nval, Rval, Ngrids, Hdecay(2,isUvel,:)) CASE ('HdecayM(isVvel)') Npts=load_r(Nval, Rval, Ngrids, Hdecay(2,isVvel,:)) CASE ('HdecayM(isTvar)') Npts=load_r(Nval, Rval, MT, Ngrids, Rtracer) DO ng=1,Ngrids DO itrc=1,NT(ng) Hdecay(2,isTvar(itrc),ng)=Rtracer(itrc,ng) END DO END DO CASE ('VdecayM(isUvel)') Npts=load_r(Nval, Rval, Ngrids, Vdecay(2,isUvel,:)) CASE ('VdecayM(isVvel)') Npts=load_r(Nval, Rval, Ngrids, Vdecay(2,isVvel,:)) CASE ('VdecayM(isTvar)') Npts=load_r(Nval, Rval, MT, Ngrids, Rtracer) DO ng=1,Ngrids DO itrc=1,NT(ng) Vdecay(2,isTvar(itrc),ng)=Rtracer(itrc,ng) END DO END DO CASE ('TdecayM(isFsur)') Npts=load_r(Nval, Rval, Ngrids, Tdecay(isFsur,:)) CASE ('TdecayM(isUbar)') Npts=load_r(Nval, Rval, Ngrids, Tdecay(isUbar,:)) CASE ('TdecayM(isVbar)') Npts=load_r(Nval, Rval, Ngrids, Tdecay(isVbar,:)) CASE ('TdecayM(isUvel)') Npts=load_r(Nval, Rval, Ngrids, Tdecay(isUvel,:)) CASE ('TdecayM(isVvel)') Npts=load_r(Nval, Rval, Ngrids, Tdecay(isVvel,:)) CASE ('TdecayM(isTvar)') Npts=load_r(Nval, Rval, MT, Ngrids, Rtracer) DO ng=1,Ngrids DO itrc=1,NT(ng) Tdecay(isTvar(itrc),ng)=Rtracer(itrc,ng) END DO END DO CASE ('HdecayI(isFsur)') Npts=load_r(Nval, Rval, Ngrids, Hdecay(1,isFsur,:)) CASE ('HdecayI(isUbar)') Npts=load_r(Nval, Rval, Ngrids, Hdecay(1,isUbar,:)) CASE ('HdecayI(isVbar)') Npts=load_r(Nval, Rval, Ngrids, Hdecay(1,isVbar,:)) CASE ('HdecayI(isUvel)') Npts=load_r(Nval, Rval, Ngrids, Hdecay(1,isUvel,:)) CASE ('HdecayI(isVvel)') Npts=load_r(Nval, Rval, Ngrids, Hdecay(1,isVvel,:)) CASE ('HdecayI(isTvar)') Npts=load_r(Nval, Rval, MT, Ngrids, Rtracer) DO ng=1,Ngrids DO itrc=1,NT(ng) Hdecay(1,isTvar(itrc),ng)=Rtracer(itrc,ng) END DO END DO CASE ('VdecayI(isUvel)') Npts=load_r(Nval, Rval, Ngrids, Vdecay(1,isUvel,:)) CASE ('VdecayI(isVvel)') Npts=load_r(Nval, Rval, Ngrids, Vdecay(1,isVvel,:)) CASE ('VdecayI(isTvar)') Npts=load_r(Nval, Rval, MT, Ngrids, Rtracer) DO ng=1,Ngrids DO itrc=1,NT(ng) Vdecay(1,isTvar(itrc),ng)=Rtracer(itrc,ng) END DO END DO CASE ('STDnameI') label='STD - initial conditions standard deviation' Npts=load_s1d(Nval, Cval, Cdim, line, label, igrid, & & Ngrids, Nfiles, 4, 1, inp_lib, STD) CASE ('STDnameM') label='STD - model error standard deviation' Npts=load_s1d(Nval, Cval, Cdim, line, label, igrid, & & Ngrids, Nfiles, 4, 2, inp_lib, STD) CASE ('STDnameB') label='STD - boundary conditions standard deviation' Npts=load_s1d(Nval, Cval, Cdim, line, label, igrid, & & Ngrids, Nfiles, 4, 3, inp_lib, STD) CASE ('STDnameF') label='STD - surface forcing standard deviation' Npts=load_s1d(Nval, Cval, Cdim, line, label, igrid, & & Ngrids, Nfiles, 4, 4, inp_lib, STD) CASE ('NRMnameI') label='NRM - initial conditions normalization' Npts=load_s1d(Nval, Cval, Cdim, line, label, igrid, & & Ngrids, Nfiles, 4, 1, inp_lib, NRM) CASE ('NRMnameM') label='NRM - model error normalization' Npts=load_s1d(Nval, Cval, Cdim, line, label, igrid, & & Ngrids, Nfiles, 4, 2, inp_lib, NRM) CASE ('NRMnameB') label='NRM - boundary conditions normalization' Npts=load_s1d(Nval, Cval, Cdim, line, label, igrid, & & Ngrids, Nfiles, 4, 3, inp_lib, NRM) CASE ('NRMnameF') label='NRM - surface forcing normalization' Npts=load_s1d(Nval, Cval, Cdim, line, label, igrid, & & Ngrids, Nfiles, 4, 4, inp_lib, NRM) CASE ('OBSname') label='OBS - data assimilation observations' Npts=load_s1d(Nval, Cval, Cdim, line, label, igrid, & & Ngrids, Nfiles, inp_lib, OBS) CASE ('HSSname') label='HSS - Hessian eigenvectors' Npts=load_s1d(Nval, Cval, Cdim, line, label, igrid, & & Ngrids, Nfiles, out_lib, HSS) CASE ('LCZname') label='LCZ - Lanczos vectors' Npts=load_s1d(Nval, Cval, Cdim, line, label, igrid, & & Ngrids, Nfiles, out_lib, LCZ) CASE ('LZEname') label='LZE - Time-evolved Lanczos vectors' Npts=load_s1d(Nval, Cval, Cdim, line, label, igrid, & & Ngrids, Nfiles, out_lib, LZE) CASE ('MODname') label='DAV - 4D-Var data assimilation variables' Npts=load_s1d(Nval, Cval, Cdim, line, label, igrid, & & Ngrids, Nfiles, out_lib, DAV) CASE ('ERRname') label='ERR - 4D-Var posterior error covariance' Npts=load_s1d(Nval, Cval, Cdim, line, label, igrid, & & Ngrids, Nfiles, out_lib, ERR) END SELECT END IF END DO 10 IF (Master) WRITE (out,50) line exit_flag=4 RETURN 20 CONTINUE ! !----------------------------------------------------------------------- ! Report input parameters. !----------------------------------------------------------------------- ! IF (Master.and.Lwrite) THEN DO ng=1,Ngrids WRITE (out,60) ng WRITE (out,100) HevecErr, 'HevecErr', & & 'Accuracy required for eigenvectors.' WRITE (out,70) LhotStart, 'LhotStart', & & 'Switch for hot start of subsequent outer loops.' WRITE (out,70) Lprecond, 'Lprecond', & & 'Switch for conjugate gradient preconditioning.' WRITE (out,70) Lritz, 'Lritz', & & 'Switch for Ritz limited-memory preconditioning.' IF (Lprecond.and.(NritzEV.gt.0)) THEN WRITE (out,80) NritzEV, 'NritzEV', & & 'Number of preconditioning eigenpairs to use.' END IF WRITE (out,170) LdefNRM(1:4,ng), 'LdefNRM', & & 'Switch to create a normalization NetCDF file.' WRITE (out,170) LwrtNRM(1:4,ng), 'LwrtNRM', & & 'Switch to write out normalization factors.' IF (ANY(LwrtNRM(:,ng))) THEN IF (Nmethod(ng).eq.0) THEN WRITE (out,80) Nmethod(ng), 'Nmethod', & & 'Correlation normalization method: Exact.' ELSE IF (Nmethod(ng).eq.1) THEN WRITE (out,80) Nmethod(ng), 'Nmethod', & & 'Correlation normalization method: Randomization.' WRITE (out,80) Rscheme(ng), 'Rscheme', & & 'Random number generation scheme' WRITE (out,80) Nrandom, 'Nrandom', & & 'Number of iterations for randomization.' END IF END IF IF (ANY(LwrtNRM(:,ng))) THEN WRITE (out,70) Cnorm(2,isFsur), 'CnormM(isFsur)', & & 'Compute model 2D RHO-normalization factors.' WRITE (out,70) Cnorm(2,isUbar), 'CnormM(isUbar)', & & 'Compute model 2D U-normalization factors.' WRITE (out,70) Cnorm(2,isVbar), 'CnormM(isVbar)', & & 'Compute model 2D V-normalization factors.' WRITE (out,70) Cnorm(2,isUvel), 'CnormM(isUvel)', & & 'Compute model 3D U-normalization factors.' WRITE (out,70) Cnorm(2,isVvel), 'CnormM(isVvel)', & & 'Compute model 3D V-normalization factors.' DO itrc=1,NT(ng) WRITE (out,110) Cnorm(2,isTvar(itrc)), 'CnormM(isTvar)', & & 'Compute model normalization factors for tracer ', & & itrc, TRIM(Vname(1,idTvar(itrc))) END DO END IF IF (ANY(LwrtNRM(:,ng))) THEN WRITE (out,70) Cnorm(1,isFsur), 'CnormI(isFsur)', & & 'Compute initial 2D RHO-normalization factors.' WRITE (out,70) Cnorm(1,isUbar), 'CnormI(isUbar)', & & 'Compute initial 2D U-normalization factors.' WRITE (out,70) Cnorm(1,isVbar), 'CnormI(isVbar)', & & 'Compute initial 2D V-normalization factors.' WRITE (out,70) Cnorm(1,isUvel), 'CnormI(isUvel)', & & 'Compute initial 3D U-normalization factors.' WRITE (out,70) Cnorm(1,isVvel), 'CnormI(isVvel)', & & 'Compute initial 3D V-normalization factors.' DO itrc=1,NT(ng) WRITE (out,110) Cnorm(1,isTvar(itrc)), 'CnormI(isTvar)', & & 'Compute initial normalization factors for tracer ', & & itrc, TRIM(Vname(1,idTvar(itrc))) END DO END IF WRITE (out,100) Hgamma(1), 'Hgamma', & & 'Horizontal diffusion factor, initial conditions.' WRITE (out,100) Hgamma(2), 'HgammaM', & & 'Horizontal diffusion factor, model error.' WRITE (out,100) Vgamma(1), 'Vgamma', & & 'Vertical diffusion factor, initial conditions.' WRITE (out,100) Vgamma(2), 'VgammaM', & & 'Vertical diffusion factor, model error.' IF (nADJ(ng).lt.ntimes(ng)) THEN WRITE (out,120) Hdecay(2,isFsur,ng), 'HdecayM(isFsur)', & & 'Model decorrelation H-scale (m), free-surface.' WRITE (out,120) Hdecay(2,isUbar,ng), 'HdecayM(isUbar)', & & 'Model decorrelation H-scale (m), 2D U-momentum.' WRITE (out,120) Hdecay(2,isVbar,ng), 'HdecayM(isVbar)', & & 'Model decorrelation H-scale (m), 2D V-momentum.' WRITE (out,120) Hdecay(2,isUvel,ng), 'HdecayM(isUvel)', & & 'Model decorrelation H-scale (m), 3D U-momentum.' WRITE (out,120) Hdecay(2,isVvel,ng), 'HdecayM(isVvel)', & & 'Model decorrelation H-scale (m), 3D V-momentum.' DO itrc=1,NT(ng) WRITE (out,130) Hdecay(2,isTvar(itrc),ng), & & 'HdecayM(idTvar)', & & 'Model decorrelation H-scale (m), ', & & TRIM(Vname(1,idTvar(itrc))) END DO WRITE (out,120) Vdecay(2,isUvel,ng), 'VdecayM(isUvel)', & & 'Model decorrelation V-scale (m), 3D U-momentum.' WRITE (out,120) Vdecay(2,isVvel,ng), 'VdecayM(isVvel)', & & 'Model decorrelation V-scale (m), 3D V-momentum.' DO itrc=1,NT(ng) WRITE (out,130) Vdecay(2,isTvar(itrc),ng), & & 'VdecayM(idTvar)', & & 'Model decorrelation V-scale (m), ', & & TRIM(Vname(1,idTvar(itrc))) END DO END IF WRITE (out,120) Hdecay(1,isFsur,ng), 'HdecayI(isFsur)', & & 'Initial decorrelation H-scale (m), free-surface.' WRITE (out,120) Hdecay(1,isUbar,ng), 'HdecayI(isUbar)', & & 'Initial decorrelation H-scale (m), 2D U-momentum.' WRITE (out,120) Hdecay(1,isVbar,ng), 'HdecayI(isVbar)', & & 'Initial decorrelation H-scale (m), 2D V-momentum.' WRITE (out,120) Hdecay(1,isUvel,ng), 'HdecayI(isUvel)', & & 'Initial decorrelation H-scale (m), 3D U-momentum.' WRITE (out,120) Hdecay(1,isVvel,ng), 'HdecayI(isVvel)', & & 'Initial decorrelation H-scale (m), 3D V-momentum.' DO itrc=1,NT(ng) WRITE (out,130) Hdecay(1,isTvar(itrc),ng), & & 'HdecayI(idTvar)', & & 'Initial decorrelation H-scale (m), ', & & TRIM(Vname(1,idTvar(itrc))) END DO WRITE (out,120) Vdecay(1,isUvel,ng), 'VdecayI(isUvel)', & & 'Initial decorrelation V-scale (m), 3D U-momentum.' WRITE (out,120) Vdecay(1,isVvel,ng), 'VdecayI(isVvel)', & & 'Initial decorrelation V-scale (m), 3D V-momentum.' DO itrc=1,NT(ng) WRITE (out,130) Vdecay(1,isTvar(itrc),ng), & & 'VdecayI(idTvar)', & & 'Initial decorrelation V-scale (m), ', & & TRIM(Vname(1,idTvar(itrc))) END DO IF (NextraObs.gt.0) THEN WRITE (out,80) NextraObs, 'NextraObs', & & 'Number of extra-observations classes to process.' DO i=1,NextraObs WRITE (out,85) ExtraIndex(i), 'ExtraIndex', i, & & 'Extra-observation type index for: '// & & TRIM(ExtraName(i)) END DO END IF END DO END IF ! !----------------------------------------------------------------------- ! Report input files and check availability of input files. !----------------------------------------------------------------------- ! IF (Master.and.Lwrite) THEN WRITE (out,150) WRITE (out,160) ' Assimilation Parameters File: ', & & TRIM(aparnam) END IF ! DO ng=1,Ngrids fname=STD(2,ng)%name IF (NSA.eq.2) THEN IF (.not.find_file(ng, out, fname, 'STDnameM')) THEN IF (FoundError(exit_flag, NoError, 1291, MyFile)) RETURN ELSE IF (Master.and.Lwrite) WRITE (out,160) & & ' Model STD File: ', TRIM(fname) END IF END IF fname=STD(1,ng)%name IF (.not.find_file(ng, out, fname, 'STDnameI')) THEN IF (FoundError(exit_flag, NoError, 1300, MyFile)) RETURN ELSE IF (Master.and.Lwrite) WRITE (out,160) & & ' Initial Conditions STD File: ', TRIM(fname) END IF fname=NRM(2,ng)%name IF (Master.and.Lwrite) WRITE (out,160) & & ' Model Norm File: ', TRIM(fname) IF (.not.LdefNRM(2,ng)) THEN IF (.not.find_file(ng, out, fname, 'NRMnameM')) THEN IF (FoundError(exit_flag, NoError, 1332, MyFile)) RETURN END IF END IF fname=NRM(1,ng)%name IF (Master.and.Lwrite) WRITE (out,160) & & ' Initial Conditions Norm File: ', TRIM(fname) IF (.not.LdefNRM(1,ng)) THEN IF (.not.find_file(ng, out, fname, 'NRMnameI')) THEN IF (FoundError(exit_flag, NoError, 1360, MyFile)) RETURN END IF END IF fname=OBS(ng)%name IF (.not.find_file(ng, out, fname, 'OBSname')) THEN IF (FoundError(exit_flag, NoError, 1409, MyFile)) RETURN ELSE IF (Master.and.Lwrite) WRITE (out,160) & & ' Input observations File: ', TRIM(fname) END IF IF (Master.and.Lwrite) WRITE (out,160) & & ' Input/Output Lanczos File: ', TRIM(LCZ(ng)%name) IF (Master.and.Lwrite) WRITE (out,160) & & ' Input/Output Hessian File: ', TRIM(HSS(ng)%name) IF (Master.and.Lwrite) WRITE (out,160) & & ' Output 4D-Var File: ', TRIM(DAV(ng)%name) END DO ! ! Stop if activating pre-conditioning for the RBLanczos minimization ! algorithm. It does not work yet. ! DO ng=1,Ngrids IF (Lprecond) THEN IF (Master) THEN WRITE (out,230) 'Lprecond', Lprecond, & & 'pre-conditioning does not work yet with ' // & & uppercase('rpcg') // ' algorithm.', & & 'Set Lprecond to F in ' // TRIM(aparnam) END IF exit_flag=5 RETURN END IF END DO ! 50 FORMAT (/,' READ_AssPar - Error while processing line: ',/,a) 55 FORMAT (/,' READ_AssPar - ',a,i4,2x,i4,/,15x,a) 60 FORMAT (/,/,' Assimilation Parameters, Grid: ',i2.2, & & /, ' =================================',/) 70 FORMAT (10x,l1,2x,a,t30,a) 80 FORMAT (1x,i10,2x,a,t30,a) 85 FORMAT (1x,i10,2x,a,'(',i2.2,')',t30,a) 90 FORMAT (1x,i10,2x,a,t30,a,/,t32,a) 95 FORMAT (1x,a,2x,a,t30,a) 100 FORMAT (1p,e11.4,2x,a,t30,a) 110 FORMAT (10x,l1,2x,a,t30,a,i2.2,':',1x,a) 120 FORMAT (f11.3,2x,a,t30,a) 130 FORMAT (f11.3,2x,a,t30,a,a,'.') 150 FORMAT (/,' Input/Output Assimilation Files:',/) 160 FORMAT (2x,a,a) 170 FORMAT (3x,4(1x,l1),2x,a,t30,a) 180 FORMAT (3x,4(1x,l1),2x,a,t30,a,i2.2,':',1x,a) 190 FORMAT (1p,e11.4,2x,a,t30,a,a,'.') 200 FORMAT (1p,e11.4,2x,a,'(',i2.2,')',t30,a,':',1x,i6) 210 FORMAT (/,' READ_ASSPAR - variable info not yet loaded, ', a) 220 FORMAT (/,' READ_ASSPAR - Grid ', i2.2, & & ', could not find input file: ',a) 230 FORMAT (/,' READ_ASSPAR - Illegal parameter, ', a, ' = ', 1x, l1, & & /,15x,a,/,15x,a) 240 FORMAT (/,' READ_ASSPAR - Illegal parameter, ', a, ' = ', 1x, i2, & & /,15x,a) RETURN END SUBROUTINE read_AssPar