- Timestamp:
- 2015-02-05T17:29:55+01:00 (9 years ago)
- Location:
- branches/2014/dev_r4650_UKMO14.12_STAND_ALONE_OBSOPER/NEMOGCM/NEMO/SAO_SRC
- Files:
-
- 2 deleted
- 5 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2014/dev_r4650_UKMO14.12_STAND_ALONE_OBSOPER/NEMOGCM/NEMO/SAO_SRC/nemogcm.F90
r5042 r5063 97 97 !! 98 98 !!---------------------------------------------------------------------- 99 !! Class 4 output stream switch100 USE obs_fbm, ONLY: ln_cl4101 99 !! Initialise NEMO 102 100 CALL nemo_init 103 101 !! Initialise Stand Alone Observation operator data 104 CALL sao_data_init( ln_cl4 ) 105 !! Loop over various model counterparts 106 DO jimatch = 1, cl4_match_len 107 !! Initialise obs_oper 108 CALL dia_obs_init 109 !! Interpolate to observation space 110 CALL sao_interp 111 !! Pipe to output files 112 CALL dia_obs_wri 113 !! Reset the obs_oper between 114 CALL dia_obs_dealloc 115 END DO 102 CALL sao_data_init 103 !! Initialise obs_oper 104 CALL dia_obs_init 105 !! Interpolate to observation space 106 CALL sao_interp 107 !! Pipe to output files 108 CALL dia_obs_wri 109 !! Reset the obs_oper between 110 CALL dia_obs_dealloc 116 111 !! Safely stop MPI 117 112 IF(lk_mpp) CALL mppstop ! end mpp communications -
branches/2014/dev_r4650_UKMO14.12_STAND_ALONE_OBSOPER/NEMOGCM/NEMO/SAO_SRC/sao.nml
r4846 r5063 1 2 1 !---------------------------------------------------------------------- 3 2 ! namsao Stand Alone Observation operator namelist … … 10 9 nn_sao_idx = 1 11 10 / 12 !----------------------------------------------------------------------13 ! namcl4 Class 4 obs_oper namelist14 !----------------------------------------------------------------------15 ! cl4_date the verfication date of the class 4 file16 ! cl4_vars the name of the variable in the class 4 file17 ! cl4_sys the forecast system being used e.g. FOAM18 ! cl4_cfg the model configuration being used e.g. amm719 ! cl4_vn the version number e.g. 12.020 ! cl4_prefix prefix which denotes the output file21 ! cl4_contact email address of file creator22 ! cl4_inst institution related to the data within the file23 ! cl4_leadtime lead time axis of class 4 file24 ! cl4_fcst_idx output file forecast index25 ! cl4_fcst_len output file forecast dimension length26 ! cl4_match_len number of match types27 !28 &namcl429 cl4_leadtime = 1230 cl4_fcst_idx = 131 cl4_fcst_len = 132 cl4_match_len = 133 cl4_date = '20130101'34 cl4_vars = 'forecast'35 cl4_sys = 'FOAM'36 cl4_cfg = 'amm7'37 cl4_vn = '12.0'38 cl4_prefix = 'class4'39 cl4_contact = ''40 cl4_inst = 'UK Met Office'41 /42 43 44 -
branches/2014/dev_r4650_UKMO14.12_STAND_ALONE_OBSOPER/NEMOGCM/NEMO/SAO_SRC/sao_data.F90
r4849 r5063 12 12 INTEGER, PARAMETER :: MaxNumFiles = 1000 13 13 14 !! Class 4 file settings15 INTEGER :: &16 & cl4_fcst_idx(MaxNumFiles), & !: forecast indices17 & cl4_match_len, & !: number of match types18 & cl4_fcst_len !: number of forecast days19 CHARACTER(len=lc) :: &20 & cl4_vars(MaxNumFiles), & !: class 4 variables21 & cl4_sys, & !: class 4 system22 & cl4_cfg, & !: class 4 configuration23 & cl4_date, & !: class 4 date24 & cl4_vn, & !: class 4 version25 & cl4_prefix, & !: class 4 prefix26 & cl4_contact, & !: class 4 contact27 & cl4_inst !: class 4 institute28 REAL :: cl4_modjuld !: model Julian day29 REAL :: &30 & cl4_leadtime(MaxNumFiles) !: Lead time data31 32 14 !! Stand Alone Observation operator settings 33 15 CHARACTER(len=lc) :: & 34 16 & sao_files(MaxNumFiles) !: model files 35 17 INTEGER :: & 36 & jifile, & !: current file list index37 18 & n_files, & !: number of files 38 & jimatch, & !: current match39 19 & nn_sao_idx(MaxNumFiles), & !: time_counter indices 40 20 & nn_sao_freq !: read frequency in time steps 41 CHARACTER(len=128) :: &42 & alt_file !: altimeter file43 21 CONTAINS 44 SUBROUTINE sao_data_init( ld_cl4)22 SUBROUTINE sao_data_init() 45 23 !!---------------------------------------------------------------------- 46 24 !! *** SUBROUTINE sao_data_init *** … … 53 31 & jf !: file dummy loop index 54 32 LOGICAL :: lmask(MaxNumFiles) !: Logical mask used for counting 55 LOGICAL, INTENT(IN) :: ld_cl4 !: Logical class 4 on/off56 33 INTEGER :: ios 57 34 … … 59 36 NAMELIST/namsao/sao_files, nn_sao_idx, nn_sao_freq 60 37 61 ! Class 4 file specifiers62 NAMELIST/namcl4/cl4_vars, cl4_sys, cl4_cfg, cl4_date, cl4_vn, &63 & cl4_prefix, cl4_contact, cl4_inst, cl4_leadtime, &64 & cl4_fcst_idx, cl4_fcst_len, cl4_match_len65 66 38 ! Standard offline obs_oper initialisation 67 jimatch = 0 !: match-up iteration variable68 jifile = 1 !: input file iteration variable69 39 n_files = 0 !: number of files to cycle through 70 40 sao_files(:) = '' !: list of files to read in 71 41 nn_sao_idx(:) = 0 !: list of indices inside each file 72 42 nn_sao_freq = -1 !: input frequency in time steps 73 74 ! Class 4 initialisation75 cl4_leadtime(:) = 0 !: Lead time axis value for each file76 cl4_fcst_len = 0 !: Length of the forecast dimension77 cl4_match_len = 1 !: Number of match types78 cl4_fcst_idx(:) = 0 !: output file forecast index79 cl4_vars(:) = '' !: output file variable names80 cl4_sys = '' !: output file system81 cl4_cfg = '' !: output file configuration82 cl4_date = '' !: output file date string83 cl4_vn = '' !: output file version84 cl4_prefix = 'class4' !: output file prefix85 cl4_contact = '' !: output file contact details86 cl4_inst = '' !: output file institution87 43 88 44 ! Standard offline obs_oper settings … … 95 51 902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsao in configuration namelist', .TRUE. ) 96 52 97 ! Read class 4 output settings98 IF (ld_cl4) THEN99 REWIND( numnam_ref ) ! Namelist namctl in reference namelist : Control prints & Benchmark100 READ ( numnam_ref, namcl4, IOSTAT = ios, ERR = 903 )101 903 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namcl4 in reference namelist', .TRUE. )102 103 REWIND( numnam_cfg ) ! Namelist namctl in confguration namelist : Control prints & Benchmark104 READ ( numnam_cfg, namcl4, IOSTAT = ios, ERR = 904 )105 904 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namcl4 in configuration namelist', .TRUE. )106 ENDIF107 53 108 54 ! count input files … … 122 68 WRITE(numout,*) 'offline obs_oper : Initialization' 123 69 WRITE(numout,*) '~~~~~~~~~~~~~~~~~' 124 WRITE(numout,*) ' Namelist namsao : set stand alone obs_oper parameters' 70 WRITE(numout,*) ' Namelist namsao : set stand alone obs_oper parameters' 125 71 DO jf = 1, n_files 126 72 WRITE(numout,'(1X,2A)') ' Input forecast file name forecastfile = ', & … … 128 74 WRITE(numout,*) ' Input forecast file index forecastindex = ', & 129 75 nn_sao_idx(jf) 130 WRITE(numout,*) ' Output forecast leadtime index leadtimeindex = ', &131 cl4_fcst_idx(jf)132 WRITE(numout,*) ' Output forecast leadtime value leadtimevalue = ', &133 cl4_leadtime(jf)134 WRITE(numout,'(1X,2A)') ' Input class 4 variable class 4 parameter = ', &135 TRIM(cl4_vars(jf))136 76 END DO 137 WRITE(numout, '(1X,2A)') ' Input class 4 system class 4 system = ', &138 TRIM(cl4_sys)139 WRITE(numout, '(1X,2A)') ' Input class 4 config class 4 config = ', &140 TRIM(cl4_cfg)141 WRITE(numout, '(1X,2A)') ' Input class 4 date class 4 date = ', &142 TRIM(cl4_date)143 WRITE(numout, '(1X,2A)') ' Input class 4 version class 4 version = ', &144 TRIM(cl4_vn)145 WRITE(numout, '(1X,2A)') ' Input class 4 prefix class 4 prefix = ', &146 TRIM(cl4_prefix)147 WRITE(numout, '(1X,2A)') ' Input class 4 contact class 4 contact = ', &148 TRIM(cl4_contact)149 WRITE(numout, '(1X,2A)') ' Input class 4 institute class 4 institute = ', &150 TRIM(cl4_inst)151 77 END IF 152 78 -
branches/2014/dev_r4650_UKMO14.12_STAND_ALONE_OBSOPER/NEMOGCM/NEMO/SAO_SRC/sao_intp.F90
r4852 r5063 27 27 !! 2. Call dia_obs at appropriate time steps 28 28 !!---------------------------------------------------------------------- 29 INTEGER :: istp ! time step index 30 !! Loop over entire run 29 INTEGER :: & 30 & istp, & ! time step index 31 & ifile ! file index 31 32 istp = nit000 - 1 32 33 nstop = 0 34 ifile = 1 35 CALL sao_rea_dri(ifile) 33 36 DO WHILE ( istp <= nitend .AND. nstop == 0 ) 34 IF ( jifile <= n_files + 1) THEN37 IF (ifile <= n_files + 1) THEN 35 38 IF ( MOD(istp, nn_sao_freq) == nit000 ) THEN 36 !! Read next model counterpart 37 CALL sao_rea_dri(jifile) 38 jifile = jifile + 1 39 CALL sao_rea_dri(ifile) 40 ifile = ifile + 1 39 41 ENDIF 40 !! Interpolate single time step41 42 CALL dia_obs(istp) 42 43 ENDIF 43 !! Increment model step44 44 istp = istp + 1 45 45 END DO -
branches/2014/dev_r4650_UKMO14.12_STAND_ALONE_OBSOPER/NEMOGCM/NEMO/SAO_SRC/sao_read.F90
r4849 r5063 1 2 1 MODULE sao_read 3 2 !!================================================================== … … 5 4 !! Read routines : I/O for Stand Alone Observation operator 6 5 !!================================================================== 7 8 6 USE mppini 9 7 USE lib_mpp … … 28 26 !! 29 27 !! Purpose : To choose appropriate read method 30 !! Method : 28 !! Method : 31 29 !! 32 30 !! Author : A. Ryan Oct 2013 … … 36 34 & kfile !: File number 37 35 CHARACTER(len=lc) :: & 38 & cdfilename, & !: File name 39 & cmatchname !: Match name 36 & cdfilename !: File name 40 37 INTEGER :: & 41 & kindex !: File index to read 42 43 !! Filename, index and match-up kind 38 & kindex !: File index to read 39 44 40 cdfilename = TRIM(sao_files(kfile)) 45 cmatchname = TRIM(cl4_vars(kfile))46 41 kindex = nn_sao_idx(kfile) 47 48 !! Update model fields 49 !! Class 4 variables: forecast, persistence, 50 !! nrt_analysis, best_estimate 51 !! Feedback variables: empty string 52 IF ( (TRIM(cmatchname) == 'forecast') .OR. & 53 & (TRIM(cmatchname) == 'persistence') .OR. & 54 & (TRIM(cmatchname) == 'nrt_analysis') .OR. & 55 & (TRIM(cmatchname) == 'best_estimate').OR. & 56 & (TRIM(cmatchname) == '') ) THEN 57 CALL sao_read_file(TRIM(cdfilename), kindex) 58 CALL sao_read_juld(TRIM(cdfilename), kindex, cl4_modjuld) 59 ELSE IF (TRIM(cmatchname) == 'climatology') THEN 60 WRITE(numout,*) 'Interpolating climatologies' 61 ELSE IF (TRIM(cmatchname) == 'altimeter') THEN 62 CALL sao_read_altbias(TRIM(cdfilename)) 63 CALL sao_read_juld(TRIM(cdfilename), kindex, cl4_modjuld) 64 END IF 42 CALL sao_read_file(TRIM(cdfilename), kindex) 65 43 66 44 END SUBROUTINE sao_rea_dri 67 68 SUBROUTINE sao_read_altbias(filename)69 !!------------------------------------------------------------------------70 !! *** sao_read_altbias ***71 !!72 !! Purpose : To read altimeter bias and set tn,sn to missing values73 !! Method : Use subdomain indices to create start and count matrices74 !! for netcdf read.75 !!76 !! Author : A. Ryan Sept 201277 !!------------------------------------------------------------------------78 CHARACTER(len=*), INTENT(IN) :: filename79 INTEGER :: ncid, &80 & varid,&81 & istat,&82 & ntimes,&83 & tdim, &84 & xdim, &85 & ydim, &86 & zdim87 INTEGER :: ii, ij, ik88 INTEGER, DIMENSION(3) :: start_s, &89 & count_s90 REAL(fbdp), DIMENSION(:,:), ALLOCATABLE :: temp_sshn91 REAL(fbdp) :: fill_val92 93 IF (TRIM(filename) == 'nofile') THEN94 tsn(:,:,:,:) = fbrmdi95 sshn(:,:) = fbrmdi96 ELSE97 ! Open Netcdf file to find dimension id98 istat = nf90_open(trim(filename),nf90_nowrite,ncid)99 istat = nf90_inq_dimid(ncid,'x',xdim)100 istat = nf90_inq_dimid(ncid,'y',ydim)101 istat = nf90_inq_dimid(ncid,'deptht',zdim)102 istat = nf90_inq_dimid(ncid,'time',tdim)103 istat = nf90_inquire_dimension(ncid, tdim, len=ntimes)104 ! Allocate temporary temperature array105 ALLOCATE(temp_sshn(nlci,nlcj))106 ! Create start and count arrays107 start_s = (/ nimpp, njmpp, 1 /)108 count_s = (/ nlci, nlcj, 1 /)109 110 ! Altimeter bias111 istat = nf90_inq_varid(ncid,'altbias',varid)112 istat = nf90_get_att(ncid, varid, '_FillValue', fill_val)113 istat = nf90_get_var(ncid, varid, temp_sshn, start_s, count_s)114 WHERE(temp_sshn(:,:) == fill_val) temp_sshn(:,:) = fbrmdi115 116 ! Initialise tsn, sshn to fbrmdi117 tsn(:,:,:,:) = fbrmdi118 sshn(:,:) = fbrmdi119 120 ! Fill sshn with altimeter bias121 sshn(1:nlci,1:nlcj) = temp_sshn(:,:) * tmask(1:nlci,1:nlcj,1)122 123 ! Remove halo from tmask, sshn to prevent double obs counting124 IF (jpi > nlci) THEN125 tmask(nlci+1:,:,:) = 0126 sshn(nlci+1:,:) = 0127 END IF128 IF (jpj > nlcj) THEN129 tmask(:,nlcj+1:,:) = 0130 sshn(:,nlcj+1:) = 0131 END IF132 133 ! Deallocate arrays134 DEALLOCATE(temp_sshn)135 136 ! Close netcdf file137 istat = nf90_close(ncid)138 END IF139 140 END SUBROUTINE sao_read_altbias141 45 142 46 SUBROUTINE sao_read_file(filename, ifcst) … … 176 80 IF (TRIM(filename) == 'nofile') THEN 177 81 tsn(:,:,:,:) = fbrmdi 178 sshn(:,:) = fbrmdi 82 sshn(:,:) = fbrmdi 179 83 ELSE 180 84 WRITE(numout,*) "Opening :", TRIM(filename) … … 190 94 istat = nf90_inq_dimid(ncid,'time_counter',tdim) 191 95 istat = nf90_inquire_dimension(ncid, tdim, len=ntimes) 192 IF (ifcst .LE. ntimes) THEN 96 IF (ifcst .LE. ntimes) THEN 193 97 ! Allocate temporary temperature array 194 98 ALLOCATE(temp_tn(nlci,nlcj,jpk)) 195 99 ALLOCATE(temp_sn(nlci,nlcj,jpk)) 196 100 ALLOCATE(temp_sshn(nlci,nlcj)) 197 101 198 102 ! Set temp_tn, temp_sn to 0. 199 103 temp_tn(:,:,:) = fbrmdi 200 104 temp_sn(:,:,:) = fbrmdi 201 105 temp_sshn(:,:) = fbrmdi 202 106 203 107 ! Create start and count arrays 204 108 start_n = (/ nimpp, njmpp, 1, ifcst /) … … 206 110 start_s = (/ nimpp, njmpp, ifcst /) 207 111 count_s = (/ nlci, nlcj, 1 /) 208 112 209 113 ! Read information into temporary arrays 210 114 ! retrieve varid and read in temperature … … 213 117 istat = nf90_get_var(ncid, varid, temp_tn, start_n, count_n) 214 118 WHERE(temp_tn(:,:,:) == fill_val) temp_tn(:,:,:) = fbrmdi 215 119 216 120 ! retrieve varid and read in salinity 217 121 istat = nf90_inq_varid(ncid,'vosaline',varid) … … 219 123 istat = nf90_get_var(ncid, varid, temp_sn, start_n, count_n) 220 124 WHERE(temp_sn(:,:,:) == fill_val) temp_sn(:,:,:) = fbrmdi 221 125 222 126 ! retrieve varid and read in SSH 223 127 istat = nf90_inq_varid(ncid,'sossheig',varid) … … 226 130 istat = nf90_inq_varid(ncid,'altbias',varid) 227 131 END IF 228 132 229 133 istat = nf90_get_att(ncid, varid, '_FillValue', fill_val) 230 134 istat = nf90_get_var(ncid, varid, temp_sshn, start_s, count_s) 231 135 WHERE(temp_sshn(:,:) == fill_val) temp_sshn(:,:) = fbrmdi 232 136 233 137 ! Initialise tsn, sshn to fbrmdi 234 138 tsn(:,:,:,:) = fbrmdi 235 sshn(:,:) = fbrmdi 139 sshn(:,:) = fbrmdi 236 140 237 141 ! Mask out missing data index … … 256 160 ! Deallocate arrays 257 161 DEALLOCATE(temp_tn, temp_sn, temp_sshn) 258 ELSE 162 ELSE 259 163 ! Mark all as missing data 260 164 tsn(:,:,:,:) = fbrmdi 261 sshn(:,:) = fbrmdi 165 sshn(:,:) = fbrmdi 262 166 ENDIF 263 167 ! Close netcdf file … … 266 170 END IF 267 171 END SUBROUTINE sao_read_file 268 269 SUBROUTINE sao_read_juld(filename, ifcst, julian) 270 USE calendar 271 !!-------------------------------------------------------------------- 272 !! *** sao_read_juld *** 273 !! 274 !! Purpose : To read model julian day information from file 275 !! Author : A. Ryan Nov 2010 276 !!-------------------------------------------------------------------- 277 278 !! Routine arguments 279 CHARACTER(len=*), INTENT(IN) :: filename 280 INTEGER, INTENT(IN) :: ifcst 281 REAL, INTENT(OUT) :: julian !: Julian day 282 283 !! Local variables 284 INTEGER :: year, & !: Date information 285 & month, & 286 & day, & 287 & hour, & 288 & minute,& 289 & second 290 INTEGER :: istat, & !: Netcdf variables 291 & ncid, & 292 & dimid, & 293 & varid, & 294 & ntime 295 REAL,DIMENSION(:),ALLOCATABLE :: r_sec !: Remainder seconds 296 CHARACTER(len=120) :: time_str !: time string 297 298 time_str='' 299 !! Read in time_counter and remainder seconds 300 istat = nf90_open(trim(filename),nf90_nowrite,ncid) 301 istat = nf90_inq_dimid(ncid,'time_counter',dimid) 302 IF (istat /= nf90_noerr) THEN 303 istat = nf90_inq_dimid(ncid,'time',dimid) 304 ENDIF 305 istat = nf90_inquire_dimension(ncid,dimid,len=ntime) 306 istat = nf90_inq_varid(ncid,'time_counter',varid) 307 IF (istat /= nf90_noerr) THEN 308 istat = nf90_inq_dimid(ncid,'time',dimid) 309 ENDIF 310 istat = nf90_get_att(ncid,varid,'units',time_str) 311 ALLOCATE(r_sec(ntime)) 312 istat = nf90_get_var(ncid,varid, r_sec) 313 istat = nf90_close(ncid) 314 315 !! Fill yyyy-mm-dd hh:mm:ss 316 !! format(('seconds since ', I4.4,'-',I2.2,'-',I2.2,' ',I2.2,':',I2.2,':',I2.2)) 317 100 format((14x, I4.4,1x,I2.2,1x,I2.2,1x,I2.2,1x,I2.2,1x,I2.2)) 318 READ( time_str, 100 ) year, month, day, hour, minute, second 319 320 CALL ymds2ju(year, month, day, r_sec(ifcst), julian) 321 322 !! To take a comment from the ymds2ju subroutine 323 324 !- In the case of the Gregorian calendar we have chosen to use 325 !- the Lilian day numbers. This is the day counter which starts 326 !- on the 15th October 1582. 327 !- This is the day at which Pope Gregory XIII introduced the 328 !- Gregorian calendar. 329 !- Compared to the true Julian calendar, which starts some 330 !- 7980 years ago, the Lilian days are smaler and are dealt with 331 !- easily on 32 bit machines. With the true Julian days you can only 332 !- the fraction of the day in the real part to a precision of 333 !- a 1/4 of a day with 32 bits. 334 335 !! The obs operator routines prefer to calculate Julian days since 336 !! 01/01/1950 00:00:00 337 !! In order to convert to the 1950 version we must adjust by the number 338 !! of days between 15th October 1582 and 1st Jan 1950 339 340 julian = julian - 134123. 341 342 DEALLOCATE(r_sec) 343 344 END SUBROUTINE sao_read_juld 345 346 END MODULE sao_read 347 172 END MODULE sao_read
Note: See TracChangeset
for help on using the changeset viewer.