- Timestamp:
- 2020-12-12T12:31:26+01:00 (3 years ago)
- Location:
- branches/UKMO/dev_r5518_obs_oper_update_sit_SMOS/NEMOGCM
- Files:
-
- 8 edited
-
CONFIG/SHARED/namelist_ref (modified) (5 diffs)
-
NEMO/OPA_SRC/OBS/diaobs.F90 (modified) (47 diffs)
-
NEMO/OPA_SRC/OBS/obs_oper.F90 (modified) (21 diffs)
-
NEMO/OPA_SRC/OBS/obs_prep.F90 (modified) (10 diffs)
-
NEMO/OPA_SRC/OBS/obs_read_surf.F90 (modified) (14 diffs)
-
NEMO/OPA_SRC/OBS/obs_rot_vel.F90 (modified) (5 diffs)
-
NEMO/OPA_SRC/OBS/obs_surf_def.F90 (modified) (8 diffs)
-
NEMO/OPA_SRC/OBS/obs_write.F90 (modified) (11 diffs)
Legend:
- Unmodified
- Added
- Removed
-
branches/UKMO/dev_r5518_obs_oper_update_sit_SMOS/NEMOGCM/CONFIG/SHARED/namelist_ref
r12820 r14165 1202 1202 ln_fbd = .false. ! Logical switch for Sea Ice Freeboard observations 1203 1203 ln_vel3d = .false. ! Logical switch for velocity observations 1204 ln_sss = .false. ! Logical swithc for SSS observations 1204 ln_sss = .false. ! Logical switch for SSS observations 1205 ln_ssv = .false. ! Logical switch for SSV (surface velocity) observations 1205 1206 ln_slchltot = .false. ! Logical switch for surface total log10(chlorophyll) obs 1206 1207 ln_slchldia = .false. ! Logical switch for surface diatom log10(chlorophyll) obs … … 1241 1242 ln_sst_fp_indegs = .true. 1242 1243 ln_sss_fp_indegs = .true. 1244 ln_ssv_fp_indegs = .true. 1243 1245 ln_sic_fp_indegs = .true. 1244 1246 ln_sit_fp_indegs = .true. … … 1253 1255 cn_velfbfiles = 'vel_01.nc' ! Velocity feedback input observation file names 1254 1256 cn_sssfbfiles = 'sss_01.nc' ! SSS feedback input observation file names 1257 cn_ssvfbfiles = 'ssv_01.nc' ! SSV feedback input observation file names 1255 1258 cn_slchltotfbfiles = 'slchltot_01.nc' ! Surface total log10(chlorophyll) obs file names 1256 1259 cn_slchldiafbfiles = 'slchldia_01.nc' ! Surface diatom log10(chlorophyll) obs file names … … 1289 1292 rn_sss_avglamscl = 0. ! E/W diameter of SSS observation footprint (metres/degrees) 1290 1293 rn_sss_avgphiscl = 0. ! N/S diameter of SSS observation footprint (metres/degrees) 1294 rn_ssv_avglamscl = 0. ! E/W diameter of SSV observation footprint (metres/degrees) 1295 rn_ssv_avgphiscl = 0. ! N/S diameter of SSV observation footprint (metres/degrees) 1291 1296 rn_sic_avglamscl = 0. ! E/W diameter of SIC observation footprint (metres/degrees) 1292 1297 rn_sic_avgphiscl = 0. ! N/S diameter of SIC observation footprint (metres/degrees) … … 1300 1305 nn_2dint_sst = -1 ! Horizontal interpolation method for SST 1301 1306 nn_2dint_sss = -1 ! Horizontal interpolation method for SSS 1307 nn_2dint_ssv = -1 ! Horizontal interpolation method for SSV 1302 1308 nn_2dint_sic = -1 ! Horizontal interpolation method for SIC 1303 1309 nn_2dint_sit = -1 ! Horizontal interpolation method for SIT -
branches/UKMO/dev_r5518_obs_oper_update_sit_SMOS/NEMOGCM/NEMO/OPA_SRC/OBS/diaobs.F90
r12820 r14165 51 51 LOGICAL :: ln_sstnight !: Logical switch for night mean SST obs 52 52 LOGICAL :: ln_default_fp_indegs !: T=> Default obs footprint size specified in degrees, F=> in metres 53 LOGICAL :: ln_sla_fp_indegs !: T=> SLA obs footprint size specified in degrees, F=> in metres 54 LOGICAL :: ln_sst_fp_indegs !: T=> SST obs footprint size specified in degrees, F=> in metres 55 LOGICAL :: ln_sss_fp_indegs !: T=> SSS obs footprint size specified in degrees, F=> in metres 53 LOGICAL :: ln_sla_fp_indegs !: T=> SLA obs footprint size specified in degrees, F=> in metres 54 LOGICAL :: ln_sst_fp_indegs !: T=> SST obs footprint size specified in degrees, F=> in metres 55 LOGICAL :: ln_sss_fp_indegs !: T=> SSS obs footprint size specified in degrees, F=> in metres 56 LOGICAL :: ln_ssv_fp_indegs !: T=> SSV obs footprint size specified in degrees, F=> in metres 56 57 LOGICAL :: ln_sic_fp_indegs !: T=> SIC obs footprint size specified in degrees, F=> in metres 57 58 LOGICAL :: ln_sit_fp_indegs !: T=> SIT obs footprint size specified in degrees, F=> in metres … … 68 69 REAL(wp) :: rn_sss_avglamscl !: E/W diameter of SSS observation footprint 69 70 REAL(wp) :: rn_sss_avgphiscl !: N/S diameter of SSS observation footprint 71 REAL(wp) :: rn_ssv_avglamscl !: E/W diameter of SSV observation footprint 72 REAL(wp) :: rn_ssv_avgphiscl !: N/S diameter of SSV observation footprint 70 73 REAL(wp) :: rn_sic_avglamscl !: E/W diameter of SIC observation footprint 71 74 REAL(wp) :: rn_sic_avgphiscl !: N/S diameter of SIC observation footprint … … 83 86 INTEGER :: nn_2dint_sst !: SST horizontal interpolation method (-1 = default) 84 87 INTEGER :: nn_2dint_sss !: SSS horizontal interpolation method (-1 = default) 88 INTEGER :: nn_2dint_ssv !: SSV horizontal interpolation method (-1 = default) 85 89 INTEGER :: nn_2dint_sic !: SIC horizontal interpolation method (-1 = default) 86 90 INTEGER :: nn_2dint_sit !: SIT horizontal interpolation method (-1 = default) … … 174 178 & cn_velfbfiles, & ! Velocity profile input filenames 175 179 & cn_sssfbfiles, & ! Sea surface salinity input filenames 180 & cn_ssvfbfiles, & ! Sea surface velocity input filenames 176 181 & cn_slchltotfbfiles, & ! Surface total log10(chlorophyll) input filenames 177 182 & cn_slchldiafbfiles, & ! Surface diatom log10(chlorophyll) input filenames … … 213 218 LOGICAL :: ln_fbd ! Logical switch for sea ice freeboard 214 219 LOGICAL :: ln_sss ! Logical switch for sea surface salinity obs 220 LOGICAL :: ln_ssv ! Logical switch for sea surface velocity obs 215 221 LOGICAL :: ln_vel3d ! Logical switch for velocity (u,v) obs 216 222 LOGICAL :: ln_slchltot ! Logical switch for surface total log10(chlorophyll) obs … … 267 273 REAL(wp), POINTER, DIMENSION(:,:,:,:) :: & 268 274 & zmask ! Model land/sea mask associated with variables 275 REAL(wp), POINTER, DIMENSION(:,:,:) :: & 276 & zmask_surf ! Surface model land/sea mask associated with variables 269 277 270 278 271 279 NAMELIST/namobs/ln_diaobs, ln_t3d, ln_s3d, ln_sla, & 272 280 & ln_sst, ln_sic, ln_sit, ln_fbd, & 273 & ln_sss, ln_ vel3d,&281 & ln_sss, ln_ssv, ln_vel3d, & 274 282 & ln_slchltot, ln_slchldia, ln_slchlnon, & 275 283 & ln_slchldin, ln_slchlmic, ln_slchlnan, & … … 287 295 & ln_time_mean_sla_bkg, ln_default_fp_indegs, & 288 296 & ln_sla_fp_indegs, ln_sst_fp_indegs, & 289 & ln_sss_fp_indegs, ln_sic_fp_indegs, & 290 & ln_sit_fp_indegs, ln_fbd_fp_indegs, & 297 & ln_sss_fp_indegs, ln_ssv_fp_indegs, & 298 & ln_sic_fp_indegs, ln_sit_fp_indegs, & 299 & ln_fbd_fp_indegs, & 291 300 & cn_profbfiles, cn_slafbfiles, & 292 301 & cn_sstfbfiles, cn_sicfbfiles, & 293 302 & cn_sitfbfiles, cn_fbdfbfiles, & 294 & cn_velfbfiles, cn_sssfbfiles, 303 & cn_velfbfiles, cn_sssfbfiles, cn_ssvfbfiles, & 295 304 & cn_slchltotfbfiles, cn_slchldiafbfiles, & 296 305 & cn_slchlnonfbfiles, cn_slchldinfbfiles, & … … 312 321 & rn_sst_avglamscl, rn_sst_avgphiscl, & 313 322 & rn_sss_avglamscl, rn_sss_avgphiscl, & 323 & rn_ssv_avglamscl, rn_ssv_avgphiscl, & 314 324 & rn_sic_avglamscl, rn_sic_avgphiscl, & 315 325 & rn_sit_avglamscl, rn_sit_avgphiscl, & … … 317 327 & nn_1dint, nn_2dint_default, & 318 328 & nn_2dint_sla, nn_2dint_sst, & 319 & nn_2dint_sss, nn_2dint_sic, & 320 & nn_2dint_sit, nn_2dint_fbd, & 329 & nn_2dint_sss, nn_2dint_ssv, & 330 & nn_2dint_sic, nn_2dint_sit, & 331 & nn_2dint_fbd, & 321 332 & nn_msshc, rn_mdtcorr, rn_mdtcutoff, & 322 333 & nn_profdavtypes … … 335 346 cn_velfbfiles(:) = '' 336 347 cn_sssfbfiles(:) = '' 348 cn_ssvfbfiles(:) = '' 337 349 cn_slchltotfbfiles(:) = '' 338 350 cn_slchldiafbfiles(:) = '' … … 400 412 WRITE(numout,*) ' Logical switch for velocity observations ln_vel3d = ', ln_vel3d 401 413 WRITE(numout,*) ' Logical switch for SSS observations ln_sss = ', ln_sss 414 WRITE(numout,*) ' Logical switch for SSV observations ln_ssv = ', ln_ssv 402 415 WRITE(numout,*) ' Logical switch for surface total logchl obs ln_slchltot = ', ln_slchltot 403 416 WRITE(numout,*) ' Logical switch for surface diatom logchl obs ln_slchldia = ', ln_slchldia … … 435 448 WRITE(numout,*) ' Type of horizontal interpolation method for SST nn_2dint_sst = ', nn_2dint_sst 436 449 WRITE(numout,*) ' Type of horizontal interpolation method for SSS nn_2dint_sss = ', nn_2dint_sss 450 WRITE(numout,*) ' Type of horizontal interpolation method for SSV nn_2dint_ssv = ', nn_2dint_ssv 437 451 WRITE(numout,*) ' Type of horizontal interpolation method for SIC nn_2dint_sic = ', nn_2dint_sic 438 452 WRITE(numout,*) ' Type of horizontal interpolation method for SIT nn_2dint_sit = ', nn_2dint_sit … … 477 491 & ln_pchltot, ln_pno3, ln_psi4, ln_ppo4, & 478 492 & ln_pdic, ln_palk, ln_pph, ln_po2 /) ) 479 nsurftypes = COUNT( (/ln_sla, ln_sst, ln_sic, ln_sit, ln_fbd, ln_sss, & 493 nsurftypes = COUNT( (/ln_sla, ln_sst, ln_sic, ln_sit, ln_fbd, & 494 & ln_sss, ln_ssv, & 480 495 & ln_slchltot, ln_slchldia, ln_slchlnon, ln_slchldin, & 481 496 & ln_slchlmic, ln_slchlnan, ln_slchlpic, ln_schltot, & … … 611 626 cobstypessurf(jtype) = 'sss' 612 627 clsurffiles(jtype,:) = cn_sssfbfiles 628 ENDIF 629 IF (ln_ssv) THEN 630 jtype = jtype + 1 631 cobstypessurf(jtype) = 'ssv' 632 clsurffiles(jtype,:) = cn_ssvfbfiles 613 633 ENDIF 614 634 IF (ln_slchltot) THEN … … 751 771 ztype_avgphiscl = rn_sss_avgphiscl 752 772 ltype_fp_indegs = ln_sss_fp_indegs 773 ltype_night = .FALSE. 774 ELSE IF ( TRIM(cobstypessurf(jtype)) == 'ssv' ) THEN 775 IF ( nn_2dint_ssv == -1 ) THEN 776 n2dint_type = nn_2dint_default 777 ELSE 778 n2dint_type = nn_2dint_ssv 779 ENDIF 780 ztype_avglamscl = rn_ssv_avglamscl 781 ztype_avgphiscl = rn_ssv_avgphiscl 782 ltype_fp_indegs = ln_ssv_fp_indegs 753 783 ltype_night = .FALSE. 754 784 ELSE … … 934 964 nvarssurf(jtype) = 1 935 965 nextrsurf(jtype) = 2 966 ELSE IF ( TRIM(cobstypessurf(jtype)) == 'ssv' ) THEN 967 nvarssurf(jtype) = 2 968 nextrsurf(jtype) = 0 936 969 ELSE 937 970 nvarssurf(jtype) = 1 … … 940 973 941 974 ALLOCATE( clvars( nvarssurf(jtype) ) ) 942 975 CALL wrk_alloc( jpi, jpj, nvarssurf(jtype), zglam ) 976 CALL wrk_alloc( jpi, jpj, nvarssurf(jtype), zgphi ) 977 CALL wrk_alloc( jpi, jpj, nvarssurf(jtype), zmask_surf ) 978 979 IF ( TRIM(cobstypessurf(jtype)) == 'ssv' ) THEN 980 zglam(:,:,1) = glamu(:,:) 981 zglam(:,:,2) = glamv(:,:) 982 zgphi(:,:,1) = gphiu(:,:) 983 zgphi(:,:,2) = gphiv(:,:) 984 zmask_surf(:,:,1) = umask(:,:,1) 985 zmask_surf(:,:,2) = vmask(:,:,1) 986 ELSE 987 DO jvar = 1, nvarssurf(jtype) 988 zglam(:,:,jvar) = glamt(:,:) 989 zgphi(:,:,jvar) = gphit(:,:) 990 zmask_surf(:,:,jvar) = tmask(:,:,1) 991 END DO 992 ENDIF 993 943 994 IF ( TRIM(cobstypessurf(jtype)) == 'sla' ) THEN 944 995 clvars(1) = 'SLA' … … 955 1006 ELSE IF ( TRIM(cobstypessurf(jtype)) == 'sss' ) THEN 956 1007 clvars(1) = 'SSS' 1008 ELSE IF ( TRIM(cobstypessurf(jtype)) == 'ssv' ) THEN 1009 clvars(1) = 'UVEL' 1010 clvars(2) = 'VVEL' 957 1011 ELSE IF ( TRIM(cobstypessurf(jtype)) == 'slchltot' ) THEN 958 1012 clvars(1) = 'SLCHLTOT' … … 994 1048 & llnightav(jtype), ltype_clim, ln_time_mean_sla_bkg, clvars ) 995 1049 996 CALL obs_pre_surf( surfdata(jtype), surfdataqc(jtype), ln_nea, ln_bound_reject, ln_seaicetypes ) 1050 CALL obs_pre_surf( surfdata(jtype), surfdataqc(jtype), & 1051 & jpi, jpj, & 1052 & zmask_surf, zglam, zgphi, & 1053 & ln_nea, ln_bound_reject, ln_seaicetypes ) 997 1054 998 1055 IF ( TRIM(cobstypessurf(jtype)) == 'sla' ) THEN … … 1024 1081 1025 1082 DEALLOCATE( clvars ) 1083 CALL wrk_dealloc( jpi, jpj, nvarssurf(jtype), zglam ) 1084 CALL wrk_dealloc( jpi, jpj, nvarssurf(jtype), zgphi ) 1085 CALL wrk_dealloc( jpi, jpj, nvarssurf(jtype), zmask_surf ) 1026 1086 1027 1087 END DO … … 1104 1164 #elif defined key_fabm 1105 1165 USE par_fabm ! FABM parameters 1106 USE fabm, ONLY: &1107 & fabm_get_interior_diagnostic_data1108 1166 #endif 1109 1167 #if defined key_spm … … 1128 1186 REAL(wp), POINTER, DIMENSION(:,:,:,:) :: & 1129 1187 & zprofmask ! Mask associated with zprofvar 1130 REAL(wp), POINTER, DIMENSION(:,: ) :: &1188 REAL(wp), POINTER, DIMENSION(:,:,:) :: & 1131 1189 & zsurfvar, & ! Model values equivalent to surface ob. 1132 1190 & zsurfclim, & ! Climatology values for variables in a surface ob. 1133 1191 & zsurfmask ! Mask associated with surface variable 1134 1192 REAL(wp), POINTER, DIMENSION(:,:,:) :: & 1135 & zglam, & ! Model longitudes for prof variables1136 & zgphi ! Model latitudes for prof variables1193 & zglam, & ! Model longitudes 1194 & zgphi ! Model latitudes 1137 1195 LOGICAL :: llog10 ! Perform log10 transform of variable 1138 1196 #if defined key_fabm … … 1316 1374 #elif defined key_fabm 1317 1375 ! Alkalinity from ERSEM 1318 zprofvar(:,:,:,1) = fabm_get_interior_diagnostic_data(model,jp_fabm_o3ta)1376 zprofvar(:,:,:,1) = model%get_interior_diagnostic_data(jp_fabm_o3ta) 1319 1377 #else 1320 1378 CALL ctl_stop( ' Trying to run palk observation operator', & … … 1331 1389 #elif defined key_fabm 1332 1390 ! pH from ERSEM 1333 zprofvar(:,:,:,1) = trn(:,:,:,jp_fabm_m1+jp_fabm_o3ph)1391 zprofvar(:,:,:,1) = model%get_interior_diagnostic_data(jp_fabm_o3ph) 1334 1392 #else 1335 1393 CALL ctl_stop( ' Trying to run pph observation operator', & … … 1380 1438 1381 1439 IF ( nsurftypes > 0 ) THEN 1382 1383 !Allocate local work arrays 1384 CALL wrk_alloc( jpi, jpj, zsurfvar ) 1385 CALL wrk_alloc( jpi, jpj, zsurfclim ) 1386 CALL wrk_alloc( jpi, jpj, zsurfmask ) 1440 1387 1441 #if defined key_fabm 1388 1442 CALL wrk_alloc( jpi, jpj, jpk, fabm_3d ) … … 1391 1445 DO jtype = 1, nsurftypes 1392 1446 1447 !Allocate local work arrays 1448 CALL wrk_alloc( jpi, jpj, surfdataqc(jtype)%nvar, zsurfvar ) 1449 CALL wrk_alloc( jpi, jpj, surfdataqc(jtype)%nvar, zsurfclim ) 1450 CALL wrk_alloc( jpi, jpj, surfdataqc(jtype)%nvar, zsurfmask ) 1451 CALL wrk_alloc( jpi, jpj, surfdataqc(jtype)%nvar, zglam ) 1452 CALL wrk_alloc( jpi, jpj, surfdataqc(jtype)%nvar, zgphi ) 1453 1393 1454 !Defaults which might be changed 1394 zsurfmask(:,:) = tmask(:,:,1) 1395 zsurfclim(:,:) = 0._wp 1455 DO jvar = 1, surfdataqc(jtype)%nvar 1456 zsurfmask(:,:,jvar) = tmask(:,:,1) 1457 zsurfclim(:,:,jvar) = 0._wp 1458 zglam(:,:,jvar) = glamt(:,:) 1459 zgphi(:,:,jvar) = gphit(:,:) 1460 END DO 1396 1461 llog10 = .FALSE. 1397 1462 1398 1463 SELECT CASE ( TRIM(cobstypessurf(jtype)) ) 1399 1464 CASE('sst') 1400 zsurfvar(:,: ) = tsn(:,:,1,jp_tem)1401 IF ( ln_output_clim ) zsurfclim(:,: ) = tclim(:,:,1)1465 zsurfvar(:,:,1) = tsn(:,:,1,jp_tem) 1466 IF ( ln_output_clim ) zsurfclim(:,:,1) = tclim(:,:,1) 1402 1467 CASE('sla') 1403 zsurfvar(:,: ) = sshn(:,:)1468 zsurfvar(:,:,1) = sshn(:,:) 1404 1469 CASE('sss') 1405 zsurfvar(:,:) = tsn(:,:,1,jp_sal) 1406 IF ( ln_output_clim ) zsurfclim(:,:) = sclim(:,:,1) 1470 zsurfvar(:,:,1) = tsn(:,:,1,jp_sal) 1471 IF ( ln_output_clim ) zsurfclim(:,:,1) = sclim(:,:,1) 1472 CASE('ssv') 1473 zsurfvar(:,:,1) = un(:,:,1) 1474 zsurfvar(:,:,2) = vn(:,:,1) 1475 zsurfmask(:,:,1) = umask(:,:,1) 1476 zsurfmask(:,:,2) = vmask(:,:,1) 1477 zglam(:,:,1) = glamu(:,:) 1478 zglam(:,:,2) = glamv(:,:) 1479 zgphi(:,:,1) = gphiu(:,:) 1480 zgphi(:,:,2) = gphiv(:,:) 1407 1481 CASE('sic') 1408 1482 IF ( kstp == 0 ) THEN … … 1415 1489 ELSE 1416 1490 #if defined key_cice 1417 zsurfvar(:,: ) = fr_i(:,:)1491 zsurfvar(:,:,1) = fr_i(:,:) 1418 1492 #elif defined key_lim2 || defined key_lim3 1419 zsurfvar(:,: ) = 1._wp - frld(:,:)1493 zsurfvar(:,:,1) = 1._wp - frld(:,:) 1420 1494 #else 1421 1495 CALL ctl_stop( ' Trying to run sea-ice concentration observation operator', & … … 1434 1508 ELSE 1435 1509 #if defined key_cice 1436 zsurfvar(:,: ) = thick_i(:,:)1510 zsurfvar(:,:,1) = thick_i(:,:) 1437 1511 #elif defined key_lim2 || defined key_lim3 1438 1512 CALL ctl_stop( ' No sea-ice thickness observation operator defined for LIM model' ) … … 1465 1539 #if defined key_hadocc 1466 1540 ! Surface chlorophyll from HadOCC 1467 zsurfvar(:,: ) = HADOCC_CHL(:,:,1)1541 zsurfvar(:,:,1) = HADOCC_CHL(:,:,1) 1468 1542 #elif defined key_medusa 1469 1543 ! Add non-diatom and diatom surface chlorophyll from MEDUSA 1470 zsurfvar(:,: ) = trn(:,:,1,jpchn) + trn(:,:,1,jpchd)1544 zsurfvar(:,:,1) = trn(:,:,1,jpchn) + trn(:,:,1,jpchd) 1471 1545 #elif defined key_fabm 1472 1546 ! Add all surface chlorophyll groups from ERSEM 1473 zsurfvar(:,: ) = trn(:,:,1,jp_fabm_m1+jp_fabm_chl1) + trn(:,:,1,jp_fabm_m1+jp_fabm_chl2) + &1547 zsurfvar(:,:,1) = trn(:,:,1,jp_fabm_m1+jp_fabm_chl1) + trn(:,:,1,jp_fabm_m1+jp_fabm_chl2) + & 1474 1548 & trn(:,:,1,jp_fabm_m1+jp_fabm_chl3) + trn(:,:,1,jp_fabm_m1+jp_fabm_chl4) 1475 1549 #else … … 1485 1559 #elif defined key_medusa 1486 1560 ! Diatom surface chlorophyll from MEDUSA 1487 zsurfvar(:,: ) = trn(:,:,1,jpchd)1561 zsurfvar(:,:,1) = trn(:,:,1,jpchd) 1488 1562 #elif defined key_fabm 1489 1563 ! Diatom surface chlorophyll from ERSEM 1490 zsurfvar(:,: ) = trn(:,:,1,jp_fabm_m1+jp_fabm_chl1)1564 zsurfvar(:,:,1) = trn(:,:,1,jp_fabm_m1+jp_fabm_chl1) 1491 1565 #else 1492 1566 CALL ctl_stop( ' Trying to run slchldia observation operator', & … … 1501 1575 #elif defined key_medusa 1502 1576 ! Non-diatom surface chlorophyll from MEDUSA 1503 zsurfvar(:,: ) = trn(:,:,1,jpchn)1577 zsurfvar(:,:,1) = trn(:,:,1,jpchn) 1504 1578 #elif defined key_fabm 1505 1579 ! Add all non-diatom surface chlorophyll groups from ERSEM 1506 zsurfvar(:,: ) = trn(:,:,1,jp_fabm_m1+jp_fabm_chl2) + &1580 zsurfvar(:,:,1) = trn(:,:,1,jp_fabm_m1+jp_fabm_chl2) + & 1507 1581 & trn(:,:,1,jp_fabm_m1+jp_fabm_chl3) + trn(:,:,1,jp_fabm_m1+jp_fabm_chl4) 1508 1582 #else … … 1521 1595 #elif defined key_fabm 1522 1596 ! Dinoflagellate surface chlorophyll from ERSEM 1523 zsurfvar(:,: ) = trn(:,:,1,jp_fabm_m1+jp_fabm_chl4)1597 zsurfvar(:,:,1) = trn(:,:,1,jp_fabm_m1+jp_fabm_chl4) 1524 1598 #else 1525 1599 CALL ctl_stop( ' Trying to run slchldin observation operator', & … … 1537 1611 #elif defined key_fabm 1538 1612 ! Add diatom and dinoflagellate surface chlorophyll from ERSEM 1539 zsurfvar(:,: ) = trn(:,:,1,jp_fabm_m1+jp_fabm_chl1) + trn(:,:,1,jp_fabm_m1+jp_fabm_chl4)1613 zsurfvar(:,:,1) = trn(:,:,1,jp_fabm_m1+jp_fabm_chl1) + trn(:,:,1,jp_fabm_m1+jp_fabm_chl4) 1540 1614 #else 1541 1615 CALL ctl_stop( ' Trying to run slchlmic observation operator', & … … 1553 1627 #elif defined key_fabm 1554 1628 ! Nanophytoplankton surface chlorophyll from ERSEM 1555 zsurfvar(:,: ) = trn(:,:,1,jp_fabm_m1+jp_fabm_chl2)1629 zsurfvar(:,:,1) = trn(:,:,1,jp_fabm_m1+jp_fabm_chl2) 1556 1630 #else 1557 1631 CALL ctl_stop( ' Trying to run slchlnan observation operator', & … … 1569 1643 #elif defined key_fabm 1570 1644 ! Picophytoplankton surface chlorophyll from ERSEM 1571 zsurfvar(:,: ) = trn(:,:,1,jp_fabm_m1+jp_fabm_chl3)1645 zsurfvar(:,:,1) = trn(:,:,1,jp_fabm_m1+jp_fabm_chl3) 1572 1646 #else 1573 1647 CALL ctl_stop( ' Trying to run slchlpic observation operator', & … … 1579 1653 #if defined key_hadocc 1580 1654 ! Surface chlorophyll from HadOCC 1581 zsurfvar(:,: ) = HADOCC_CHL(:,:,1)1655 zsurfvar(:,:,1) = HADOCC_CHL(:,:,1) 1582 1656 #elif defined key_medusa 1583 1657 ! Add non-diatom and diatom surface chlorophyll from MEDUSA 1584 zsurfvar(:,: ) = trn(:,:,1,jpchn) + trn(:,:,1,jpchd)1658 zsurfvar(:,:,1) = trn(:,:,1,jpchn) + trn(:,:,1,jpchd) 1585 1659 #elif defined key_fabm 1586 1660 ! Add all surface chlorophyll groups from ERSEM 1587 zsurfvar(:,: ) = trn(:,:,1,jp_fabm_m1+jp_fabm_chl1) + trn(:,:,1,jp_fabm_m1+jp_fabm_chl2) + &1588 & trn(:,:,1,jp_fabm_m1+jp_fabm_chl3) + trn(:,:,1,jp_fabm_m1+jp_fabm_chl4)1661 zsurfvar(:,:,1) = trn(:,:,1,jp_fabm_m1+jp_fabm_chl1) + trn(:,:,1,jp_fabm_m1+jp_fabm_chl2) + & 1662 & trn(:,:,1,jp_fabm_m1+jp_fabm_chl3) + trn(:,:,1,jp_fabm_m1+jp_fabm_chl4) 1589 1663 #else 1590 1664 CALL ctl_stop( ' Trying to run schltot observation operator', & … … 1595 1669 #if defined key_hadocc 1596 1670 ! Surface phytoplankton nitrogen from HadOCC multiplied by C:N ratio 1597 zsurfvar(:,: ) = trn(:,:,1,jp_had_phy) * c2n_p1671 zsurfvar(:,:,1) = trn(:,:,1,jp_had_phy) * c2n_p 1598 1672 #elif defined key_medusa 1599 1673 ! Add non-diatom and diatom surface phytoplankton nitrogen from MEDUSA 1600 1674 ! multiplied by C:N ratio for each 1601 zsurfvar(:,: ) = (trn(:,:,1,jpphn) * xthetapn) + (trn(:,:,1,jpphd) * xthetapd)1675 zsurfvar(:,:,1) = (trn(:,:,1,jpphn) * xthetapn) + (trn(:,:,1,jpphd) * xthetapd) 1602 1676 #elif defined key_fabm 1603 1677 ! Add all surface phytoplankton carbon groups from ERSEM 1604 zsurfvar(:,: ) = trn(:,:,1,jp_fabm_m1+jp_fabm_p1c) + trn(:,:,1,jp_fabm_m1+jp_fabm_p2c) + &1678 zsurfvar(:,:,1) = trn(:,:,1,jp_fabm_m1+jp_fabm_p1c) + trn(:,:,1,jp_fabm_m1+jp_fabm_p2c) + & 1605 1679 & trn(:,:,1,jp_fabm_m1+jp_fabm_p3c) + trn(:,:,1,jp_fabm_m1+jp_fabm_p4c) 1606 1680 #else … … 1616 1690 #elif defined key_medusa 1617 1691 ! Diatom surface phytoplankton nitrogen from MEDUSA multiplied by C:N ratio 1618 zsurfvar(:,: ) = trn(:,:,1,jpphd) * xthetapd1692 zsurfvar(:,:,1) = trn(:,:,1,jpphd) * xthetapd 1619 1693 #elif defined key_fabm 1620 1694 ! Diatom surface phytoplankton carbon from ERSEM 1621 zsurfvar(:,: ) = trn(:,:,1,jp_fabm_m1+jp_fabm_p1c)1695 zsurfvar(:,:,1) = trn(:,:,1,jp_fabm_m1+jp_fabm_p1c) 1622 1696 #else 1623 1697 CALL ctl_stop( ' Trying to run slphydia observation operator', & … … 1632 1706 #elif defined key_medusa 1633 1707 ! Non-diatom surface phytoplankton nitrogen from MEDUSA multiplied by C:N ratio 1634 zsurfvar(:,: ) = trn(:,:,1,jpphn) * xthetapn1708 zsurfvar(:,:,1) = trn(:,:,1,jpphn) * xthetapn 1635 1709 #elif defined key_fabm 1636 1710 ! Add all non-diatom surface phytoplankton carbon groups from ERSEM 1637 zsurfvar(:,: ) = trn(:,:,1,jp_fabm_m1+jp_fabm_p2c) + &1638 & trn(:,:,1,jp_fabm_m1+jp_fabm_p3c) + trn(:,:,1,jp_fabm_m1+jp_fabm_p4c)1711 zsurfvar(:,:,1) = trn(:,:,1,jp_fabm_m1+jp_fabm_p2c) + & 1712 & trn(:,:,1,jp_fabm_m1+jp_fabm_p3c) + trn(:,:,1,jp_fabm_m1+jp_fabm_p4c) 1639 1713 #else 1640 1714 CALL ctl_stop( ' Trying to run slphynon observation operator', & … … 1645 1719 CASE('sspm') 1646 1720 #if defined key_spm 1647 zsurfvar(:,: ) = 0.01721 zsurfvar(:,:,1) = 0.0 1648 1722 DO jn = 1, jp_spm 1649 zsurfvar(:,: ) = zsurfvar(:,:) + trn(:,:,1,jn) ! sum SPM sizes1723 zsurfvar(:,:,1) = zsurfvar(:,:,1) + trn(:,:,1,jn) ! sum SPM sizes 1650 1724 END DO 1651 1725 #else … … 1662 1736 & ' but MEDUSA does not explicitly simulate Kd490' ) 1663 1737 #elif defined key_fabm 1664 ! light_xEPS diagnostic variable 1665 fabm_3d(:,:,:) = fabm_get_interior_diagnostic_data(model, jp_fabm_xeps) 1666 zsurfvar(:,:) = fabm_3d(:,:,1) 1738 ! light_Kd_band3 diagnostic variable if using spectral optical model 1739 ! light_xEPS diagnostic variable if using standard ERSEM light model 1740 IF ( jp_fabm_kd490 /= -1 ) THEN 1741 fabm_3d(:,:,:) = model%get_interior_diagnostic_data(jp_fabm_kd490) 1742 ELSEIF ( jp_fabm_xeps /= -1 ) THEN 1743 fabm_3d(:,:,:) = model%get_interior_diagnostic_data(jp_fabm_xeps) 1744 ELSE 1745 CALL ctl_stop( ' Trying to run skd490 observation operator', & 1746 & ' but cannot access Kd490 from ERSEM' ) 1747 ENDIF 1748 zsurfvar(:,:,1) = fabm_3d(:,:,1) 1667 1749 #else 1668 1750 CALL ctl_stop( ' Trying to run skd490 observation operator', & … … 1672 1754 CASE('sfco2') 1673 1755 #if defined key_hadocc 1674 zsurfvar(:,: ) = HADOCC_FCO2(:,:) ! fCO2 from HadOCC1756 zsurfvar(:,:,1) = HADOCC_FCO2(:,:) ! fCO2 from HadOCC 1675 1757 IF ( ( MINVAL( HADOCC_FCO2 ) == HADOCC_FILL_FLT ) .AND. & 1676 1758 & ( MAXVAL( HADOCC_FCO2 ) == HADOCC_FILL_FLT ) ) THEN 1677 1759 zsurfvar(:,:) = obfillflt 1678 zsurfmask(:,: ) = 01760 zsurfmask(:,:,1) = 0 1679 1761 CALL ctl_warn( ' HadOCC fCO2 values masked out for observation operator', & 1680 1762 & ' as HADOCC_FCO2(:,:) == HADOCC_FILL_FLT' ) 1681 1763 ENDIF 1682 1764 #elif defined key_medusa && defined key_roam 1683 zsurfvar(:,: ) = f2_fco2w(:,:)1765 zsurfvar(:,:,1) = f2_fco2w(:,:) 1684 1766 #elif defined key_fabm 1685 1767 ! First, get pCO2 from FABM 1686 fabm_3d(:,:,:) = fabm_get_interior_diagnostic_data(model,jp_fabm_o3pc)1687 zsurfvar(:,: ) = fabm_3d(:,:,1)1768 fabm_3d(:,:,:) = model%get_interior_diagnostic_data(jp_fabm_o3pc) 1769 zsurfvar(:,:,1) = fabm_3d(:,:,1) 1688 1770 ! Now, convert pCO2 to fCO2, based on SST in K. This follows the standard methodology of: 1689 1771 ! Pierrot et al. (2009), Recommendations for autonomous underway pCO2 measuring systems … … 1699 1781 ! because xCO2 in atm is ~0, and so this term will only affect the result to the 3rd decimal 1700 1782 ! place for typical values, and xCO2 would need to be approximated from pCO2 anyway. 1701 zsurfvar(:,: ) = zsurfvar(:,:) * EXP((-1636.75 + &1702 & 12.0408 * (tsn(:,:,1,jp_tem)+rt0) - &1703 & 0.0327957 * (tsn(:,:,1,jp_tem)+rt0)*(tsn(:,:,1,jp_tem)+rt0) + &1704 & 0.0000316528 * (tsn(:,:,1,jp_tem)+rt0)*(tsn(:,:,1,jp_tem)+rt0)*(tsn(:,:,1,jp_tem)+rt0) + &1705 & 2.0 * (57.7 - 0.118 * (tsn(:,:,1,jp_tem)+rt0))) / &1706 & (82.0578 * (tsn(:,:,1,jp_tem)+rt0)))1783 zsurfvar(:,:,1) = zsurfvar(:,:,1) * EXP((-1636.75 + & 1784 & 12.0408 * (tsn(:,:,1,jp_tem)+rt0) - & 1785 & 0.0327957 * (tsn(:,:,1,jp_tem)+rt0)*(tsn(:,:,1,jp_tem)+rt0) + & 1786 & 0.0000316528 * (tsn(:,:,1,jp_tem)+rt0)*(tsn(:,:,1,jp_tem)+rt0)*(tsn(:,:,1,jp_tem)+rt0) + & 1787 & 2.0 * (57.7 - 0.118 * (tsn(:,:,1,jp_tem)+rt0))) / & 1788 & (82.0578 * (tsn(:,:,1,jp_tem)+rt0))) 1707 1789 #else 1708 1790 CALL ctl_stop( ' Trying to run sfco2 observation operator', & … … 1712 1794 CASE('spco2') 1713 1795 #if defined key_hadocc 1714 zsurfvar(:,: ) = HADOCC_PCO2(:,:) ! pCO2 from HadOCC1796 zsurfvar(:,:,1) = HADOCC_PCO2(:,:) ! pCO2 from HadOCC 1715 1797 IF ( ( MINVAL( HADOCC_PCO2 ) == HADOCC_FILL_FLT ) .AND. & 1716 1798 & ( MAXVAL( HADOCC_PCO2 ) == HADOCC_FILL_FLT ) ) THEN 1717 zsurfvar(:,: ) = obfillflt1718 zsurfmask(:,: ) = 01799 zsurfvar(:,:,1) = obfillflt 1800 zsurfmask(:,:,1) = 0 1719 1801 CALL ctl_warn( ' HadOCC pCO2 values masked out for observation operator', & 1720 1802 & ' as HADOCC_PCO2(:,:) == HADOCC_FILL_FLT' ) 1721 1803 ENDIF 1722 1804 #elif defined key_medusa && defined key_roam 1723 zsurfvar(:,: ) = f2_pco2w(:,:)1724 #elif defined key_fabm 1725 fabm_3d(:,:,:) = fabm_get_interior_diagnostic_data(model,jp_fabm_o3pc)1726 zsurfvar(:,: ) = fabm_3d(:,:,1)1805 zsurfvar(:,:,1) = f2_pco2w(:,:) 1806 #elif defined key_fabm 1807 fabm_3d(:,:,:) = model%get_interior_diagnostic_data(jp_fabm_o3pc) 1808 zsurfvar(:,:,1) = fabm_3d(:,:,1) 1727 1809 #else 1728 1810 CALL ctl_stop( ' Trying to run spco2 observation operator', & … … 1739 1821 ! Take the log10 where we can, otherwise exclude 1740 1822 tiny = 1.0e-20 1741 WHERE(zsurfvar(:,: ) > tiny .AND. zsurfvar(:,:) /= obfillflt )1742 zsurfvar(:,: ) = LOG10(zsurfvar(:,:))1823 WHERE(zsurfvar(:,:,1) > tiny .AND. zsurfvar(:,:,1) /= obfillflt ) 1824 zsurfvar(:,:,1) = LOG10(zsurfvar(:,:,1)) 1743 1825 ELSEWHERE 1744 zsurfvar(:,: ) = obfillflt1745 zsurfmask(:,: ) = 01826 zsurfvar(:,:,1) = obfillflt 1827 zsurfmask(:,:,1) = 0 1746 1828 END WHERE 1747 1829 ENDIF 1748 1830 1749 IF ( TRIM(cobstypessurf(jtype)) == 'sla' .AND. & 1750 & ln_time_mean_sla_bkg ) THEN 1751 !Number of time-steps in meaning period 1752 imeanstp = NINT( ( MeanPeriodHours * 60. * 60. ) / rdt ) 1753 CALL obs_surf_opt( surfdataqc(jtype), kstp, jpi, jpj, & 1754 & nit000, idaystp, zsurfvar, & 1755 & zsurfclim, zsurfmask, & 1756 & n2dintsurf(jtype), llnightav(jtype), & 1757 & ravglamscl(jtype), ravgphiscl(jtype), & 1758 & lfpindegs(jtype), kmeanstp = imeanstp ) 1759 1760 1761 ELSE 1762 CALL obs_surf_opt( surfdataqc(jtype), kstp, jpi, jpj, & 1763 & nit000, idaystp, zsurfvar, & 1764 & zsurfclim, zsurfmask, & 1765 & n2dintsurf(jtype), llnightav(jtype), & 1766 & ravglamscl(jtype), ravgphiscl(jtype), & 1767 & lfpindegs(jtype) ) 1768 ENDIF 1831 DO jvar = 1, surfdataqc(jtype)%nvar 1832 1833 IF ( TRIM(cobstypessurf(jtype)) == 'sla' .AND. & 1834 & ln_time_mean_sla_bkg ) THEN 1835 !Number of time-steps in meaning period 1836 imeanstp = NINT( ( MeanPeriodHours * 60. * 60. ) / rdt ) 1837 CALL obs_surf_opt( surfdataqc(jtype), kstp, jpi, jpj, & 1838 & nit000, idaystp, jvar, & 1839 & zsurfvar(:,:,jvar), & 1840 & zsurfclim(:,:,jvar), & 1841 & zsurfmask(:,:,jvar), & 1842 & zglam(:,:,jvar), zgphi(:,:,jvar), & 1843 & n2dintsurf(jtype), llnightav(jtype), & 1844 & ravglamscl(jtype), ravgphiscl(jtype), & 1845 & lfpindegs(jtype), kmeanstp = imeanstp ) 1846 1847 ELSE 1848 CALL obs_surf_opt( surfdataqc(jtype), kstp, jpi, jpj, & 1849 & nit000, idaystp, jvar, & 1850 & zsurfvar(:,:,jvar), & 1851 & zsurfclim(:,:,jvar), & 1852 & zsurfmask(:,:,jvar), & 1853 & zglam(:,:,jvar), zgphi(:,:,jvar), & 1854 & n2dintsurf(jtype), llnightav(jtype), & 1855 & ravglamscl(jtype), ravgphiscl(jtype), & 1856 & lfpindegs(jtype) ) 1857 ENDIF 1858 1859 END DO 1860 1861 CALL wrk_dealloc( jpi, jpj, surfdataqc(jtype)%nvar, zsurfvar ) 1862 CALL wrk_dealloc( jpi, jpj, surfdataqc(jtype)%nvar, zsurfclim ) 1863 CALL wrk_dealloc( jpi, jpj, surfdataqc(jtype)%nvar, zsurfmask ) 1864 CALL wrk_dealloc( jpi, jpj, surfdataqc(jtype)%nvar, zglam ) 1865 CALL wrk_dealloc( jpi, jpj, surfdataqc(jtype)%nvar, zgphi ) 1769 1866 1770 1867 END DO 1771 1772 CALL wrk_dealloc( jpi, jpj, zsurfvar )1773 CALL wrk_dealloc( jpi, jpj, zsurfmask )1774 1868 #if defined key_fabm 1775 1869 CALL wrk_dealloc( jpi, jpj, jpk, fabm_3d ) … … 1827 1921 & ) 1828 1922 1829 CALL obs_rotvel ( profdataqc(jtype), nn_2dint_default, zu, zv )1923 CALL obs_rotvel_pro( profdataqc(jtype), nn_2dint_default, zu, zv ) 1830 1924 1831 1925 DO jo = 1, profdataqc(jtype)%nprof … … 1860 1954 1861 1955 DO jtype = 1, nsurftypes 1956 1957 IF ( TRIM(cobstypessurf(jtype)) == 'vel' ) THEN 1958 1959 ! For velocity data, rotate the model velocities to N/S, E/W 1960 ! using the compressed data structure. 1961 ALLOCATE( & 1962 & zu(surfdataqc(jtype)%nsurf), & 1963 & zv(surfdataqc(jtype)%nsurf) & 1964 & ) 1965 1966 CALL obs_rotvel_surf( surfdataqc(jtype), nn_2dint_default, zu, zv ) 1967 1968 DO jo = 1, surfdataqc(jtype)%nsurf 1969 surfdataqc(jtype)%rmod(jo,1) = zu(jo) 1970 surfdataqc(jtype)%rmod(jo,2) = zv(jo) 1971 END DO 1972 1973 DEALLOCATE( zu ) 1974 DEALLOCATE( zv ) 1975 1976 END IF 1977 1862 1978 1863 1979 CALL obs_surf_decompress( surfdataqc(jtype), & -
branches/UKMO/dev_r5518_obs_oper_update_sit_SMOS/NEMOGCM/NEMO/OPA_SRC/OBS/obs_oper.F90
r12610 r14165 544 544 545 545 SUBROUTINE obs_surf_opt( surfdataqc, kt, kpi, kpj, & 546 & kit000, kdaystp, psurf, pclim, psurfmask, & 546 & kit000, kdaystp, kvar, & 547 & psurf, pclim, psurfmask, & 548 & plam, pphi, & 547 549 & k2dint, ldnightav, plamscl, pphiscl, & 548 550 & lindegrees, kmeanstp ) … … 579 581 !! ! 15-02 (M. Martin) Combined routine for surface types 580 582 !! ! 17-03 (M. Martin) Added horizontal averaging options 583 !! ! 20-08 (M. Martin) Added surface velocity options 581 584 !!----------------------------------------------------------------------- 582 585 … … 595 598 ! (kit000-1 = restart time) 596 599 INTEGER, INTENT(IN) :: kdaystp ! Number of time steps per day 600 INTEGER, INTENT(IN) :: kvar ! Number of variables in surfdataqc 597 601 INTEGER, INTENT(IN) :: k2dint ! Horizontal interpolation type (see header) 598 602 INTEGER, INTENT(IN), OPTIONAL :: & … … 603 607 & pclim, & ! Climatological surface field 604 608 & psurfmask ! Land-sea mask 609 REAL(KIND=wp), INTENT(IN), DIMENSION(kpi,kpj) :: & 610 & plam, & ! Model longitudes for variable 611 & pphi ! Model latitudes for variable 605 612 LOGICAL, INTENT(IN) :: ldnightav ! Logical for averaging night-time data 606 613 REAL(KIND=wp), INTENT(IN) :: & … … 670 677 CALL obs_max_fpsize( k2dint, plamscl, pphiscl, lindegrees, psurfmask, imaxifp, imaxjfp ) 671 678 672 673 679 IF ( l_timemean ) THEN 674 680 ! Initialize time mean for first timestep … … 681 687 DO jj = 1, jpj 682 688 DO ji = 1, jpi 683 surfdataqc%vdmean(ji,jj ) = 0.0689 surfdataqc%vdmean(ji,jj,kvar) = 0.0 684 690 END DO 685 691 END DO … … 690 696 DO jj = 1, jpj 691 697 DO ji = 1, jpi 692 surfdataqc%vdmean(ji,jj ) = surfdataqc%vdmean(ji,jj) &693 & + psurf(ji,jj)698 surfdataqc%vdmean(ji,jj,kvar) = surfdataqc%vdmean(ji,jj,kvar) & 699 & + psurf(ji,jj) 694 700 END DO 695 701 END DO … … 701 707 DO jj = 1, jpj 702 708 DO ji = 1, jpi 703 surfdataqc%vdmean(ji,jj ) = surfdataqc%vdmean(ji,jj) &704 & * zmeanstp709 surfdataqc%vdmean(ji,jj,kvar) = surfdataqc%vdmean(ji,jj,kvar) & 710 & * zmeanstp 705 711 END DO 706 712 END DO … … 729 735 DO jj = 1, jpj 730 736 DO ji = 1, jpi 731 surfdataqc%vdmean(ji,jj ) = 0.0737 surfdataqc%vdmean(ji,jj,kvar) = 0.0 732 738 zmeanday(ji,jj) = 0.0 733 739 icount_night(ji,jj) = 0 … … 743 749 DO ji = 1, jpi 744 750 ! Increment the temperature field for computing night mean and counter 745 surfdataqc%vdmean(ji,jj ) = surfdataqc%vdmean(ji,jj) &746 & + psurf(ji,jj) * REAL( imask_night(ji,jj) )751 surfdataqc%vdmean(ji,jj,kvar) = surfdataqc%vdmean(ji,jj,kvar) & 752 & + psurf(ji,jj) * REAL( imask_night(ji,jj) ) 747 753 zmeanday(ji,jj) = zmeanday(ji,jj) + psurf(ji,jj) 748 754 icount_night(ji,jj) = icount_night(ji,jj) + imask_night(ji,jj) … … 758 764 ! Test if "no night" point 759 765 IF ( icount_night(ji,jj) > 0 ) THEN 760 surfdataqc%vdmean(ji,jj ) = surfdataqc%vdmean(ji,jj) &761 & / REAL( icount_night(ji,jj) )766 surfdataqc%vdmean(ji,jj,kvar) = surfdataqc%vdmean(ji,jj,kvar) & 767 & / REAL( icount_night(ji,jj) ) 762 768 ELSE 763 769 !At locations where there is no night (e.g. poles), 764 770 ! calculate daily mean instead of night-time mean. 765 surfdataqc%vdmean(ji,jj ) = zmeanday(ji,jj) * zdaystp771 surfdataqc%vdmean(ji,jj,kvar) = zmeanday(ji,jj) * zdaystp 766 772 ENDIF 767 773 END DO … … 772 778 773 779 ! Get the data for interpolation 774 775 780 ALLOCATE( & 776 & zweig(imaxifp,imaxjfp,1), & 777 & igrdi(imaxifp,imaxjfp,isurf), & 778 & igrdj(imaxifp,imaxjfp,isurf), & 779 & zglam(imaxifp,imaxjfp,isurf), & 780 & zgphi(imaxifp,imaxjfp,isurf), & 781 & zmask(imaxifp,imaxjfp,isurf), & 782 & zsurf(imaxifp,imaxjfp,isurf), & 783 & zsurftmp(imaxifp,imaxjfp,isurf), & 784 & zglamf(imaxifp+1,imaxjfp+1,isurf), & 785 & zgphif(imaxifp+1,imaxjfp+1,isurf), & 786 & igrdip1(imaxifp+1,imaxjfp+1,isurf), & 787 & igrdjp1(imaxifp+1,imaxjfp+1,isurf) & 781 & zweig(imaxifp,imaxjfp,1), & 782 & igrdi(imaxifp,imaxjfp,isurf), & 783 & igrdj(imaxifp,imaxjfp,isurf), & 784 & zglam(imaxifp,imaxjfp,isurf), & 785 & zgphi(imaxifp,imaxjfp,isurf), & 786 & zmask(imaxifp,imaxjfp,isurf), & 787 & zsurf(imaxifp,imaxjfp,isurf), & 788 & zsurftmp(imaxifp,imaxjfp,isurf) & 788 789 & ) 790 791 IF ( k2dint > 4 ) THEN 792 ALLOCATE( & 793 & zglamf(imaxifp+1,imaxjfp+1,isurf), & 794 & zgphif(imaxifp+1,imaxjfp+1,isurf), & 795 & igrdip1(imaxifp+1,imaxjfp+1,isurf), & 796 & igrdjp1(imaxifp+1,imaxjfp+1,isurf) & 797 & ) 798 ENDIF 789 799 790 800 IF ( surfdataqc%lclim ) ALLOCATE( zclim(imaxifp,imaxjfp,isurf) ) … … 793 803 iobs = jobs - surfdataqc%nsurfup 794 804 DO ji = 0, imaxifp 795 imodi = surfdataqc%mi(jobs ) - int(imaxifp/2) + ji - 1805 imodi = surfdataqc%mi(jobs,kvar) - int(imaxifp/2) + ji - 1 796 806 797 807 !Deal with wrap around in longitude … … 800 810 801 811 DO jj = 0, imaxjfp 802 imodj = surfdataqc%mj(jobs ) - int(imaxjfp/2) + jj - 1812 imodj = surfdataqc%mj(jobs,kvar) - int(imaxjfp/2) + jj - 1 803 813 !If model values are out of the domain to the north/south then 804 814 !set them to be the edge of the domain … … 806 816 IF ( imodj > jpjglo ) imodj = jpjglo 807 817 808 igrdip1(ji+1,jj+1,iobs) = imodi 809 igrdjp1(ji+1,jj+1,iobs) = imodj 818 IF ( k2dint > 4 ) THEN 819 igrdip1(ji+1,jj+1,iobs) = imodi 820 igrdjp1(ji+1,jj+1,iobs) = imodj 821 ENDIF 810 822 811 823 IF ( ji >= 1 .AND. jj >= 1 ) THEN … … 819 831 820 832 CALL obs_int_comm_2d( imaxifp, imaxjfp, isurf, kpi, kpj, & 821 & igrdi, igrdj, glamt, zglam )833 & igrdi, igrdj, plam, zglam ) 822 834 CALL obs_int_comm_2d( imaxifp, imaxjfp, isurf, kpi, kpj, & 823 & igrdi, igrdj, gphit, zgphi )835 & igrdi, igrdj, pphi, zgphi ) 824 836 CALL obs_int_comm_2d( imaxifp, imaxjfp, isurf, kpi, kpj, & 825 837 & igrdi, igrdj, psurfmask, zmask ) … … 831 843 IF (lwp) WRITE(numout,*)' Interpolating the time mean values on time step: ',kt 832 844 CALL obs_int_comm_2d( imaxifp, imaxjfp, isurf, kpi, kpj, & 833 & igrdi, igrdj, surfdataqc%vdmean(:,: ), zsurfm )845 & igrdi, igrdj, surfdataqc%vdmean(:,:,kvar), zsurfm ) 834 846 ENDIF 835 847 ELSE … … 858 870 859 871 CALL obs_int_comm_2d( imaxifp,imaxjfp, isurf, kpi, kpj, igrdi, igrdj, & 860 & surfdataqc%vdmean(:,:), zsurfm )872 & surfdataqc%vdmean(:,:,kvar), zsurfm ) 861 873 862 874 ENDIF … … 937 949 surfdataqc%rmod(jobs,1) = surfdataqc%rext(jobs,1) - surfdataqc%rext(jobs,2) 938 950 ELSE 939 surfdataqc%rmod(jobs, 1) = zext(1)951 surfdataqc%rmod(jobs,kvar) = zext(1) 940 952 ENDIF 941 953 … … 985 997 & zmask, & 986 998 & zsurf, & 987 & zsurftmp, & 988 & zglamf, & 989 & zgphif, & 990 & igrdip1,& 991 & igrdjp1 & 999 & zsurftmp & 992 1000 & ) 993 1001 1002 IF ( k2dint > 4 ) THEN 1003 DEALLOCATE( & 1004 & zglamf, & 1005 & zgphif, & 1006 & igrdip1,& 1007 & igrdjp1 & 1008 & ) 1009 ENDIF 1010 994 1011 IF ( surfdataqc%lclim ) DEALLOCATE( zclim ) 995 1012 … … 1001 1018 ENDIF 1002 1019 1003 surfdataqc%nsurfup = surfdataqc%nsurfup + isurf 1020 IF ( kvar == surfdataqc%nvar ) THEN 1021 surfdataqc%nsurfup = surfdataqc%nsurfup + isurf 1022 ENDIF 1004 1023 1005 1024 END SUBROUTINE obs_surf_opt -
branches/UKMO/dev_r5518_obs_oper_update_sit_SMOS/NEMOGCM/NEMO/OPA_SRC/OBS/obs_prep.F90
r12610 r14165 53 53 CONTAINS 54 54 55 SUBROUTINE obs_pre_surf( surfdata, surfdataqc, ld_nea, ld_bound_reject, & 56 ld_seaicetypes, kqc_cutoff ) 57 !!---------------------------------------------------------------------- 58 !! *** ROUTINE obs_pre_sla *** 55 SUBROUTINE obs_pre_surf( surfdata, surfdataqc, & 56 & kpi, kpj, & 57 & zmask, pglam, pgphi, & 58 & ld_nea, ld_bound_reject, & 59 & ld_seaicetypes, kqc_cutoff ) 60 !!---------------------------------------------------------------------- 61 !! *** ROUTINE obs_pre_surf *** 59 62 !! 60 63 !! ** Purpose : First level check and screening of surface observations … … 82 85 TYPE(obs_surf), INTENT(INOUT) :: surfdata ! Full set of surface data 83 86 TYPE(obs_surf), INTENT(INOUT) :: surfdataqc ! Subset of surface data not failing screening 87 INTEGER, INTENT(IN) :: kpi, kpj ! Local domain sizes 88 REAL(wp), INTENT(IN), DIMENSION(kpi,kpj,surfdata%nvar) :: & 89 & zmask 90 REAL(wp), INTENT(IN), DIMENSION(kpi,kpj,surfdata%nvar) :: & 91 & pglam, & 92 & pgphi 84 93 LOGICAL, INTENT(IN) :: ld_nea ! Switch for rejecting observation near land 85 94 LOGICAL, INTENT(IN) :: ld_bound_reject ! Switch for rejecting obs near the boundary … … 94 103 INTEGER :: imin0 95 104 INTEGER :: icycle ! Current assimilation cycle 96 ! Counters for observations that 97 INTEGER :: iotdobs ! - outside time domain 98 INTEGER :: iosdsobs ! - outside space domain 99 INTEGER :: ilansobs ! - within a model land cell 100 INTEGER :: inlasobs ! - close to land 101 INTEGER :: igrdobs ! - fail the grid search 102 INTEGER :: ibdysobs ! - close to open boundary 103 ! Global counters for observations that 104 INTEGER :: iotdobsmpp ! - outside time domain 105 INTEGER :: iosdsobsmpp ! - outside space domain 106 INTEGER :: ilansobsmpp ! - within a model land cell 107 INTEGER :: inlasobsmpp ! - close to land 108 INTEGER :: igrdobsmpp ! - fail the grid search 109 INTEGER :: ibdysobsmpp ! - close to open boundary 105 ! Counters for observations that are 106 INTEGER :: iotdobs ! - outside time domain 107 INTEGER, DIMENSION(surfdata%nvar) :: iosdsobs ! - outside space domain 108 INTEGER, DIMENSION(surfdata%nvar) :: ilansobs ! - within a model land cell 109 INTEGER, DIMENSION(surfdata%nvar) :: inlasobs ! - close to land 110 INTEGER, DIMENSION(surfdata%nvar) :: ibdysobs ! - close to open boundary 111 INTEGER :: igrdobs ! - fail the grid search 112 ! Global counters for observations that 113 INTEGER :: iotdobsmpp ! - outside time domain 114 INTEGER, DIMENSION(surfdata%nvar) :: iosdsobsmpp ! - outside space domain 115 INTEGER, DIMENSION(surfdata%nvar) :: ilansobsmpp ! - within a model land cell 116 INTEGER, DIMENSION(surfdata%nvar) :: inlasobsmpp ! - close to land 117 INTEGER, DIMENSION(surfdata%nvar) :: ibdysobsmpp ! - close to open boundary 118 INTEGER :: igrdobsmpp ! - fail the grid search 119 110 120 LOGICAL, DIMENSION(:), ALLOCATABLE :: & 111 121 & llvalid ! SLA data selection 112 INTEGER :: jobs ! Obs. loop variable 122 INTEGER :: jobs ! Obs. loop counter 123 INTEGER :: jvar ! Variable loop counter 113 124 INTEGER :: jstp ! Time loop variable 114 125 INTEGER :: inrc ! Time index variable 115 126 CHARACTER(LEN=256) :: cout1 ! Diagnostic output line 127 116 128 IF(lwp) WRITE(numout,*)'obs_pre_surf : Preparing the surface observations...' 117 129 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' … … 130 142 iotdobs = 0 131 143 igrdobs = 0 132 iosdsobs = 0133 ilansobs = 0134 inlasobs = 0135 ibdysobs = 0144 iosdsobs(:) = 0 145 ilansobs(:) = 0 146 inlasobs(:) = 0 147 ibdysobs(:) = 0 136 148 137 149 ! Set QC cutoff to optional value if provided … … 162 174 ! Check for surface data failing the grid search 163 175 ! ----------------------------------------------------------------------- 164 165 CALL obs_coo_grd( surfdata%nsurf, surfdata%mi, surfdata%mj, & 166 & surfdata%nqc, igrdobs ) 167 176 DO jvar = 1, surfdata%nvar 177 CALL obs_coo_grd( surfdata%nsurf, surfdata%mi(:,jvar), surfdata%mj(:,jvar), & 178 & surfdata%nqc, igrdobs ) 179 END DO 180 168 181 CALL obs_mpp_sum_integer( igrdobs, igrdobsmpp ) 169 182 … … 172 185 ! ----------------------------------------------------------------------- 173 186 174 CALL obs_coo_spc_2d( surfdata%nsurf, & 175 & jpi, jpj, & 176 & surfdata%mi, surfdata%mj, & 177 & surfdata%rlam, surfdata%rphi, & 178 & glamt, gphit, & 179 & tmask(:,:,1), surfdata%nqc, & 180 & iosdsobs, ilansobs, & 181 & inlasobs, ld_nea, & 182 & ibdysobs, ld_bound_reject, & 183 & iqc_cutoff ) 184 185 CALL obs_mpp_sum_integer( iosdsobs, iosdsobsmpp ) 186 CALL obs_mpp_sum_integer( ilansobs, ilansobsmpp ) 187 CALL obs_mpp_sum_integer( inlasobs, inlasobsmpp ) 188 CALL obs_mpp_sum_integer( ibdysobs, ibdysobsmpp ) 189 187 DO jvar = 1, surfdata%nvar 188 CALL obs_coo_spc_2d( surfdata%nsurf, & 189 & jpi, jpj, & 190 & surfdata%mi(:,jvar), surfdata%mj(:,jvar), & 191 & surfdata%rlam, surfdata%rphi, & 192 & pglam(:,:,jvar), pgphi(:,:,jvar), & 193 & zmask(:,:,jvar), surfdata%nqc(:), & 194 & iosdsobs(jvar), ilansobs(jvar), & 195 & inlasobs(jvar), ld_nea, & 196 & ibdysobs(jvar), ld_bound_reject, & 197 & iqc_cutoff ) 198 CALL obs_mpp_sum_integer( iosdsobs(jvar), iosdsobsmpp(jvar) ) 199 CALL obs_mpp_sum_integer( ilansobs(jvar), ilansobsmpp(jvar) ) 200 CALL obs_mpp_sum_integer( inlasobs(jvar), inlasobsmpp(jvar) ) 201 CALL obs_mpp_sum_integer( ibdysobs(jvar), ibdysobsmpp(jvar) ) 202 END DO 203 190 204 ! ----------------------------------------------------------------------- 191 205 ! Copy useful data from the surfdata data structure to … … 216 230 217 231 IF(lwp) THEN 232 233 DO jvar = 1, surfdataqc%nvar 234 IF ( jvar == 1 ) THEN 235 cout1=TRIM(surfdataqc%cvars(1)) 236 ELSE 237 WRITE(cout1,'(A,A1,A)') TRIM(cout1), '/', TRIM(surfdataqc%cvars(jvar)) 238 ENDIF 239 END DO 240 218 241 WRITE(numout,*) 219 WRITE(numout,*) ' '// surfdataqc%cvars(1)//' data outside time domain = ', &242 WRITE(numout,*) ' '//TRIM(cout1)//' data outside time domain = ', & 220 243 & iotdobsmpp 221 WRITE(numout,*) ' Remaining '// surfdataqc%cvars(1)//' data that failed grid search = ', &244 WRITE(numout,*) ' Remaining '//TRIM(cout1)//' data that failed grid search = ', & 222 245 & igrdobsmpp 223 WRITE(numout,*) ' Remaining '//surfdataqc%cvars(1)//' data outside space domain = ', & 224 & iosdsobsmpp 225 WRITE(numout,*) ' Remaining '//surfdataqc%cvars(1)//' data at land points = ', & 226 & ilansobsmpp 227 IF (ld_nea) THEN 228 WRITE(numout,*) ' Remaining '//surfdataqc%cvars(1)//' data near land points (removed) = ', & 229 & inlasobsmpp 230 ELSE 231 WRITE(numout,*) ' Remaining '//surfdataqc%cvars(1)//' data near land points (kept) = ', & 232 & inlasobsmpp 233 ENDIF 234 WRITE(numout,*) ' Remaining '//surfdataqc%cvars(1)//' data near open boundary (removed) = ', & 235 & ibdysobsmpp 236 WRITE(numout,*) ' '//surfdataqc%cvars(1)//' data accepted = ', & 237 & surfdataqc%nsurfmpp 246 247 DO jvar = 1, surfdataqc%nvar 248 WRITE(numout,*) ' Remaining '//surfdataqc%cvars(jvar)//' data outside space domain = ', & 249 & iosdsobsmpp(jvar) 250 WRITE(numout,*) ' Remaining '//surfdataqc%cvars(jvar)//' data at land points = ', & 251 & ilansobsmpp(jvar) 252 IF (ld_nea) THEN 253 WRITE(numout,*) ' Remaining '//surfdataqc%cvars(jvar)//' data near land points (removed) = ', & 254 & inlasobsmpp(jvar) 255 ELSE 256 WRITE(numout,*) ' Remaining '//surfdataqc%cvars(jvar)//' data near land points (kept) = ', & 257 & inlasobsmpp(jvar) 258 ENDIF 259 WRITE(numout,*) ' Remaining '//surfdataqc%cvars(jvar)//' data near open boundary (removed) = ', & 260 & ibdysobsmpp(jvar) 261 END DO 262 WRITE(numout,*) ' '//TRIM(cout1)//' data accepted = ', & 263 & surfdataqc%nsurfmpp 238 264 239 265 WRITE(numout,*) 240 266 WRITE(numout,*) ' Number of observations per time step :' 241 267 WRITE(numout,*) 242 WRITE(numout,'(10X,A,10X,A)')'Time step', surfdataqc%cvars(1)268 WRITE(numout,'(10X,A,10X,A)')'Time step',TRIM(cout1) 243 269 WRITE(numout,'(10X,A,5X,A)')'---------','-----------------' 244 270 CALL FLUSH(numout) … … 445 471 446 472 IF ( TRIM(profdata%cvars(1)) == 'UVEL' ) THEN 447 CALL obs_uv_rej ( profdata, iuvchku, iuvchkv, iqc_cutoff )473 CALL obs_uv_rej_pro( profdata, iuvchku, iuvchkv, iqc_cutoff ) 448 474 CALL obs_mpp_sum_integer( iuvchku, iuvchkumpp ) 449 475 CALL obs_mpp_sum_integer( iuvchkv, iuvchkvmpp ) … … 1457 1483 END SUBROUTINE obs_pro_rej 1458 1484 1459 SUBROUTINE obs_uv_rej ( profdata, knumu, knumv, kqc_cutoff )1460 !!---------------------------------------------------------------------- 1461 !! *** ROUTINE obs_uv_rej ***1485 SUBROUTINE obs_uv_rej_pro( profdata, knumu, knumv, kqc_cutoff ) 1486 !!---------------------------------------------------------------------- 1487 !! *** ROUTINE obs_uv_rej_pro *** 1462 1488 !! 1463 1489 !! ** Purpose : Reject u if v is rejected and vice versa … … 1513 1539 END DO 1514 1540 1515 END SUBROUTINE obs_uv_rej 1541 END SUBROUTINE obs_uv_rej_pro 1516 1542 1517 1543 END MODULE obs_prep -
branches/UKMO/dev_r5518_obs_oper_update_sit_SMOS/NEMOGCM/NEMO/OPA_SRC/OBS/obs_read_surf.F90
r12610 r14165 85 85 INTEGER :: jj 86 86 INTEGER :: jk 87 INTEGER :: jvar 87 88 INTEGER :: iflag 88 89 INTEGER :: inobf … … 107 108 & ityp, & 108 109 & itypmpp 109 INTEGER, DIMENSION(: ), ALLOCATABLE :: &110 INTEGER, DIMENSION(:,:), ALLOCATABLE :: & 110 111 & iobsi, & 111 112 & iobsj, & 112 & iproc, & 113 & iproc 114 INTEGER, DIMENSION(:), ALLOCATABLE :: & 113 115 & iindx, & 114 116 & ifileidx, & … … 126 128 TYPE(obfbdata), POINTER, DIMENSION(:) :: & 127 129 & inpfiles 128 130 CHARACTER(LEN=256) :: cout1 ! Diagnostic output line 131 129 132 ! Local initialization 130 133 iobs = 0 … … 278 281 279 282 IF ( inpfiles(jj)%nobs > 0 ) THEN 280 inpfiles(jj)%iproc = -1281 inpfiles(jj)%iobsi = -1282 inpfiles(jj)%iobsj = -1283 inpfiles(jj)%iproc(:,:) = -1 284 inpfiles(jj)%iobsi(:,:) = -1 285 inpfiles(jj)%iobsj(:,:) = -1 283 286 ENDIF 284 287 … … 295 298 END DO 296 299 ENDIF 297 300 298 301 inowin = 0 299 302 DO ji = 1, inpfiles(jj)%nobs … … 305 308 ALLOCATE( zlam(inowin) ) 306 309 ALLOCATE( zphi(inowin) ) 307 ALLOCATE( iobsi(inowin ) )308 ALLOCATE( iobsj(inowin ) )309 ALLOCATE( iproc(inowin ) )310 ALLOCATE( iobsi(inowin,kvars) ) 311 ALLOCATE( iobsj(inowin,kvars) ) 312 ALLOCATE( iproc(inowin,kvars) ) 310 313 inowin = 0 311 314 DO ji = 1, inpfiles(jj)%nobs … … 318 321 END DO 319 322 320 CALL obs_grid_search( inowin, zlam, zphi, iobsi, iobsj, iproc, 'T' ) 323 ! Assume anything other than velocity is on T grid 324 IF ( TRIM( inpfiles(jj)%cname(1) ) == 'UVEL' ) THEN 325 CALL obs_grid_search( inowin, zlam, zphi, iobsi(:,1), iobsj(:,1), & 326 & iproc(:,1), 'U' ) 327 CALL obs_grid_search( inowin, zlam, zphi, iobsi(:,2), iobsj(:,2), & 328 & iproc(:,2), 'V' ) 329 ELSE 330 CALL obs_grid_search( inowin, zlam, zphi, iobsi(:,1), iobsj(:,1), & 331 & iproc(:,1), 'T' ) 332 IF ( kvars > 1 ) THEN 333 DO jvar = 2, kvars 334 iobsi(:,jvar) = iobsi(:,1) 335 iobsj(:,jvar) = iobsj(:,1) 336 iproc(:,jvar) = iproc(:,1) 337 END DO 338 ENDIF 339 ENDIF 321 340 322 341 inowin = 0 … … 325 344 & ( inpfiles(jj)%ptim(ji) <= djulend(jj) ) ) THEN 326 345 inowin = inowin + 1 327 inpfiles(jj)%iproc(ji,1) = iproc(inowin) 328 inpfiles(jj)%iobsi(ji,1) = iobsi(inowin) 329 inpfiles(jj)%iobsj(ji,1) = iobsj(inowin) 346 DO jvar = 1, kvars 347 inpfiles(jj)%iproc(ji,jvar) = iproc(inowin,jvar) 348 inpfiles(jj)%iobsi(ji,jvar) = iobsi(inowin,jvar) 349 inpfiles(jj)%iobsj(ji,jvar) = iobsj(inowin,jvar) 350 END DO 330 351 ENDIF 331 352 END DO … … 448 469 449 470 ! Coordinate search parameters 450 surfdata%mi (iobs) = inpfiles(jj)%iobsi(ji,1) 451 surfdata%mj (iobs) = inpfiles(jj)%iobsj(ji,1) 452 471 DO jvar = 1, kvars 472 surfdata%mi(iobs,jvar) = inpfiles(jj)%iobsi(ji,jvar) 473 surfdata%mj(iobs,jvar) = inpfiles(jj)%iobsj(ji,jvar) 474 END DO 475 453 476 ! WMO number 454 477 surfdata%cwmo(iobs) = inpfiles(jj)%cdwmo(ji) … … 480 503 481 504 ! Observed value 482 surfdata%robs(iobs,1) = inpfiles(jj)%pob(1,ji,1) 505 DO jvar = 1, kvars 506 surfdata%robs(iobs,jvar) = inpfiles(jj)%pob(1,ji,jvar) 507 END DO 483 508 IF ( TRIM(surfdata%cvars(1)) == 'FBD' ) THEN 484 509 surfdata%rext(iobs,1) = inpfiles(jj)%pob(1,ji,1) … … 488 513 ! Model and MDT is set to fbrmdi unless read from file 489 514 IF ( ldmod ) THEN 490 surfdata%rmod(iobs,1) = inpfiles(jj)%padd(1,ji,1,1) 515 DO jvar = 1, kvars 516 surfdata%rmod(iobs,jvar) = inpfiles(jj)%padd(1,ji,1,jvar) 517 END DO 491 518 IF ( TRIM(surfdata%cvars(1)) == 'SLA' ) THEN 492 519 surfdata%rext(iobs,1) = inpfiles(jj)%padd(1,ji,2,1) … … 494 521 ENDIF 495 522 ELSE 496 surfdata%rmod(iobs,1) = fbrmdi 523 DO jvar = 1, kvars 524 surfdata%rmod(iobs,jvar) = fbrmdi 525 END DO 497 526 IF ( TRIM(surfdata%cvars(1)) == 'SLA' ) surfdata%rext(iobs,:) = fbrmdi 498 527 ENDIF 499 528 500 529 ! Initialise climatology if set 501 IF ( surfdata%lclim ) surfdata%rclm(iobs,1) = fbrmdi 530 IF ( surfdata%lclim ) THEN 531 DO jvar = 1, kvars 532 surfdata%rclm(iobs,jvar) = fbrmdi 533 END DO 534 ENDIF 502 535 503 536 ! STD (obs error standard deviation) read from file and passed through obs operator … … 521 554 !----------------------------------------------------------------------- 522 555 IF (lwp) THEN 523 556 DO jvar = 1, surfdata%nvar 557 IF ( jvar == 1 ) THEN 558 cout1=TRIM(surfdata%cvars(1)) 559 ELSE 560 WRITE(cout1,'(A,A1,A)') TRIM(cout1), '/', TRIM(surfdata%cvars(jvar)) 561 ENDIF 562 END DO 563 524 564 WRITE(numout,*) 525 WRITE(numout,'(1X,A)')TRIM( surfdata%cvars(1))//' data'565 WRITE(numout,'(1X,A)')TRIM( cout1 )//' data' 526 566 WRITE(numout,'(1X,A)')'--------------' 527 567 DO jj = 1,8 … … 533 573 & '---------------------------------------------------------------' 534 574 WRITE(numout,'(1X,A,I8)') & 535 & 'Total data for variable '//TRIM( surfdata%cvars(1))// &575 & 'Total data for variable '//TRIM( cout1 )// & 536 576 & ' = ', iobsmpp 537 577 WRITE(numout,'(1X,A)') & -
branches/UKMO/dev_r5518_obs_oper_update_sit_SMOS/NEMOGCM/NEMO/OPA_SRC/OBS/obs_rot_vel.F90
r7992 r14165 6 6 7 7 !!---------------------------------------------------------------------- 8 !! obs_rotvel : Rotate velocity data into N-S,E-W directorions 8 !! obs_rotvel_pro : Rotate profile velocity data into N-S,E-W directions 9 !! obs_rotvel_surf : Rotate surface velocity data into N-S,E-W directions 9 10 !!---------------------------------------------------------------------- 10 11 !! * Modules used … … 17 18 USE obs_utils ! For error handling 18 19 USE obs_profiles_def ! Profile definitions 20 USE obs_surf_def ! Surface definitions 19 21 USE obs_inter_h2d ! Horizontal interpolation 20 22 USE obs_inter_sup ! MPP support routines for interpolation … … 27 29 PRIVATE 28 30 29 PUBLIC obs_rotvel ! Rotate the observations 30 31 PUBLIC obs_rotvel_pro, & ! Rotate the profile velocity observations 32 & obs_rotvel_surf ! Rotate the surface velocity observations 33 31 34 !!---------------------------------------------------------------------- 32 35 !! NEMO/OPA 3.3 , NEMO Consortium (2010) … … 37 40 CONTAINS 38 41 39 SUBROUTINE obs_rotvel ( profdata, k2dint, pu, pv )42 SUBROUTINE obs_rotvel_pro( profdata, k2dint, pu, pv ) 40 43 !!--------------------------------------------------------------------- 41 44 !! … … 228 231 CALL wrk_dealloc(jpi,jpj,zsingu,zcosgu,zsingv,zcosgv) 229 232 230 END SUBROUTINE obs_rotvel 233 END SUBROUTINE obs_rotvel_pro 234 235 SUBROUTINE obs_rotvel_surf( surfdata, k2dint, pu, pv ) 236 !!--------------------------------------------------------------------- 237 !! 238 !! *** ROUTINE obs_rotvel_surf *** 239 !! 240 !! ** Purpose : Rotate surface velocity data into N-S,E-W directorions 241 !! 242 !! ** Method : Interpolation of geo2ocean coefficients on U,V grid 243 !! to observation point followed by a similar computations 244 !! as in geo2ocean. 245 !! 246 !! ** Action : Review if there is a better way to do this. 247 !! 248 !! References : 249 !! 250 !! History : 251 !! ! : 2009-02 (K. Mogensen) : New routine 252 !!---------------------------------------------------------------------- 253 !! * Modules used 254 !! * Arguments 255 TYPE(obs_surf), INTENT(INOUT) :: surfdata ! Surface data to be read 256 INTEGER, INTENT(IN) :: k2dint ! Horizontal interpolation methed 257 REAL(wp), DIMENSION(*) :: & 258 & pu, & 259 & pv 260 !! * Local declarations 261 REAL(wp), DIMENSION(2,2,1) :: zweig 262 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: & 263 & zmasku, & 264 & zmaskv, & 265 & zcoslu, & 266 & zsinlu, & 267 & zcoslv, & 268 & zsinlv, & 269 & zglamu, & 270 & zgphiu, & 271 & zglamv, & 272 & zgphiv 273 REAL(wp), DIMENSION(1) :: & 274 & zsinu, & 275 & zcosu, & 276 & zsinv, & 277 & zcosv 278 REAL(wp) :: zsin 279 REAL(wp) :: zcos 280 REAL(wp), DIMENSION(1) :: zobsmask 281 REAL(wp), POINTER, DIMENSION(:,:) :: zsingu,zcosgu,zsingv,zcosgv 282 INTEGER, DIMENSION(:,:,:), ALLOCATABLE :: & 283 & igrdiu, & 284 & igrdju, & 285 & igrdiv, & 286 & igrdjv 287 INTEGER :: ji 288 INTEGER :: jk 289 290 CALL wrk_alloc(jpi,jpj,zsingu,zcosgu,zsingv,zcosgv) 291 292 !----------------------------------------------------------------------- 293 ! Allocate data for message parsing and interpolation 294 !----------------------------------------------------------------------- 295 296 ALLOCATE( & 297 & igrdiu(2,2,surfdata%nsurf), & 298 & igrdju(2,2,surfdata%nsurf), & 299 & zglamu(2,2,surfdata%nsurf), & 300 & zgphiu(2,2,surfdata%nsurf), & 301 & zmasku(2,2,surfdata%nsurf), & 302 & zcoslu(2,2,surfdata%nsurf), & 303 & zsinlu(2,2,surfdata%nsurf), & 304 & igrdiv(2,2,surfdata%nsurf), & 305 & igrdjv(2,2,surfdata%nsurf), & 306 & zglamv(2,2,surfdata%nsurf), & 307 & zgphiv(2,2,surfdata%nsurf), & 308 & zmaskv(2,2,surfdata%nsurf), & 309 & zcoslv(2,2,surfdata%nsurf), & 310 & zsinlv(2,2,surfdata%nsurf) & 311 & ) 312 313 !----------------------------------------------------------------------- 314 ! Receive the angles on the U and V grids. 315 !----------------------------------------------------------------------- 316 317 CALL obs_rot( zsingu, zcosgu, zsingv, zcosgv ) 318 319 DO ji = 1, surfdata%nsurf 320 igrdiu(1,1,ji) = surfdata%mi(ji,1)-1 321 igrdju(1,1,ji) = surfdata%mj(ji,1)-1 322 igrdiu(1,2,ji) = surfdata%mi(ji,1)-1 323 igrdju(1,2,ji) = surfdata%mj(ji,1) 324 igrdiu(2,1,ji) = surfdata%mi(ji,1) 325 igrdju(2,1,ji) = surfdata%mj(ji,1)-1 326 igrdiu(2,2,ji) = surfdata%mi(ji,1) 327 igrdju(2,2,ji) = surfdata%mj(ji,1) 328 igrdiv(1,1,ji) = surfdata%mi(ji,2)-1 329 igrdjv(1,1,ji) = surfdata%mj(ji,2)-1 330 igrdiv(1,2,ji) = surfdata%mi(ji,2)-1 331 igrdjv(1,2,ji) = surfdata%mj(ji,2) 332 igrdiv(2,1,ji) = surfdata%mi(ji,2) 333 igrdjv(2,1,ji) = surfdata%mj(ji,2)-1 334 igrdiv(2,2,ji) = surfdata%mi(ji,2) 335 igrdjv(2,2,ji) = surfdata%mj(ji,2) 336 END DO 337 338 CALL obs_int_comm_2d( 2, 2, surfdata%nsurf, jpi, jpj, igrdiu, igrdju, & 339 & glamu, zglamu ) 340 CALL obs_int_comm_2d( 2, 2, surfdata%nsurf, jpi, jpj, igrdiu, igrdju, & 341 & gphiu, zgphiu ) 342 CALL obs_int_comm_2d( 2, 2, surfdata%nsurf, jpi, jpj, igrdiu, igrdju, & 343 & umask(:,:,1), zmasku ) 344 CALL obs_int_comm_2d( 2, 2, surfdata%nsurf, jpi, jpj, igrdiu, igrdju, & 345 & zsingu, zsinlu ) 346 CALL obs_int_comm_2d( 2, 2, surfdata%nsurf, jpi, jpj, igrdiu, igrdju, & 347 & zcosgu, zcoslu ) 348 CALL obs_int_comm_2d( 2, 2, surfdata%nsurf, jpi, jpj, igrdiv, igrdjv, & 349 & glamv, zglamv ) 350 CALL obs_int_comm_2d( 2, 2, surfdata%nsurf, jpi, jpj, igrdiv, igrdjv, & 351 & gphiv, zgphiv ) 352 CALL obs_int_comm_2d( 2, 2, surfdata%nsurf, jpi, jpj, igrdiv, igrdjv, & 353 & vmask(:,:,1), zmaskv ) 354 CALL obs_int_comm_2d( 2, 2, surfdata%nsurf, jpi, jpj, igrdiv, igrdjv, & 355 & zsingv, zsinlv ) 356 CALL obs_int_comm_2d( 2, 2, surfdata%nsurf, jpi, jpj, igrdiv, igrdjv, & 357 & zcosgv, zcoslv ) 358 359 DO ji = 1, surfdata%nsurf 360 361 CALL obs_int_h2d_init( 1, 1, k2dint, & 362 & surfdata%rlam(ji), surfdata%rphi(ji), & 363 & zglamu(:,:,ji), zgphiu(:,:,ji), & 364 & zmasku(:,:,ji), zweig, zobsmask ) 365 366 CALL obs_int_h2d( 1, 1, zweig, zsinlu(:,:,ji), zsinu ) 367 368 CALL obs_int_h2d( 1, 1, zweig, zcoslu(:,:,ji), zcosu ) 369 370 CALL obs_int_h2d_init( 1, 1, k2dint, & 371 & surfdata%rlam(ji), surfdata%rphi(ji), & 372 & zglamv(:,:,ji), zgphiv(:,:,ji), & 373 & zmaskv(:,:,ji), zweig, zobsmask ) 374 375 CALL obs_int_h2d( 1, 1, zweig, zsinlv(:,:,ji), zsinv ) 376 377 CALL obs_int_h2d( 1, 1, zweig, zcoslv(:,:,ji), zcosv ) 378 379 ! Assume that the angle at observation point is the 380 ! mean of u and v cosines/sines 381 382 zcos = 0.5_wp * ( zcosu(1) + zcosv(1) ) 383 zsin = 0.5_wp * ( zsinu(1) + zsinv(1) ) 384 385 IF ( ( surfdata%rmod(ji,1) /= fbrmdi ) .AND. & 386 & ( surfdata%rmod(ji,2) /= fbrmdi ) ) THEN 387 pu(ji) = surfdata%rmod(ji,1) * zcos - & 388 & surfdata%rmod(ji,2) * zsin 389 pv(ji) = surfdata%rmod(ji,2) * zcos + & 390 & surfdata%rmod(ji,1) * zsin 391 ELSE 392 pu(ji) = fbrmdi 393 pv(ji) = fbrmdi 394 ENDIF 395 396 397 END DO 398 399 DEALLOCATE( & 400 & igrdiu, & 401 & igrdju, & 402 & zglamu, & 403 & zgphiu, & 404 & zmasku, & 405 & zcoslu, & 406 & zsinlu, & 407 & igrdiv, & 408 & igrdjv, & 409 & zglamv, & 410 & zgphiv, & 411 & zmaskv, & 412 & zcoslv, & 413 & zsinlv & 414 & ) 415 416 CALL wrk_dealloc(jpi,jpj,zsingu,zcosgu,zsingv,zcosgv) 417 418 END SUBROUTINE obs_rotvel_surf 231 419 232 420 END MODULE obs_rot_vel -
branches/UKMO/dev_r5518_obs_oper_update_sit_SMOS/NEMOGCM/NEMO/OPA_SRC/OBS/obs_surf_def.F90
r11468 r14165 57 57 58 58 INTEGER, POINTER, DIMENSION(:) :: & 59 & mi, & !: i-th grid coord. for interpolating to surface observation60 & mj, & !: j-th grid coord. for interpolating to surface observation61 59 & mt, & !: time record number for gridded data 62 60 & nsidx,& !: Surface observation number … … 71 69 & ntyp !: Type of surface observation product 72 70 71 INTEGER, POINTER, DIMENSION(:,:) :: & 72 & mi, & !: i-th grid coord. for interpolating to surface observation 73 & mj !: j-th grid coord. for interpolating to surface observation 74 75 73 76 CHARACTER(len=8), POINTER, DIMENSION(:) :: & 74 77 & cvars !: Variable names … … 92 95 & rext !: Extra fields interpolated to observation points 93 96 94 REAL(KIND=wp), POINTER, DIMENSION(:,: ) :: &97 REAL(KIND=wp), POINTER, DIMENSION(:,:,:) :: & 95 98 & vdmean !: Time averaged of model field 96 99 … … 176 179 177 180 ALLOCATE( & 178 & surf%mi(ksurf), &179 & surf%mj(ksurf), &180 181 & surf%mt(ksurf), & 181 182 & surf%nsidx(ksurf), & … … 201 202 202 203 ALLOCATE( & 204 & surf%mi(ksurf,kvar), & 205 & surf%mj(ksurf,kvar), & 203 206 & surf%robs(ksurf,kvar), & 204 & surf%rmod(ksurf,kvar) &207 & surf%rmod(ksurf,kvar) & 205 208 & ) 206 209 … … 230 233 231 234 ALLOCATE( & 232 & surf%vdmean(kpi,kpj ) &235 & surf%vdmean(kpi,kpj,kvar) & 233 236 & ) 234 237 … … 405 408 insurf = insurf + 1 406 409 407 newsurf%mi(insurf ) = surf%mi(ji)408 newsurf%mj(insurf ) = surf%mj(ji)410 newsurf%mi(insurf,:) = surf%mi(ji,:) 411 newsurf%mj(insurf,:) = surf%mj(ji,:) 409 412 newsurf%mt(insurf) = surf%mt(ji) 410 413 newsurf%nsidx(insurf) = surf%nsidx(ji) … … 496 499 jj=surf%nsind(ji) 497 500 498 oldsurf%mi(jj ) = surf%mi(ji)499 oldsurf%mj(jj ) = surf%mj(ji)501 oldsurf%mi(jj,:) = surf%mi(ji,:) 502 oldsurf%mj(jj,:) = surf%mj(ji,:) 500 503 oldsurf%mt(jj) = surf%mt(ji) 501 504 oldsurf%nsidx(jj) = surf%nsidx(ji) -
branches/UKMO/dev_r5518_obs_oper_update_sit_SMOS/NEMOGCM/NEMO/OPA_SRC/OBS/obs_write.F90
r12820 r14165 349 349 fbdata%ivlqc(ik,jo,jvar) = IBSET(profdata%var(jvar)%nvqc(jk),2) 350 350 fbdata%ivlqcf(1,ik,jo,jvar) = profdata%var(jvar)%nvqcf(1,jk) 351 fbdata%ivlqcf(2,ik,jo,jvar) = IAND(profdata%var(jvar)%nvqc(jk),b'0000 0000 11111111')351 fbdata%ivlqcf(2,ik,jo,jvar) = IAND(profdata%var(jvar)%nvqc(jk),b'0000000011111111') 352 352 ELSE 353 353 fbdata%ivlqc(ik,jo,jvar) = profdata%var(jvar)%nvqc(jk) … … 437 437 CHARACTER(LEN=10) :: clfiletype 438 438 CHARACTER(LEN=12), PARAMETER :: cpname = 'obs_wri_surf' 439 CHARACTER(LEN=ilenlong) :: cllongname ! Long name of variable440 CHARACTER(LEN=ilenunit) :: clunits ! Units of variable441 CHARACTER(LEN=ilengrid) :: clgrid ! Grid of variable439 CHARACTER(LEN=ilenlong), DIMENSION(surfdata%nvar) :: cllongname ! Long name of variable 440 CHARACTER(LEN=ilenunit), DIMENSION(surfdata%nvar) :: clunits ! Units of variable 441 CHARACTER(LEN=ilengrid), DIMENSION(surfdata%nvar) :: clgrid ! Grid of variable 442 442 INTEGER :: jo 443 443 INTEGER :: ja 444 444 INTEGER :: je 445 INTEGER :: jv 445 446 INTEGER :: iadd 446 447 INTEGER :: iext … … 522 523 CASE('SST') 523 524 524 clfiletype = 'sstfb'525 cllongname = 'Sea surface temperature'526 clunits = 'Degree centigrade'527 clgrid = 'T'525 clfiletype = 'sstfb' 526 cllongname(1) = 'Sea surface temperature' 527 clunits(1) = 'Degree centigrade' 528 clgrid(1) = 'T' 528 529 529 530 CASE('ICECONC') 530 531 531 clfiletype = 'sicfb'532 cllongname = 'Sea ice concentration'533 clunits = 'Fraction'534 clgrid = 'T'532 clfiletype = 'sicfb' 533 cllongname(1) = 'Sea ice concentration' 534 clunits(1) = 'Fraction' 535 clgrid(1) = 'T' 535 536 536 537 CASE('SIT') 537 538 538 clfiletype = 'sitfb'539 cllongname = 'Sea ice thickness'540 clunits = 'm'541 clgrid = 'T'539 clfiletype = 'sitfb' 540 cllongname(1) = 'Sea ice thickness' 541 clunits(1) = 'm' 542 clgrid(1) = 'T' 542 543 543 544 CASE('FBD') … … 552 553 CASE('SSS') 553 554 554 clfiletype = 'sssfb' 555 cllongname = 'Sea surface salinity' 556 clunits = 'psu' 557 clgrid = 'T' 555 clfiletype = 'sssfb' 556 cllongname(1) = 'Sea surface salinity' 557 clunits(1) = 'psu' 558 clgrid(1) = 'T' 559 560 CASE('UVEL') 561 562 clfiletype = 'ssvfb' 563 cllongname(1) = 'Eastward sea surface velocity' 564 clunits(1) = 'm s-1' 565 clgrid(1) = 'U' 566 cllongname(2) = 'Northward sea surface velocity' 567 clunits(2) = 'm s-1' 568 clgrid(2) = 'V' 558 569 559 570 CASE('SLCHLTOT') 560 571 561 clfiletype = 'slchltotfb'562 cllongname = 'Surface total log10(chlorophyll)'563 clunits = 'log10(mg/m3)'564 clgrid = 'T'572 clfiletype = 'slchltotfb' 573 cllongname(1) = 'Surface total log10(chlorophyll)' 574 clunits(1) = 'log10(mg/m3)' 575 clgrid(1) = 'T' 565 576 566 577 CASE('SLCHLDIA') 567 578 568 clfiletype = 'slchldiafb'569 cllongname = 'Surface diatom log10(chlorophyll)'570 clunits = 'log10(mg/m3)'571 clgrid = 'T'579 clfiletype = 'slchldiafb' 580 cllongname(1) = 'Surface diatom log10(chlorophyll)' 581 clunits(1) = 'log10(mg/m3)' 582 clgrid(1) = 'T' 572 583 573 584 CASE('SLCHLNON') 574 585 575 clfiletype = 'slchlnonfb'576 cllongname = 'Surface non-diatom log10(chlorophyll)'577 clunits = 'log10(mg/m3)'578 clgrid = 'T'586 clfiletype = 'slchlnonfb' 587 cllongname(1) = 'Surface non-diatom log10(chlorophyll)' 588 clunits(1) = 'log10(mg/m3)' 589 clgrid(1) = 'T' 579 590 580 591 CASE('SLCHLDIN') 581 592 582 clfiletype = 'slchldinfb'583 cllongname = 'Surface dinoflagellate log10(chlorophyll)'584 clunits = 'log10(mg/m3)'585 clgrid = 'T'593 clfiletype = 'slchldinfb' 594 cllongname(1) = 'Surface dinoflagellate log10(chlorophyll)' 595 clunits(1) = 'log10(mg/m3)' 596 clgrid(1) = 'T' 586 597 587 598 CASE('SLCHLMIC') 588 599 589 clfiletype = 'slchlmicfb'590 cllongname = 'Surface microphytoplankton log10(chlorophyll)'591 clunits = 'log10(mg/m3)'592 clgrid = 'T'600 clfiletype = 'slchlmicfb' 601 cllongname(1) = 'Surface microphytoplankton log10(chlorophyll)' 602 clunits(1) = 'log10(mg/m3)' 603 clgrid(1) = 'T' 593 604 594 605 CASE('SLCHLNAN') 595 606 596 clfiletype = 'slchlnanfb'597 cllongname = 'Surface nanophytoplankton log10(chlorophyll)'598 clunits = 'log10(mg/m3)'599 clgrid = 'T'607 clfiletype = 'slchlnanfb' 608 cllongname(1) = 'Surface nanophytoplankton log10(chlorophyll)' 609 clunits(1) = 'log10(mg/m3)' 610 clgrid(1) = 'T' 600 611 601 612 CASE('SLCHLPIC') 602 613 603 clfiletype = 'slchlpicfb'604 cllongname = 'Surface picophytoplankton log10(chlorophyll)'605 clunits = 'log10(mg/m3)'606 clgrid = 'T'614 clfiletype = 'slchlpicfb' 615 cllongname(1) = 'Surface picophytoplankton log10(chlorophyll)' 616 clunits(1) = 'log10(mg/m3)' 617 clgrid(1) = 'T' 607 618 608 619 CASE('SCHLTOT') … … 615 626 CASE('SLPHYTOT') 616 627 617 clfiletype = 'slphytotfb'618 cllongname = 'Surface total log10(phytoplankton carbon)'619 clunits = 'log10(mmolC/m3)'620 clgrid = 'T'628 clfiletype = 'slphytotfb' 629 cllongname(1) = 'Surface total log10(phytoplankton carbon)' 630 clunits(1) = 'log10(mmolC/m3)' 631 clgrid(1) = 'T' 621 632 622 633 CASE('SLPHYDIA') 623 634 624 clfiletype = 'slphydiafb'625 cllongname = 'Surface diatom log10(phytoplankton carbon)'626 clunits = 'log10(mmolC/m3)'627 clgrid = 'T'635 clfiletype = 'slphydiafb' 636 cllongname(1) = 'Surface diatom log10(phytoplankton carbon)' 637 clunits(1) = 'log10(mmolC/m3)' 638 clgrid(1) = 'T' 628 639 629 640 CASE('SLPHYNON') 630 641 631 clfiletype = 'slphynonfb'632 cllongname = 'Surface non-diatom log10(phytoplankton carbon)'633 clunits = 'log10(mmolC/m3)'634 clgrid = 'T'642 clfiletype = 'slphynonfb' 643 cllongname(1) = 'Surface non-diatom log10(phytoplankton carbon)' 644 clunits(1) = 'log10(mmolC/m3)' 645 clgrid(1) = 'T' 635 646 636 647 CASE('SSPM') 637 648 638 clfiletype = 'sspmfb'639 cllongname = 'Surface suspended particulate matter'640 clunits = 'g/m3'641 clgrid = 'T'649 clfiletype = 'sspmfb' 650 cllongname(1) = 'Surface suspended particulate matter' 651 clunits(1) = 'g/m3' 652 clgrid(1) = 'T' 642 653 643 654 CASE('SKD490') 644 655 645 clfiletype = 'skd490fb'646 cllongname = 'Surface attenuation coefficient of downwelling radiation at 490 nm'647 clunits = 'm-1'648 clgrid = 'T'656 clfiletype = 'skd490fb' 657 cllongname(1) = 'Surface attenuation coefficient of downwelling radiation at 490 nm' 658 clunits(1) = 'm-1' 659 clgrid(1) = 'T' 649 660 650 661 CASE('SFCO2') 651 662 652 clfiletype = 'sfco2fb'653 cllongname = 'Surface fugacity of carbon dioxide'654 clunits = 'uatm'655 clgrid = 'T'663 clfiletype = 'sfco2fb' 664 cllongname(1) = 'Surface fugacity of carbon dioxide' 665 clunits(1) = 'uatm' 666 clgrid(1) = 'T' 656 667 657 668 CASE('SPCO2') 658 669 659 clfiletype = 'spco2fb'660 cllongname = 'Surface partial pressure of carbon dioxide'661 clunits = 'uatm'662 clgrid = 'T'670 clfiletype = 'spco2fb' 671 cllongname(1) = 'Surface partial pressure of carbon dioxide' 672 clunits(1) = 'uatm' 673 clgrid(1) = 'T' 663 674 664 675 CASE DEFAULT … … 673 684 IF ( TRIM(surfdata%cvars(1)) /= 'SLA' ) THEN 674 685 675 CALL alloc_obfbdata( fbdata, 1, surfdata%nsurf, 1, &686 CALL alloc_obfbdata( fbdata, surfdata%nvar, surfdata%nsurf, 1, & 676 687 & 1 + iadd_std + iadd_clm + iadd, iext, .TRUE. ) 677 688 678 fbdata%cname(1) = surfdata%cvars(1) 679 fbdata%coblong(1) = cllongname 680 fbdata%cobunit(1) = clunits 689 DO jv = 1, surfdata%nvar 690 fbdata%cname(jv) = surfdata%cvars(jv) 691 fbdata%coblong(jv) = cllongname(jv) 692 fbdata%cobunit(jv) = clunits(jv) 693 END DO 681 694 DO je = 1, iext 682 695 fbdata%cextname(je) = pext%cdname(je) … … 684 697 fbdata%cextunit(je) = pext%cdunit(je,1) 685 698 END DO 686 IF ( TRIM(surfdata%cvars(1)) == 'ICECONC' ) THEN 687 fbdata%caddlong(1,1) = 'Model interpolated ICE' 688 ELSE 689 fbdata%caddlong(1,1) = 'Model interpolated ' // TRIM(surfdata%cvars(1)) 690 ENDIF 691 fbdata%caddunit(1,1) = clunits 692 fbdata%cgrid(1) = clgrid 699 DO jv = 1, surfdata%nvar 700 IF ( TRIM(surfdata%cvars(1)) == 'ICECONC' ) THEN 701 fbdata%caddlong(1,jv) = 'Model interpolated ICE' 702 ELSE 703 fbdata%caddlong(1,jv) = 'Model interpolated ' // TRIM(surfdata%cvars(jv)) 704 ENDIF 705 fbdata%caddunit(1,jv) = clunits(jv) 706 fbdata%cgrid(jv) = clgrid(jv) 707 END DO 693 708 DO ja = 1, iadd 694 709 fbdata%caddname(1+iadd_mdt+iadd_std+iadd_clm+ja) = padd%cdname(ja) 695 fbdata%caddlong(1+iadd_mdt+iadd_std+iadd_clm+ja,1) = padd%cdlong(ja,1) 696 fbdata%caddunit(1+iadd_mdt+iadd_std+iadd_clm+ja,1) = padd%cdunit(ja,1) 697 END DO 698 699 ENDIF 700 710 DO jv = 1, surfdata%nvar 711 fbdata%caddlong(1+iadd_mdt+iadd_std+iadd_clm+ja,jv) = padd%cdlong(ja,jv) 712 fbdata%caddunit(1+iadd_mdt+iadd_std+iadd_clm+ja,jv) = padd%cdunit(ja,jv) 713 END DO 714 END DO 715 ENDIF 716 701 717 fbdata%caddname(1) = 'Hx' 702 718 IF ( indx_std /= -1 ) THEN 703 719 fbdata%caddname(1+iadd_mdt+iadd_std) = surfdata%cext(indx_std) 704 fbdata%caddlong(1+iadd_mdt+iadd_std,1) = 'Obs error standard deviation' 705 fbdata%caddunit(1+iadd_mdt+iadd_std,1) = fbdata%cobunit(1) 720 DO jv = 1, surfdata%nvar 721 fbdata%caddlong(1+iadd_mdt+iadd_std,1) = 'Obs error standard deviation' 722 fbdata%caddunit(1+iadd_mdt+iadd_std,1) = fbdata%cobunit(1) 723 END DO 706 724 ENDIF 707 725 708 726 IF ( surfdata%lclim ) THEN 709 727 fbdata%caddname(1+iadd_mdt+iadd_std+iadd_clm) = 'CLM' 710 fbdata%caddlong(1+iadd_mdt+iadd_std+iadd_clm,1) = 'Climatology' 711 fbdata%caddunit(1+iadd_mdt+iadd_std+iadd_clm,1) = fbdata%cobunit(1) 712 ENDIF 713 714 728 DO jv = 1, surfdata%nvar 729 fbdata%caddlong(1+iadd_mdt+iadd_std+iadd_clm,jv) = 'Climatology' 730 fbdata%caddunit(1+iadd_mdt+iadd_std+iadd_clm,jv) = fbdata%cobunit(1) 731 END DO 732 ENDIF 733 715 734 WRITE(clfname, FMT="(A,'_fdbk_',I4.4,'.nc')") TRIM(clfiletype), nproc 716 735 … … 719 738 WRITE(numout,*)'obs_wri_surf :' 720 739 WRITE(numout,*)'~~~~~~~~~~~~~' 721 WRITE(numout,*)'Writing '//TRIM( surfdata%cvars(1))//' feedback file : ',TRIM(clfname)740 WRITE(numout,*)'Writing '//TRIM(clfiletype)//' feedback file : ',TRIM(clfname) 722 741 ENDIF 723 742 … … 733 752 fbdata%ioqc(jo) = 4 734 753 fbdata%ioqcf(1,jo) = 0 735 fbdata%ioqcf(2,jo) = IAND(surfdata%nqc(jo),b'0000 0000 11111111')754 fbdata%ioqcf(2,jo) = IAND(surfdata%nqc(jo),b'0000000011111111') 736 755 ELSE 737 756 fbdata%ioqc(jo) = surfdata%nqc(jo) … … 745 764 fbdata%kindex(jo) = surfdata%nsfil(jo) 746 765 IF (ln_grid_global) THEN 747 fbdata%iobsi(jo,1) = surfdata%mi(jo) 748 fbdata%iobsj(jo,1) = surfdata%mj(jo) 766 DO jv = 1, surfdata%nvar 767 fbdata%iobsi(jo,jv) = surfdata%mi(jo,jv) 768 fbdata%iobsj(jo,jv) = surfdata%mj(jo,jv) 769 END DO 749 770 ELSE 750 fbdata%iobsi(jo,1) = mig(surfdata%mi(jo)) 751 fbdata%iobsj(jo,1) = mjg(surfdata%mj(jo)) 771 DO jv = 1, surfdata%nvar 772 fbdata%iobsi(jo,jv) = mig(surfdata%mi(jo,jv)) 773 fbdata%iobsj(jo,jv) = mjg(surfdata%mj(jo,jv)) 774 END DO 752 775 ENDIF 753 776 CALL greg2jul( 0, & … … 759 782 & fbdata%ptim(jo), & 760 783 & krefdate = 19500101 ) 761 762 fbdata%pob(1,jo,1) = surfdata%robs(jo,1) 784 763 785 fbdata%pdep(1,jo) = 0.0 764 786 fbdata%idqc(1,jo) = 0 765 787 fbdata%idqcf(:,1,jo) = 0 766 IF ( surfdata%nqc(jo) > 255 ) THEN 767 fbdata%ivqc(jo,1) = 4 768 fbdata%ivlqc(1,jo,1) = 4 769 fbdata%ivlqcf(1,1,jo,1) = 0 770 fbdata%ivlqcf(2,1,jo,1) = IAND(surfdata%nqc(jo),b'0000 0000 1111 1111') 771 ELSE 772 fbdata%ivqc(jo,1) = surfdata%nqc(jo) 773 fbdata%ivlqc(1,jo,1) = surfdata%nqc(jo) 774 fbdata%ivlqcf(:,1,jo,1) = 0 775 ENDIF 776 fbdata%iobsk(1,jo,1) = 0 777 778 ! Additional variables. 779 ! Hx is always the first additional variable 780 fbdata%padd(1,jo,1,1) = surfdata%rmod(jo,1) 781 ! MDT is output as an additional variable if SLA obs type 782 IF ( TRIM(surfdata%cvars(1)) == 'SLA' ) THEN 783 fbdata%padd(1,jo,1+iadd_mdt,1) = surfdata%rext(jo,1) 784 ENDIF 785 ! STD is output as an additional variable if available 786 IF ( indx_std /= -1 ) THEN 787 fbdata%padd(1,jo,1+iadd_mdt+iadd_std,1) = surfdata%rext(jo,indx_std) 788 ENDIF 789 ! CLM is output as an additional variable if available 790 IF ( surfdata%lclim ) THEN 791 fbdata%padd(1,jo,1+iadd_mdt+iadd_std+iadd_clm,1) = surfdata%rclm(jo,1) 792 ENDIF 793 ! Then other additional variables are output 794 DO ja = 1, iadd 795 fbdata%padd(1,jo,1+iadd_mdt+iadd_std+iadd_clm+ja,1) = & 796 & surfdata%rext(jo,padd%ipoint(ja)) 797 END DO 788 DO jv = 1, surfdata%nvar 789 fbdata%pob(1,jo,jv) = surfdata%robs(jo,jv) 798 790 791 IF ( surfdata%nqc(jo) > 255 ) THEN 792 fbdata%ivqc(jo,jv) = 4 793 fbdata%ivlqc(1,jo,jv) = 4 794 fbdata%ivlqcf(1,1,jo,jv) = 0 795 fbdata%ivlqcf(2,1,jo,jv) = IAND(surfdata%nqc(jo),b'0000000011111111') 796 ELSE 797 fbdata%ivqc(jo,jv) = surfdata%nqc(jo) 798 fbdata%ivlqc(1,jo,jv) = surfdata%nqc(jo) 799 fbdata%ivlqcf(:,1,jo,jv) = 0 800 ENDIF 801 fbdata%iobsk(1,jo,jv) = 0 802 803 804 ! Additional variables. 805 ! Hx is always the first additional variable 806 fbdata%padd(1,jo,1,jv) = surfdata%rmod(jo,jv) 807 ! MDT is output as an additional variable if SLA obs type 808 IF ( TRIM(surfdata%cvars(jv)) == 'SLA' ) THEN 809 fbdata%padd(1,jo,1+iadd_mdt,jv) = surfdata%rext(jo,1) 810 ENDIF 811 ! STD is output as an additional variable if available 812 IF ( indx_std /= -1 ) THEN 813 fbdata%padd(1,jo,1+iadd_mdt+iadd_std,jv) = surfdata%rext(jo,indx_std) 814 ENDIF 815 ! CLM is output as an additional variable if available 816 IF ( surfdata%lclim ) THEN 817 fbdata%padd(1,jo,1+iadd_mdt+iadd_std+iadd_clm,jv) = surfdata%rclm(jo,jv) 818 ENDIF 819 ! Then other additional variables are output 820 DO ja = 1, iadd 821 fbdata%padd(1,jo,1+iadd_mdt+iadd_std+iadd_clm+ja,jv) = & 822 & surfdata%rext(jo,padd%ipoint(ja)) 823 END DO 824 END DO 799 825 ! Extra variables 800 826 IF ( TRIM(surfdata%cvars(1)) == 'SLA' ) fbdata%pext(1,jo,1) = surfdata%rext(jo,2)
Note: See TracChangeset
for help on using the changeset viewer.
