New URL for NEMO forge! http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.

obs_prep.F90 in branches/UKMO/dev_r5518_obs_oper_update/NEMOGCM/NEMO/OPA_SRC/OBS – NEMO

Context Navigation

source: branches/UKMO/dev_r5518_obs_oper_update/NEMOGCM/NEMO/OPA_SRC/OBS/obs_prep.F90 @ 13393

Last change on this file since 13393 was 13393, checked in by mattmartin, 4 years ago
Include surface velocity observation operator. See ticket https://code.metoffice.gov.uk/trac/utils/ticket/376.
File size: 63.0 KB

Rev	Line
[2128]	1	MODULE obs_prep
	2	!!=====================================================================
	3	!! * MODULE obs_prep *
	4	!! Observation diagnostics: Prepare observation arrays: screening,
	5	!! sorting, coordinate search
	6	!!=====================================================================
	7
	8	!!---------------------------------------------------------------------
[7992]	9	!! obs_pre_prof : First level check and screening of profile observations
	10	!! obs_pre_surf : First level check and screening of surface observations
	11	!! obs_scr : Basic screening of the observations
	12	!! obs_coo_tim : Compute number of time steps to the observation time
	13	!! obs_sor : Sort the observation arrays
[2128]	14	!!---------------------------------------------------------------------
	15	!! * Modules used
	16	USE par_kind, ONLY : & ! Precision variables
	17	& wp
[12506]	18	USE dom_oce ! ocean space and time domain
[2128]	19	USE in_out_manager ! I/O manager
	20	USE obs_profiles_def ! Definitions for storage arrays for profiles
	21	USE obs_surf_def ! Definitions for storage arrays for surface data
	22	USE obs_mpp, ONLY : & ! MPP support routines for observation diagnostics
	23	& obs_mpp_sum_integer, &
	24	& obs_mpp_sum_integers
	25	USE obs_inter_sup ! Interpolation support
	26	USE obs_oper ! Observation operators
[7992]	27	#if defined key_bdy
	28	USE bdy_oce, ONLY : & ! Boundary information
	29	idx_bdy, nb_bdy
	30	#endif
[2715]	31	USE lib_mpp, ONLY : &
	32	& ctl_warn, ctl_stop
[2128]	33
	34	IMPLICIT NONE
	35
	36	!! * Routine accessibility
	37	PRIVATE
	38
	39	PUBLIC &
[7992]	40	& obs_pre_prof, & ! First level check and screening of profile obs
	41	& obs_pre_surf, & ! First level check and screening of surface obs
	42	& calc_month_len ! Calculate the number of days in the months of a year
[2128]	43
[2287]	44	!!----------------------------------------------------------------------
	45	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
	46	!! $Id$
	47	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
	48	!!----------------------------------------------------------------------
	49
[7992]	50	!! * Substitutions
	51	# include "domzgr_substitute.h90"
	52
[2128]	53	CONTAINS
	54
[13393]	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 )
[2128]	60	!!----------------------------------------------------------------------
[13393]	61	!! * ROUTINE obs_pre_surf *
[2128]	62	!!
[7992]	63	!! ** Purpose : First level check and screening of surface observations
[2128]	64	!!
[7992]	65	!! ** Method : First level check and screening of surface observations
[2128]	66	!!
	67	!! ** Action :
	68	!!
	69	!! References :
	70	!!
	71	!! History :
	72	!! ! 2007-03 (A. Weaver, K. Mogensen) Original
	73	!! ! 2007-06 (K. Mogensen et al) Reject obs. near land.
[7992]	74	!! ! 2015-02 (M. Martin) Combined routine for surface types.
[2128]	75	!!----------------------------------------------------------------------
	76	!! * Modules used
	77	USE domstp ! Domain: set the time-step
	78	USE par_oce ! Ocean parameters
	79	USE dom_oce, ONLY : & ! Geographical information
	80	& glamt, &
	81	& gphit, &
	82	& tmask, &
	83	& nproc
	84	!! * Arguments
[7992]	85	TYPE(obs_surf), INTENT(INOUT) :: surfdata ! Full set of surface data
	86	TYPE(obs_surf), INTENT(INOUT) :: surfdataqc ! Subset of surface data not failing screening
[13393]	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
[7992]	93	LOGICAL, INTENT(IN) :: ld_nea ! Switch for rejecting observation near land
	94	LOGICAL, INTENT(IN) :: ld_bound_reject ! Switch for rejecting obs near the boundary
[12610]	95	LOGICAL, INTENT(IN) :: ld_seaicetypes ! Switch to indicate sea ice data
	96	INTEGER, INTENT(IN), OPTIONAL :: kqc_cutoff ! cut off for QC value
[2128]	97	!! * Local declarations
[7992]	98	INTEGER :: iqc_cutoff = 255 ! cut off for QC value
[2128]	99	INTEGER :: iyea0 ! Initial date
	100	INTEGER :: imon0 ! - (year, month, day, hour, minute)
	101	INTEGER :: iday0
	102	INTEGER :: ihou0
	103	INTEGER :: imin0
	104	INTEGER :: icycle ! Current assimilation cycle
[13393]	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
[2128]	120	LOGICAL, DIMENSION(:), ALLOCATABLE :: &
	121	& llvalid ! SLA data selection
[13393]	122	INTEGER :: jobs ! Obs. loop counter
	123	INTEGER :: jvar ! Variable loop counter
[2128]	124	INTEGER :: jstp ! Time loop variable
	125	INTEGER :: inrc ! Time index variable
[13393]	126	CHARACTER(LEN=256) :: cout1 ! Diagnostic output line
	127
[7992]	128	IF(lwp) WRITE(numout,*)'obs_pre_surf : Preparing the surface observations...'
	129	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
	130
[2128]	131	! Initial date initialization (year, month, day, hour, minute)
	132	iyea0 = ndate0 / 10000
	133	imon0 = ( ndate0 - iyea0 * 10000 ) / 100
	134	iday0 = ndate0 - iyea0 * 10000 - imon0 * 100
	135	ihou0 = 0
	136	imin0 = 0
	137
	138	icycle = no ! Assimilation cycle
	139
	140	! Diagnotics counters for various failures.
	141
	142	iotdobs = 0
	143	igrdobs = 0
[13393]	144	iosdsobs(:) = 0
	145	ilansobs(:) = 0
	146	inlasobs(:) = 0
	147	ibdysobs(:) = 0
[2128]	148
[7992]	149	! Set QC cutoff to optional value if provided
	150	IF ( PRESENT(kqc_cutoff) ) iqc_cutoff=kqc_cutoff
	151
[2128]	152	! -----------------------------------------------------------------------
[7992]	153	! Find time coordinate for surface data
[2128]	154	! -----------------------------------------------------------------------
	155
[12610]	156	IF ( ld_seaicetypes ) THEN
	157	CALL obs_coo_tim( icycle, &
	158	& iyea0, imon0, iday0, ihou0, imin0, &
	159	& surfdata%nsurf, surfdata%nyea, surfdata%nmon, &
	160	& surfdata%nday, surfdata%nhou, surfdata%nmin, &
	161	& surfdata%nqc, surfdata%mstp, iotdobs, &
	162	& ld_seaicetypes = ld_seaicetypes )
	163	ELSE
	164	CALL obs_coo_tim( icycle, &
	165	& iyea0, imon0, iday0, ihou0, imin0, &
	166	& surfdata%nsurf, surfdata%nyea, surfdata%nmon, &
	167	& surfdata%nday, surfdata%nhou, surfdata%nmin, &
	168	& surfdata%nqc, surfdata%mstp, iotdobs )
	169	ENDIF
[2128]	170
	171	CALL obs_mpp_sum_integer( iotdobs, iotdobsmpp )
	172
	173	! -----------------------------------------------------------------------
[7992]	174	! Check for surface data failing the grid search
[2128]	175	! -----------------------------------------------------------------------
[13393]	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
[2128]	181	CALL obs_mpp_sum_integer( igrdobs, igrdobsmpp )
	182
	183	! -----------------------------------------------------------------------
	184	! Check for land points.
	185	! -----------------------------------------------------------------------
	186
[13393]	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
[2128]	204	! -----------------------------------------------------------------------
[7992]	205	! Copy useful data from the surfdata data structure to
	206	! the surfdataqc data structure
[2128]	207	! -----------------------------------------------------------------------
	208
	209	! Allocate the selection arrays
	210
[7992]	211	ALLOCATE( llvalid(surfdata%nsurf) )
[2128]	212
[7992]	213	! We want all data which has qc flags <= iqc_cutoff
[2128]	214
[7992]	215	llvalid(:) = ( surfdata%nqc(:) <= iqc_cutoff )
[2128]	216
	217	! The actual copying
	218
[7992]	219	CALL obs_surf_compress( surfdata, surfdataqc, .TRUE., numout, &
[2128]	220	& lvalid=llvalid )
	221
	222	! Dellocate the selection arrays
	223	DEALLOCATE( llvalid )
	224
	225	! -----------------------------------------------------------------------
	226	! Print information about what observations are left after qc
	227	! -----------------------------------------------------------------------
	228
	229	! Update the total observation counter array
	230
	231	IF(lwp) THEN
[13393]	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
[2128]	241	WRITE(numout,*)
[13393]	242	WRITE(numout,*) ' '//TRIM(cout1)//' data outside time domain = ', &
[2128]	243	& iotdobsmpp
[13393]	244	WRITE(numout,*) ' Remaining '//TRIM(cout1)//' data that failed grid search = ', &
[2128]	245	& igrdobsmpp
	246
[13393]	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
	264
[2128]	265	WRITE(numout,*)
	266	WRITE(numout,*) ' Number of observations per time step :'
	267	WRITE(numout,*)
[13393]	268	WRITE(numout,'(10X,A,10X,A)')'Time step',TRIM(cout1)
[7992]	269	WRITE(numout,'(10X,A,5X,A)')'---------','-----------------'
	270	CALL FLUSH(numout)
[2128]	271	ENDIF
	272
[7992]	273	DO jobs = 1, surfdataqc%nsurf
	274	inrc = surfdataqc%mstp(jobs) + 2 - nit000
	275	surfdataqc%nsstp(inrc) = surfdataqc%nsstp(inrc) + 1
[2128]	276	END DO
	277
[7992]	278	CALL obs_mpp_sum_integers( surfdataqc%nsstp, surfdataqc%nsstpmpp, &
[2128]	279	& nitend - nit000 + 2 )
	280
	281	IF ( lwp ) THEN
	282	DO jstp = nit000 - 1, nitend
	283	inrc = jstp - nit000 + 2
[7992]	284	WRITE(numout,1999) jstp, surfdataqc%nsstpmpp(inrc)
	285	CALL FLUSH(numout)
[2128]	286	END DO
	287	ENDIF
	288
	289	1999 FORMAT(10X,I9,5X,I17)
	290
[7992]	291	END SUBROUTINE obs_pre_surf
[2128]	292
	293
[9306]	294	SUBROUTINE obs_pre_prof( profdata, prodatqc, ld_var, &
[7992]	295	& kpi, kpj, kpk, &
[9306]	296	& zmask, pglam, pgphi, &
[7992]	297	& ld_nea, ld_bound_reject, kdailyavtypes, kqc_cutoff )
[2128]	298
[7992]	299	!!----------------------------------------------------------------------
	300	!! * ROUTINE obs_pre_prof *
[2128]	301	!!
[7992]	302	!! ** Purpose : First level check and screening of profiles
[2128]	303	!!
[7992]	304	!! ** Method : First level check and screening of profiles
[2128]	305	!!
	306	!! History :
	307	!! ! 2007-06 (K. Mogensen) original : T and S profile data
	308	!! ! 2008-09 (M. Valdivieso) : TAO velocity data
[7992]	309	!! ! 2009-01 (K. Mogensen) : New feedback stricture
	310	!! ! 2015-02 (M. Martin) : Combined profile routine.
[2128]	311	!!
	312	!!----------------------------------------------------------------------
	313	!! * Modules used
	314	USE domstp ! Domain: set the time-step
	315	USE par_oce ! Ocean parameters
	316	USE dom_oce, ONLY : & ! Geographical information
[4292]	317	& gdept_1d, &
[2128]	318	& nproc
[7992]	319
[2128]	320	!! * Arguments
	321	TYPE(obs_prof), INTENT(INOUT) :: profdata ! Full set of profile data
	322	TYPE(obs_prof), INTENT(INOUT) :: prodatqc ! Subset of profile data not failing screening
[9306]	323	LOGICAL, DIMENSION(profdata%nvar), INTENT(IN) :: &
	324	& ld_var ! Observed variables switches
[7992]	325	LOGICAL, INTENT(IN) :: ld_nea ! Switch for rejecting observation near land
	326	LOGICAL, INTENT(IN) :: ld_bound_reject ! Switch for rejecting observations near the boundary
	327	INTEGER, INTENT(IN) :: kpi, kpj, kpk ! Local domain sizes
	328	INTEGER, DIMENSION(imaxavtypes), OPTIONAL :: &
	329	& kdailyavtypes ! Types for daily averages
[9306]	330	REAL(wp), INTENT(IN), DIMENSION(kpi,kpj,kpk,profdata%nvar) :: &
	331	& zmask
	332	REAL(wp), INTENT(IN), DIMENSION(kpi,kpj,profdata%nvar) :: &
	333	& pglam, &
	334	& pgphi
[7992]	335	INTEGER, INTENT(IN), OPTIONAL :: kqc_cutoff ! cut off for QC value
	336
[2128]	337	!! * Local declarations
[7992]	338	INTEGER :: iqc_cutoff = 255 ! cut off for QC value
[2128]	339	INTEGER :: iyea0 ! Initial date
	340	INTEGER :: imon0 ! - (year, month, day, hour, minute)
	341	INTEGER :: iday0
	342	INTEGER :: ihou0
	343	INTEGER :: imin0
	344	INTEGER :: icycle ! Current assimilation cycle
[9306]	345	! Counters for observations that are
	346	INTEGER :: iotdobs ! - outside time domain
	347	INTEGER, DIMENSION(profdata%nvar) :: iosdvobs ! - outside space domain
	348	INTEGER, DIMENSION(profdata%nvar) :: ilanvobs ! - within a model land cell
	349	INTEGER, DIMENSION(profdata%nvar) :: inlavobs ! - close to land
	350	INTEGER, DIMENSION(profdata%nvar) :: ibdyvobs ! - boundary
	351	INTEGER :: igrdobs ! - fail the grid search
	352	INTEGER :: iuvchku ! - reject UVEL if VVEL rejected
	353	INTEGER :: iuvchkv ! - reject VVEL if UVEL rejected
	354	! Global counters for observations that are
	355	INTEGER :: iotdobsmpp ! - outside time domain
	356	INTEGER, DIMENSION(profdata%nvar) :: iosdvobsmpp ! - outside space domain
	357	INTEGER, DIMENSION(profdata%nvar) :: ilanvobsmpp ! - within a model land cell
	358	INTEGER, DIMENSION(profdata%nvar) :: inlavobsmpp ! - close to land
	359	INTEGER, DIMENSION(profdata%nvar) :: ibdyvobsmpp ! - boundary
	360	INTEGER :: igrdobsmpp ! - fail the grid search
	361	INTEGER :: iuvchkumpp ! - reject UVEL if VVEL rejected
	362	INTEGER :: iuvchkvmpp ! - reject VVEL if UVEL rejected
[2128]	363	TYPE(obs_prof_valid) :: llvalid ! Profile selection
	364	TYPE(obs_prof_valid), DIMENSION(profdata%nvar) :: &
[9306]	365	& llvvalid ! vars selection
[2128]	366	INTEGER :: jvar ! Variable loop variable
	367	INTEGER :: jobs ! Obs. loop variable
	368	INTEGER :: jstp ! Time loop variable
	369	INTEGER :: inrc ! Time index variable
[9306]	370	CHARACTER(LEN=256) :: cout1 ! Diagnostic output line
	371	CHARACTER(LEN=256) :: cout2 ! Diagnostic output line
[2128]	372
[7992]	373	IF(lwp) WRITE(numout,*)'obs_pre_prof: Preparing the profile data...'
	374	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
[2128]	375
	376	! Initial date initialization (year, month, day, hour, minute)
	377	iyea0 = ndate0 / 10000
	378	imon0 = ( ndate0 - iyea0 * 10000 ) / 100
	379	iday0 = ndate0 - iyea0 * 10000 - imon0 * 100
	380	ihou0 = 0
	381	imin0 = 0
	382
	383	icycle = no ! Assimilation cycle
	384
[9306]	385	! Diagnostics counters for various failures.
[2128]	386
[9306]	387	iotdobs = 0
	388	igrdobs = 0
	389	iosdvobs(:) = 0
	390	ilanvobs(:) = 0
	391	inlavobs(:) = 0
	392	ibdyvobs(:) = 0
	393	iuvchku = 0
	394	iuvchkv = 0
[2128]	395
[7992]	396
	397	! Set QC cutoff to optional value if provided
	398	IF ( PRESENT(kqc_cutoff) ) iqc_cutoff=kqc_cutoff
	399
[2128]	400	! -----------------------------------------------------------------------
	401	! Find time coordinate for profiles
	402	! -----------------------------------------------------------------------
	403
[7992]	404	IF ( PRESENT(kdailyavtypes) ) THEN
	405	CALL obs_coo_tim_prof( icycle, &
	406	& iyea0, imon0, iday0, ihou0, imin0, &
	407	& profdata%nprof, profdata%nyea, profdata%nmon, &
	408	& profdata%nday, profdata%nhou, profdata%nmin, &
	409	& profdata%ntyp, profdata%nqc, profdata%mstp, &
	410	& iotdobs, kdailyavtypes = kdailyavtypes, &
	411	& kqc_cutoff = iqc_cutoff )
	412	ELSE
	413	CALL obs_coo_tim_prof( icycle, &
	414	& iyea0, imon0, iday0, ihou0, imin0, &
	415	& profdata%nprof, profdata%nyea, profdata%nmon, &
	416	& profdata%nday, profdata%nhou, profdata%nmin, &
	417	& profdata%ntyp, profdata%nqc, profdata%mstp, &
	418	& iotdobs, kqc_cutoff = iqc_cutoff )
	419	ENDIF
	420
[2128]	421	CALL obs_mpp_sum_integer( iotdobs, iotdobsmpp )
	422
	423	! -----------------------------------------------------------------------
	424	! Check for profiles failing the grid search
	425	! -----------------------------------------------------------------------
	426
[9306]	427	DO jvar = 1, profdata%nvar
	428	CALL obs_coo_grd( profdata%nprof, profdata%mi(:,jvar), profdata%mj(:,jvar), &
	429	& profdata%nqc, igrdobs )
	430	END DO
[2128]	431
	432	CALL obs_mpp_sum_integer( igrdobs, igrdobsmpp )
	433
	434	! -----------------------------------------------------------------------
[7992]	435	! Reject all observations for profiles with nqc > iqc_cutoff
[2128]	436	! -----------------------------------------------------------------------
	437
[7992]	438	CALL obs_pro_rej( profdata, kqc_cutoff = iqc_cutoff )
[2128]	439
	440	! -----------------------------------------------------------------------
	441	! Check for land points. This includes points below the model
	442	! bathymetry so this is done for every point in the profile
	443	! -----------------------------------------------------------------------
	444
[9306]	445	DO jvar = 1, profdata%nvar
	446	CALL obs_coo_spc_3d( profdata%nprof, profdata%nvprot(jvar), &
	447	& profdata%npvsta(:,jvar), profdata%npvend(:,jvar), &
	448	& jpi, jpj, &
	449	& jpk, &
	450	& profdata%mi, profdata%mj, &
	451	& profdata%var(jvar)%mvk, &
	452	& profdata%rlam, profdata%rphi, &
	453	& profdata%var(jvar)%vdep, &
	454	& pglam(:,:,jvar), pgphi(:,:,jvar), &
	455	& gdept_1d, zmask(:,:,:,jvar), &
	456	& profdata%nqc, profdata%var(jvar)%nvqc, &
	457	& iosdvobs(jvar), ilanvobs(jvar), &
	458	& inlavobs(jvar), ld_nea, &
	459	& ibdyvobs(jvar), ld_bound_reject, &
	460	& iqc_cutoff )
[2128]	461
[9306]	462	CALL obs_mpp_sum_integer( iosdvobs(jvar), iosdvobsmpp(jvar) )
	463	CALL obs_mpp_sum_integer( ilanvobs(jvar), ilanvobsmpp(jvar) )
	464	CALL obs_mpp_sum_integer( inlavobs(jvar), inlavobsmpp(jvar) )
	465	CALL obs_mpp_sum_integer( ibdyvobs(jvar), ibdyvobsmpp(jvar) )
	466	END DO
[2128]	467
	468	! -----------------------------------------------------------------------
	469	! Reject u if v is rejected and vice versa
	470	! -----------------------------------------------------------------------
	471
[7992]	472	IF ( TRIM(profdata%cvars(1)) == 'UVEL' ) THEN
[13393]	473	CALL obs_uv_rej_pro( profdata, iuvchku, iuvchkv, iqc_cutoff )
[7992]	474	CALL obs_mpp_sum_integer( iuvchku, iuvchkumpp )
	475	CALL obs_mpp_sum_integer( iuvchkv, iuvchkvmpp )
	476	ENDIF
[2128]	477
	478	! -----------------------------------------------------------------------
	479	! Copy useful data from the profdata data structure to
	480	! the prodatqc data structure
	481	! -----------------------------------------------------------------------
	482
	483	! Allocate the selection arrays
	484
	485	ALLOCATE( llvalid%luse(profdata%nprof) )
	486	DO jvar = 1,profdata%nvar
	487	ALLOCATE( llvvalid(jvar)%luse(profdata%nvprot(jvar)) )
	488	END DO
	489
[7992]	490	! We want all data which has qc flags <= iqc_cutoff
[2128]	491
[7992]	492	llvalid%luse(:) = ( profdata%nqc(:) <= iqc_cutoff )
[2128]	493	DO jvar = 1,profdata%nvar
[7992]	494	llvvalid(jvar)%luse(:) = ( profdata%var(jvar)%nvqc(:) <= iqc_cutoff )
[2128]	495	END DO
	496
	497	! The actual copying
	498
	499	CALL obs_prof_compress( profdata, prodatqc, .TRUE., numout, &
	500	& lvalid=llvalid, lvvalid=llvvalid )
	501
	502	! Dellocate the selection arrays
	503	DEALLOCATE( llvalid%luse )
	504	DO jvar = 1,profdata%nvar
	505	DEALLOCATE( llvvalid(jvar)%luse )
	506	END DO
	507
	508	! -----------------------------------------------------------------------
	509	! Print information about what observations are left after qc
	510	! -----------------------------------------------------------------------
	511
	512	! Update the total observation counter array
	513
	514	IF(lwp) THEN
[7992]	515
[2128]	516	WRITE(numout,*)
[9306]	517	WRITE(numout,*) ' Profiles outside time domain = ', &
[2128]	518	& iotdobsmpp
[9306]	519	WRITE(numout,*) ' Remaining profiles that failed grid search = ', &
[2128]	520	& igrdobsmpp
[9306]	521	DO jvar = 1, profdata%nvar
	522	WRITE(numout,*) ' Remaining '//prodatqc%cvars(jvar)//' data outside space domain = ', &
	523	& iosdvobsmpp(jvar)
	524	WRITE(numout,*) ' Remaining '//prodatqc%cvars(jvar)//' data at land points = ', &
	525	& ilanvobsmpp(jvar)
	526	IF (ld_nea) THEN
	527	WRITE(numout,*) ' Remaining '//prodatqc%cvars(jvar)//' data near land points (removed) = ',&
	528	& inlavobsmpp(jvar)
	529	ELSE
	530	WRITE(numout,*) ' Remaining '//prodatqc%cvars(jvar)//' data near land points (kept) = ',&
	531	& inlavobsmpp(jvar)
	532	ENDIF
	533	IF ( TRIM(profdata%cvars(jvar)) == 'UVEL' ) THEN
	534	WRITE(numout,*) ' U observation rejected since V rejected = ', &
	535	& iuvchku
	536	ELSE IF ( TRIM(profdata%cvars(jvar)) == 'VVEL' ) THEN
	537	WRITE(numout,*) ' V observation rejected since U rejected = ', &
	538	& iuvchkv
	539	ENDIF
	540	WRITE(numout,*) ' Remaining '//prodatqc%cvars(jvar)//' data near open boundary (removed) = ',&
	541	& ibdyvobsmpp(jvar)
	542	WRITE(numout,*) ' '//prodatqc%cvars(jvar)//' data accepted = ', &
	543	& prodatqc%nvprotmpp(jvar)
	544	END DO
[2128]	545
	546	WRITE(numout,*)
	547	WRITE(numout,*) ' Number of observations per time step :'
	548	WRITE(numout,*)
[9306]	549	WRITE(cout1,'(10X,A9,5X,A8)') 'Time step', 'Profiles'
	550	WRITE(cout2,'(10X,A9,5X,A8)') '---------', '--------'
	551	DO jvar = 1, prodatqc%nvar
	552	WRITE(cout1,'(A,5X,A11)') TRIM(cout1), TRIM(prodatqc%cvars(jvar))
	553	WRITE(cout2,'(A,5X,A11)') TRIM(cout2), '-----------'
	554	END DO
	555	WRITE(numout,*) cout1
	556	WRITE(numout,*) cout2
[2128]	557	ENDIF
	558
	559	DO jobs = 1, prodatqc%nprof
	560	inrc = prodatqc%mstp(jobs) + 2 - nit000
	561	prodatqc%npstp(inrc) = prodatqc%npstp(inrc) + 1
	562	DO jvar = 1, prodatqc%nvar
	563	IF ( prodatqc%npvend(jobs,jvar) > 0 ) THEN
	564	prodatqc%nvstp(inrc,jvar) = prodatqc%nvstp(inrc,jvar) + &
	565	& ( prodatqc%npvend(jobs,jvar) - &
	566	& prodatqc%npvsta(jobs,jvar) + 1 )
	567	ENDIF
	568	END DO
	569	END DO
	570
	571
	572	CALL obs_mpp_sum_integers( prodatqc%npstp, prodatqc%npstpmpp, &
	573	& nitend - nit000 + 2 )
	574	DO jvar = 1, prodatqc%nvar
	575	CALL obs_mpp_sum_integers( prodatqc%nvstp(:,jvar), &
	576	& prodatqc%nvstpmpp(:,jvar), &
	577	& nitend - nit000 + 2 )
	578	END DO
	579
	580	IF ( lwp ) THEN
	581	DO jstp = nit000 - 1, nitend
	582	inrc = jstp - nit000 + 2
[9306]	583	WRITE(cout1,'(10X,I9,5X,I8)') jstp, prodatqc%npstpmpp(inrc)
	584	DO jvar = 1, prodatqc%nvar
	585	WRITE(cout1,'(A,5X,I11)') TRIM(cout1), prodatqc%nvstpmpp(inrc,jvar)
	586	END DO
	587	WRITE(numout,*) cout1
[2128]	588	END DO
	589	ENDIF
	590
[7992]	591	END SUBROUTINE obs_pre_prof
[2128]	592
	593	SUBROUTINE obs_coo_tim( kcycle, &
	594	& kyea0, kmon0, kday0, khou0, kmin0, &
	595	& kobsno, &
	596	& kobsyea, kobsmon, kobsday, kobshou, kobsmin, &
[12610]	597	& kobsqc, kobsstp, kotdobs, ld_seaicetypes )
[2128]	598	!!----------------------------------------------------------------------
	599	!! * ROUTINE obs_coo_tim *
	600	!!
	601	!! ** Purpose : Compute the number of time steps to the observation time.
	602	!!
	603	!! ** Method : For time coordinates ( yea_obs, mon_obs, day_obs,
	604	!! hou_obs, min_obs ), this routine locates the time step
	605	!! that is closest to this time.
	606	!!
	607	!! ** Action :
	608	!!
	609	!! References :
	610	!!
	611	!! History :
	612	!! ! 1997-07 (A. Weaver) Original
	613	!! ! 2006-08 (A. Weaver) NEMOVAR migration
	614	!! ! 2006-10 (A. Weaver) Cleanup
	615	!! ! 2007-01 (K. Mogensen) Rewritten with loop
	616	!! ! 2010-05 (D. Lea) Fix in leap year calculation for NEMO vn3.2
	617	!!----------------------------------------------------------------------
	618	!! * Modules used
	619	USE dom_oce, ONLY : & ! Geographical information
	620	& rdt
	621	USE phycst, ONLY : & ! Physical constants
	622	& rday, &
	623	& rmmss, &
	624	& rhhmm
	625	!! * Arguments
	626	INTEGER, INTENT(IN) :: kcycle ! Current cycle
	627	INTEGER, INTENT(IN) :: kyea0 ! Initial date coordinates
	628	INTEGER, INTENT(IN) :: kmon0
	629	INTEGER, INTENT(IN) :: kday0
	630	INTEGER, INTENT(IN) :: khou0
	631	INTEGER, INTENT(IN) :: kmin0
	632	INTEGER, INTENT(IN) :: kobsno ! Number of observations
	633	INTEGER, INTENT(INOUT) :: kotdobs ! Number of observations failing time check
	634	INTEGER, DIMENSION(kobsno), INTENT(IN ) :: &
	635	& kobsyea, & ! Observation time coordinates
	636	& kobsmon, &
	637	& kobsday, &
	638	& kobshou, &
	639	& kobsmin
	640	INTEGER, DIMENSION(kobsno), INTENT(INOUT) :: &
	641	& kobsqc ! Quality control flag
	642	INTEGER, DIMENSION(kobsno), INTENT(OUT) :: &
	643	& kobsstp ! Number of time steps up to the
	644	! observation time
[12610]	645	LOGICAL, OPTIONAL, INTENT(IN) :: ld_seaicetypes
[2128]	646
	647	!! * Local declarations
	648	INTEGER :: jyea
	649	INTEGER :: jmon
	650	INTEGER :: jday
	651	INTEGER :: jobs
	652	INTEGER :: iyeastr
	653	INTEGER :: iyeaend
	654	INTEGER :: imonstr
	655	INTEGER :: imonend
	656	INTEGER :: idaystr
	657	INTEGER :: idayend
	658	INTEGER :: iskip
	659	INTEGER :: idaystp
[12610]	660	INTEGER :: icecount
[2128]	661	REAL(KIND=wp) :: zminstp
	662	REAL(KIND=wp) :: zhoustp
	663	REAL(KIND=wp) :: zobsstp
	664	INTEGER, DIMENSION(12) :: imonth_len !: length in days of the months of the current year
	665
	666	!-----------------------------------------------------------------------
	667	! Initialization
	668	!-----------------------------------------------------------------------
	669
	670	! Intialize the number of time steps per day
	671	idaystp = NINT( rday / rdt )
	672
	673	!---------------------------------------------------------------------
	674	! Locate the model time coordinates for interpolation
	675	!---------------------------------------------------------------------
	676
	677	DO jobs = 1, kobsno
	678
	679	! Initialize the time step counter
	680	kobsstp(jobs) = nit000 - 1
	681
	682	! Flag if observation date is less than the initial date
	683
	684	IF ( ( kobsyea(jobs) < kyea0 ) &
	685	& .OR. ( ( kobsyea(jobs) == kyea0 ) &
	686	& .AND. ( kobsmon(jobs) < kmon0 ) ) &
	687	& .OR. ( ( kobsyea(jobs) == kyea0 ) &
	688	& .AND. ( kobsmon(jobs) == kmon0 ) &
	689	& .AND. ( kobsday(jobs) < kday0 ) ) &
	690	& .OR. ( ( kobsyea(jobs) == kyea0 ) &
	691	& .AND. ( kobsmon(jobs) == kmon0 ) &
	692	& .AND. ( kobsday(jobs) == kday0 ) &
	693	& .AND. ( kobshou(jobs) < khou0 ) ) &
	694	& .OR. ( ( kobsyea(jobs) == kyea0 ) &
	695	& .AND. ( kobsmon(jobs) == kmon0 ) &
	696	& .AND. ( kobsday(jobs) == kday0 ) &
	697	& .AND. ( kobshou(jobs) == khou0 ) &
	698	& .AND. ( kobsmin(jobs) <= kmin0 ) ) ) THEN
	699	kobsstp(jobs) = -1
[7992]	700	kobsqc(jobs) = IBSET(kobsqc(jobs),13)
[2128]	701	kotdobs = kotdobs + 1
	702	CYCLE
	703	ENDIF
	704
	705	! Compute the number of time steps to the observation day
	706	iyeastr = kyea0
	707	iyeaend = kobsyea(jobs)
	708
	709	!---------------------------------------------------------------------
	710	! Year loop
	711	!---------------------------------------------------------------------
	712	DO jyea = iyeastr, iyeaend
	713
	714	CALL calc_month_len( jyea, imonth_len )
	715
	716	imonstr = 1
	717	IF ( jyea == kyea0 ) imonstr = kmon0
	718	imonend = 12
	719	IF ( jyea == kobsyea(jobs) ) imonend = kobsmon(jobs)
	720
	721	! Month loop
	722	DO jmon = imonstr, imonend
	723
	724	idaystr = 1
	725	IF ( ( jmon == kmon0 ) &
	726	& .AND. ( jyea == kyea0 ) ) idaystr = kday0
	727	idayend = imonth_len(jmon)
	728	IF ( ( jmon == kobsmon(jobs) ) &
	729	& .AND. ( jyea == kobsyea(jobs) ) ) idayend = kobsday(jobs) - 1
	730
	731	! Day loop
	732	DO jday = idaystr, idayend
	733	kobsstp(jobs) = kobsstp(jobs) + idaystp
	734	END DO
	735
	736	END DO
	737
	738	END DO
	739
	740	! Add in the number of time steps to the observation minute
	741	zminstp = rmmss / rdt
	742	zhoustp = rhhmm * zminstp
	743
	744	zobsstp = REAL( kobsmin(jobs) - kmin0, KIND=wp ) * zminstp &
	745	& + REAL( kobshou(jobs) - khou0, KIND=wp ) * zhoustp
	746	kobsstp(jobs) = kobsstp(jobs) + NINT( zobsstp )
	747
	748	! Flag if observation step outside the time window
	749	IF ( ( kobsstp(jobs) < ( nit000 - 1 ) ) &
	750	& .OR.( kobsstp(jobs) > nitend ) ) THEN
[7992]	751	kobsqc(jobs) = IBSET(kobsqc(jobs),13)
[2128]	752	kotdobs = kotdobs + 1
	753	CYCLE
	754	ENDIF
	755
[12610]	756	! Flag sea ice observations falling on initial timestep
	757	IF ( PRESENT(ld_seaicetypes) ) THEN
	758
	759	IF ( ( kobsstp(jobs) == (nit000 - 1) ) ) THEN
	760	IF (lwp) WRITE(numout,*)( 'Sea-ice not initialised on zeroth '// &
	761	& 'time-step but SIT observation valid then, flagging '// &
	762	'in time check subroutine obs_coo_tim.' )
	763	kobsqc(jobs) = IBSET(kobsqc(jobs),13)
	764	kotdobs = kotdobs + 1
	765	CYCLE
	766	ENDIF
	767	ENDIF
	768
[2128]	769	END DO
	770
	771	END SUBROUTINE obs_coo_tim
	772
	773	SUBROUTINE calc_month_len( iyear, imonth_len )
	774	!!----------------------------------------------------------------------
	775	!! * ROUTINE calc_month_len *
	776	!!
	777	!! ** Purpose : Compute the number of days in a months given a year.
	778	!!
	779	!! ** Method :
	780	!!
	781	!! ** Action :
	782	!!
	783	!! History :
	784	!! ! 10-05 (D. Lea) New routine based on day_init
	785	!!----------------------------------------------------------------------
	786
	787	INTEGER, DIMENSION(12) :: imonth_len !: length in days of the months of the current year
	788	INTEGER :: iyear !: year
	789
	790	! length of the month of the current year (from nleapy, read in namelist)
	791	IF ( nleapy < 2 ) THEN
	792	imonth_len(:) = (/ 31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31 /)
	793	IF ( nleapy == 1 ) THEN ! we are using calendar with leap years
	794	IF ( MOD(iyear, 4) == 0 .AND. ( MOD(iyear, 400) == 0 .OR. MOD(iyear, 100) /= 0 ) ) THEN
	795	imonth_len(2) = 29
	796	ENDIF
	797	ENDIF
	798	ELSE
	799	imonth_len(:) = nleapy ! all months with nleapy days per year
	800	ENDIF
	801
	802	END SUBROUTINE
	803
	804	SUBROUTINE obs_coo_tim_prof( kcycle, &
	805	& kyea0, kmon0, kday0, khou0, kmin0, &
	806	& kobsno, &
	807	& kobsyea, kobsmon, kobsday, kobshou, kobsmin, &
	808	& ktyp, kobsqc, kobsstp, kotdobs, kdailyavtypes, &
[7992]	809	& kqc_cutoff )
[2128]	810	!!----------------------------------------------------------------------
	811	!! * ROUTINE obs_coo_tim *
	812	!!
	813	!! ** Purpose : Compute the number of time steps to the observation time.
	814	!!
	815	!! ** Method : For time coordinates ( yea_obs, mon_obs, day_obs,
	816	!! hou_obs, min_obs ), this routine locates the time step
	817	!! that is closest to this time.
	818	!!
	819	!! ** Action :
	820	!!
	821	!! References :
	822	!!
	823	!! History :
	824	!! ! 1997-07 (A. Weaver) Original
	825	!! ! 2006-08 (A. Weaver) NEMOVAR migration
	826	!! ! 2006-10 (A. Weaver) Cleanup
	827	!! ! 2007-01 (K. Mogensen) Rewritten with loop
	828	!!----------------------------------------------------------------------
	829	!! * Modules used
	830	!! * Arguments
	831	INTEGER, INTENT(IN) :: kcycle ! Current cycle
	832	INTEGER, INTENT(IN) :: kyea0 ! Initial date coordinates
	833	INTEGER, INTENT(IN) :: kmon0
	834	INTEGER, INTENT(IN) :: kday0
	835	INTEGER, INTENT(IN) :: khou0
	836	INTEGER, INTENT(IN) :: kmin0
	837	INTEGER, INTENT(IN) :: kobsno ! Number of observations
	838	INTEGER, INTENT(INOUT) :: kotdobs ! Number of observations failing time check
	839	INTEGER, DIMENSION(kobsno), INTENT(IN ) :: &
	840	& kobsyea, & ! Observation time coordinates
	841	& kobsmon, &
	842	& kobsday, &
	843	& kobshou, &
	844	& kobsmin, &
	845	& ktyp ! Observation type.
	846	INTEGER, DIMENSION(kobsno), INTENT(INOUT) :: &
	847	& kobsqc ! Quality control flag
	848	INTEGER, DIMENSION(kobsno), INTENT(OUT) :: &
	849	& kobsstp ! Number of time steps up to the
	850	! observation time
	851	INTEGER, DIMENSION(imaxavtypes), OPTIONAL :: &
	852	& kdailyavtypes ! Types for daily averages
[7992]	853	INTEGER, OPTIONAL, INTENT(IN) :: kqc_cutoff ! QC cutoff value
	854
[2128]	855	!! * Local declarations
	856	INTEGER :: jobs
[7992]	857	INTEGER :: iqc_cutoff=255
[2128]	858
	859	!-----------------------------------------------------------------------
	860	! Call standard obs_coo_tim
	861	!-----------------------------------------------------------------------
	862
	863	CALL obs_coo_tim( kcycle, &
	864	& kyea0, kmon0, kday0, khou0, kmin0, &
	865	& kobsno, &
	866	& kobsyea, kobsmon, kobsday, kobshou, kobsmin, &
	867	& kobsqc, kobsstp, kotdobs )
	868
	869	!------------------------------------------------------------------------
	870	! Always reject daily averaged data (e.g. MRB data (820)) at initial time
	871	!------------------------------------------------------------------------
	872
	873	IF ( PRESENT(kdailyavtypes) ) THEN
	874	DO jobs = 1, kobsno
	875
[7992]	876	IF ( kobsqc(jobs) <= iqc_cutoff ) THEN
[2128]	877
	878	IF ( ( kobsstp(jobs) == (nit000 - 1) ).AND.&
	879	& ( ANY (kdailyavtypes(:) == ktyp(jobs)) ) ) THEN
[7992]	880	kobsqc(jobs) = IBSET(kobsqc(jobs),13)
[2128]	881	kotdobs = kotdobs + 1
	882	CYCLE
	883	ENDIF
	884
	885	ENDIF
	886	END DO
	887	ENDIF
	888
	889
	890	END SUBROUTINE obs_coo_tim_prof
	891
	892	SUBROUTINE obs_coo_grd( kobsno, kobsi, kobsj, kobsqc, kgrdobs )
	893	!!----------------------------------------------------------------------
	894	!! * ROUTINE obs_coo_grd *
	895	!!
	896	!! ** Purpose : Verify that the grid search has not failed
	897	!!
	898	!! ** Method : The previously computed i,j indeces are checked
	899	!!
	900	!! ** Action :
	901	!!
	902	!! References :
	903	!!
	904	!! History :
	905	!! ! 2007-01 (K. Mogensen) Original
	906	!!----------------------------------------------------------------------
	907	!! * Arguments
	908	INTEGER, INTENT(IN) :: kobsno ! Number of observations
	909	INTEGER, DIMENSION(kobsno), INTENT(IN ) :: &
	910	& kobsi, & ! i,j indeces previously computed
	911	& kobsj
	912	INTEGER, INTENT(INOUT) :: kgrdobs ! Number of observations failing the check
	913	INTEGER, DIMENSION(kobsno), INTENT(INOUT) :: &
	914	& kobsqc ! Quality control flag
	915
	916	!! * Local declarations
	917	INTEGER :: jobs ! Loop variable
	918
	919	! Flag if the grid search failed
	920
	921	DO jobs = 1, kobsno
	922	IF ( ( kobsi(jobs) <= 0 ) .AND. ( kobsj(jobs) <= 0 ) ) THEN
[7992]	923	kobsqc(jobs) = IBSET(kobsqc(jobs),12)
[2128]	924	kgrdobs = kgrdobs + 1
	925	ENDIF
	926	END DO
	927
	928	END SUBROUTINE obs_coo_grd
	929
	930	SUBROUTINE obs_coo_spc_2d( kobsno, kpi, kpj, &
	931	& kobsi, kobsj, pobslam, pobsphi, &
	932	& plam, pphi, pmask, &
	933	& kobsqc, kosdobs, klanobs, &
[7992]	934	& knlaobs,ld_nea, &
	935	& kbdyobs,ld_bound_reject, &
	936	& kqc_cutoff )
[2128]	937	!!----------------------------------------------------------------------
	938	!! * ROUTINE obs_coo_spc_2d *
	939	!!
	940	!! ** Purpose : Check for points outside the domain and land points
	941	!!
	942	!! ** Method : Remove the observations that are outside the model space
	943	!! and time domain or located within model land cells.
	944	!!
	945	!! ** Action :
	946	!!
	947	!! History :
	948	!! ! 2007-03 (A. Weaver, K. Mogensen) Original
	949	!! ! 2007-06 (K. Mogensen et al) Reject obs. near land.
	950	!!----------------------------------------------------------------------
	951	!! * Modules used
	952
	953	!! * Arguments
	954	INTEGER, INTENT(IN) :: kobsno ! Total number of observations
	955	INTEGER, INTENT(IN) :: kpi ! Number of grid points in (i,j)
	956	INTEGER, INTENT(IN) :: kpj
	957	INTEGER, DIMENSION(kobsno), INTENT(IN) :: &
	958	& kobsi, & ! Observation (i,j) coordinates
	959	& kobsj
	960	REAL(KIND=wp), DIMENSION(kobsno), INTENT(IN) :: &
	961	& pobslam, & ! Observation (lon,lat) coordinates
	962	& pobsphi
	963	REAL(KIND=wp), DIMENSION(kpi,kpj), INTENT(IN) :: &
	964	& plam, pphi ! Model (lon,lat) coordinates
	965	REAL(KIND=wp), DIMENSION(kpi,kpj), INTENT(IN) :: &
	966	& pmask ! Land mask array
	967	INTEGER, DIMENSION(kobsno), INTENT(INOUT) :: &
	968	& kobsqc ! Observation quality control
[7992]	969	INTEGER, INTENT(INOUT) :: kosdobs ! Observations outside space domain
	970	INTEGER, INTENT(INOUT) :: klanobs ! Observations within a model land cell
	971	INTEGER, INTENT(INOUT) :: knlaobs ! Observations near land
	972	INTEGER, INTENT(INOUT) :: kbdyobs ! Observations near boundary
	973	LOGICAL, INTENT(IN) :: ld_nea ! Flag observations near land
	974	LOGICAL, INTENT(IN) :: ld_bound_reject ! Flag observations near open boundary
	975	INTEGER, INTENT(IN) :: kqc_cutoff ! Cutoff QC value
	976
[2128]	977	!! * Local declarations
	978	REAL(KIND=wp), DIMENSION(2,2,kobsno) :: &
	979	& zgmsk ! Grid mask
[7992]	980	#if defined key_bdy
[2128]	981	REAL(KIND=wp), DIMENSION(2,2,kobsno) :: &
[7992]	982	& zbmsk ! Boundary mask
	983	REAL(KIND=wp), DIMENSION(jpi,jpj) :: zbdymask
	984	#endif
	985	REAL(KIND=wp), DIMENSION(2,2,kobsno) :: &
[2128]	986	& zglam, & ! Model longitude at grid points
	987	& zgphi ! Model latitude at grid points
	988	INTEGER, DIMENSION(2,2,kobsno) :: &
	989	& igrdi, & ! Grid i,j
	990	& igrdj
	991	LOGICAL :: lgridobs ! Is observation on a model grid point.
	992	INTEGER :: iig, ijg ! i,j of observation on model grid point.
	993	INTEGER :: jobs, ji, jj
	994
	995	! Get grid point indices
	996
	997	DO jobs = 1, kobsno
	998
	999	! For invalid points use 2,2
	1000
[7992]	1001	IF ( kobsqc(jobs) >= kqc_cutoff ) THEN
[2128]	1002
	1003	igrdi(1,1,jobs) = 1
	1004	igrdj(1,1,jobs) = 1
	1005	igrdi(1,2,jobs) = 1
	1006	igrdj(1,2,jobs) = 2
	1007	igrdi(2,1,jobs) = 2
	1008	igrdj(2,1,jobs) = 1
	1009	igrdi(2,2,jobs) = 2
	1010	igrdj(2,2,jobs) = 2
	1011
	1012	ELSE
	1013
	1014	igrdi(1,1,jobs) = kobsi(jobs)-1
	1015	igrdj(1,1,jobs) = kobsj(jobs)-1
	1016	igrdi(1,2,jobs) = kobsi(jobs)-1
	1017	igrdj(1,2,jobs) = kobsj(jobs)
	1018	igrdi(2,1,jobs) = kobsi(jobs)
	1019	igrdj(2,1,jobs) = kobsj(jobs)-1
	1020	igrdi(2,2,jobs) = kobsi(jobs)
	1021	igrdj(2,2,jobs) = kobsj(jobs)
	1022
	1023	ENDIF
	1024
	1025	END DO
[7992]	1026
	1027	#if defined key_bdy
	1028	! Create a mask grid points in boundary rim
	1029	IF (ld_bound_reject) THEN
	1030	zbdymask(:,:) = 1.0_wp
	1031	DO ji = 1, nb_bdy
	1032	DO jj = 1, idx_bdy(ji)%nblen(1)
	1033	zbdymask(idx_bdy(ji)%nbi(jj,1),idx_bdy(ji)%nbj(jj,1)) = 0.0_wp
	1034	ENDDO
	1035	ENDDO
	1036
[8035]	1037	CALL obs_int_comm_2d( 2, 2, kobsno, kpi, kpj, igrdi, igrdj, zbdymask, zbmsk )
[7992]	1038	ENDIF
	1039	#endif
[2128]	1040
[7992]	1041	CALL obs_int_comm_2d( 2, 2, kobsno, kpi, kpj, igrdi, igrdj, pmask, zgmsk )
	1042	CALL obs_int_comm_2d( 2, 2, kobsno, kpi, kpj, igrdi, igrdj, plam, zglam )
	1043	CALL obs_int_comm_2d( 2, 2, kobsno, kpi, kpj, igrdi, igrdj, pphi, zgphi )
[2128]	1044
	1045	DO jobs = 1, kobsno
	1046
	1047	! Skip bad observations
[7992]	1048	IF ( kobsqc(jobs) >= kqc_cutoff ) CYCLE
[2128]	1049
	1050	! Flag if the observation falls outside the model spatial domain
	1051	IF ( ( pobslam(jobs) < -180. ) &
	1052	& .OR. ( pobslam(jobs) > 180. ) &
	1053	& .OR. ( pobsphi(jobs) < -90. ) &
	1054	& .OR. ( pobsphi(jobs) > 90. ) ) THEN
[7992]	1055	kobsqc(jobs) = IBSET(kobsqc(jobs),11)
[2128]	1056	kosdobs = kosdobs + 1
	1057	CYCLE
	1058	ENDIF
	1059
	1060	! Flag if the observation falls with a model land cell
	1061	IF ( SUM( zgmsk(1:2,1:2,jobs) ) == 0.0_wp ) THEN
[7992]	1062	kobsqc(jobs) = IBSET(kobsqc(jobs),10)
[2128]	1063	klanobs = klanobs + 1
	1064	CYCLE
	1065	ENDIF
	1066
	1067	! Check if this observation is on a grid point
	1068
	1069	lgridobs = .FALSE.
	1070	iig = -1
	1071	ijg = -1
	1072	DO jj = 1, 2
	1073	DO ji = 1, 2
	1074	IF ( ( ABS( zgphi(ji,jj,jobs) - pobsphi(jobs) ) < 1.0e-6_wp ) &
	1075	& .AND. &
[8222]	1076	& ( ABS( MOD( zglam(ji,jj,jobs) - pobslam(jobs),360.0) ) &
	1077	& < 1.0e-6_wp ) ) THEN
[2128]	1078	lgridobs = .TRUE.
	1079	iig = ji
	1080	ijg = jj
	1081	ENDIF
	1082	END DO
	1083	END DO
[7992]	1084
[2128]	1085	IF (lgridobs) THEN
	1086	IF (zgmsk(iig,ijg,jobs) == 0.0_wp ) THEN
[7992]	1087	kobsqc(jobs) = IBSET(kobsqc(jobs),10)
[2128]	1088	klanobs = klanobs + 1
	1089	CYCLE
	1090	ENDIF
	1091	ENDIF
[7992]	1092
	1093
[2128]	1094	! Flag if the observation falls is close to land
	1095	IF ( MINVAL( zgmsk(1:2,1:2,jobs) ) == 0.0_wp) THEN
	1096	knlaobs = knlaobs + 1
[7992]	1097	IF (ld_nea) THEN
	1098	kobsqc(jobs) = IBSET(kobsqc(jobs),9)
	1099	CYCLE
	1100	ENDIF
[2128]	1101	ENDIF
[7992]	1102
	1103	#if defined key_bdy
	1104	! Flag if the observation falls close to the boundary rim
	1105	IF (ld_bound_reject) THEN
	1106	IF ( MINVAL( zbmsk(1:2,1:2,jobs) ) == 0.0_wp ) THEN
	1107	kobsqc(jobs) = IBSET(kobsqc(jobs),8)
	1108	kbdyobs = kbdyobs + 1
	1109	CYCLE
	1110	ENDIF
	1111	! for observations on the grid...
	1112	IF (lgridobs) THEN
	1113	IF (zbmsk(iig,ijg,jobs) == 0.0_wp ) THEN
	1114	kobsqc(jobs) = IBSET(kobsqc(jobs),8)
	1115	kbdyobs = kbdyobs + 1
	1116	CYCLE
	1117	ENDIF
	1118	ENDIF
	1119	ENDIF
	1120	#endif
[2128]	1121
	1122	END DO
	1123
	1124	END SUBROUTINE obs_coo_spc_2d
	1125
	1126	SUBROUTINE obs_coo_spc_3d( kprofno, kobsno, kpstart, kpend, &
	1127	& kpi, kpj, kpk, &
	1128	& kobsi, kobsj, kobsk, &
	1129	& pobslam, pobsphi, pobsdep, &
	1130	& plam, pphi, pdep, pmask, &
	1131	& kpobsqc, kobsqc, kosdobs, &
[7992]	1132	& klanobs, knlaobs, ld_nea, &
	1133	& kbdyobs, ld_bound_reject, &
	1134	& kqc_cutoff )
[2128]	1135	!!----------------------------------------------------------------------
	1136	!! * ROUTINE obs_coo_spc_3d *
	1137	!!
	1138	!! ** Purpose : Check for points outside the domain and land points
	1139	!! Reset depth of observation above highest model level
	1140	!! to the value of highest model level
	1141	!!
	1142	!! ** Method : Remove the observations that are outside the model space
	1143	!! and time domain or located within model land cells.
	1144	!!
	1145	!! NB. T and S profile observations lying between the ocean
	1146	!! surface and the depth of the first model T point are
	1147	!! assigned a depth equal to that of the first model T pt.
	1148	!!
	1149	!! ** Action :
	1150	!!
	1151	!! History :
	1152	!! ! 2007-01 (K. Mogensen) Rewrite of parts of obs_scr
	1153	!! ! 2007-06 (K. Mogensen et al) Reject obs. near land.
	1154	!!----------------------------------------------------------------------
	1155	!! * Modules used
	1156	USE dom_oce, ONLY : & ! Geographical information
[7992]	1157	& gdepw_1d, &
	1158	& gdepw_0, &
	1159	#if defined key_vvl
	1160	& gdepw_n, &
	1161	& gdept_n, &
	1162	#endif
	1163	& ln_zco, &
	1164	& ln_zps, &
	1165	& lk_vvl
[2128]	1166
	1167	!! * Arguments
	1168	INTEGER, INTENT(IN) :: kprofno ! Number of profiles
	1169	INTEGER, INTENT(IN) :: kobsno ! Total number of observations
	1170	INTEGER, INTENT(IN) :: kpi ! Number of grid points in (i,j,k)
	1171	INTEGER, INTENT(IN) :: kpj
	1172	INTEGER, INTENT(IN) :: kpk
	1173	INTEGER, DIMENSION(kprofno), INTENT(IN) :: &
	1174	& kpstart, & ! Start of individual profiles
	1175	& kpend ! End of individual profiles
	1176	INTEGER, DIMENSION(kprofno), INTENT(IN) :: &
	1177	& kobsi, & ! Observation (i,j) coordinates
	1178	& kobsj
	1179	INTEGER, DIMENSION(kobsno), INTENT(IN) :: &
	1180	& kobsk ! Observation k coordinate
	1181	REAL(KIND=wp), DIMENSION(kprofno), INTENT(IN) :: &
	1182	& pobslam, & ! Observation (lon,lat) coordinates
	1183	& pobsphi
	1184	REAL(KIND=wp), DIMENSION(kobsno), INTENT(INOUT) :: &
	1185	& pobsdep ! Observation depths
	1186	REAL(KIND=wp), DIMENSION(kpi,kpj), INTENT(IN) :: &
	1187	& plam, pphi ! Model (lon,lat) coordinates
	1188	REAL(KIND=wp), DIMENSION(kpk), INTENT(IN) :: &
	1189	& pdep ! Model depth coordinates
	1190	REAL(KIND=wp), DIMENSION(kpi,kpj,kpk), INTENT(IN) :: &
	1191	& pmask ! Land mask array
	1192	INTEGER, DIMENSION(kprofno), INTENT(INOUT) :: &
	1193	& kpobsqc ! Profile quality control
	1194	INTEGER, DIMENSION(kobsno), INTENT(INOUT) :: &
	1195	& kobsqc ! Observation quality control
	1196	INTEGER, INTENT(INOUT) :: kosdobs ! Observations outside space domain
	1197	INTEGER, INTENT(INOUT) :: klanobs ! Observations within a model land cell
	1198	INTEGER, INTENT(INOUT) :: knlaobs ! Observations near land
[7992]	1199	INTEGER, INTENT(INOUT) :: kbdyobs ! Observations near boundary
[2128]	1200	LOGICAL, INTENT(IN) :: ld_nea ! Flag observations near land
[7992]	1201	LOGICAL, INTENT(IN) :: ld_bound_reject ! Flag observations near open boundary
	1202	INTEGER, INTENT(IN) :: kqc_cutoff ! Cutoff QC value
	1203
[2128]	1204	!! * Local declarations
	1205	REAL(KIND=wp), DIMENSION(2,2,kpk,kprofno) :: &
	1206	& zgmsk ! Grid mask
[7992]	1207	#if defined key_bdy
[2128]	1208	REAL(KIND=wp), DIMENSION(2,2,kprofno) :: &
[7992]	1209	& zbmsk ! Boundary mask
	1210	REAL(KIND=wp), DIMENSION(jpi,jpj) :: zbdymask
	1211	#endif
	1212	REAL(KIND=wp), DIMENSION(2,2,kpk,kprofno) :: &
[11461]	1213	& zgdept, &
[7992]	1214	& zgdepw
	1215	REAL(KIND=wp), DIMENSION(2,2,kprofno) :: &
[2128]	1216	& zglam, & ! Model longitude at grid points
[12506]	1217	& zgphi, & ! Model latitude at grid points
	1218	& zbathy ! Index of deepest wet level at grid points
[2128]	1219	INTEGER, DIMENSION(2,2,kprofno) :: &
	1220	& igrdi, & ! Grid i,j
	1221	& igrdj
	1222	LOGICAL :: lgridobs ! Is observation on a model grid point.
[7992]	1223	LOGICAL :: ll_next_to_land ! Is a profile next to land
[2128]	1224	INTEGER :: iig, ijg ! i,j of observation on model grid point.
	1225	INTEGER :: jobs, jobsp, jk, ji, jj
[12506]	1226	REAL(KIND=wp) :: maxdept, maxdepw
[2128]	1227
	1228	! Get grid point indices
	1229
	1230	DO jobs = 1, kprofno
	1231
	1232	! For invalid points use 2,2
	1233
[7992]	1234	IF ( kpobsqc(jobs) >= kqc_cutoff ) THEN
[2128]	1235
	1236	igrdi(1,1,jobs) = 1
	1237	igrdj(1,1,jobs) = 1
	1238	igrdi(1,2,jobs) = 1
	1239	igrdj(1,2,jobs) = 2
	1240	igrdi(2,1,jobs) = 2
	1241	igrdj(2,1,jobs) = 1
	1242	igrdi(2,2,jobs) = 2
	1243	igrdj(2,2,jobs) = 2
	1244
	1245	ELSE
	1246
	1247	igrdi(1,1,jobs) = kobsi(jobs)-1
	1248	igrdj(1,1,jobs) = kobsj(jobs)-1
	1249	igrdi(1,2,jobs) = kobsi(jobs)-1
	1250	igrdj(1,2,jobs) = kobsj(jobs)
	1251	igrdi(2,1,jobs) = kobsi(jobs)
	1252	igrdj(2,1,jobs) = kobsj(jobs)-1
	1253	igrdi(2,2,jobs) = kobsi(jobs)
	1254	igrdj(2,2,jobs) = kobsj(jobs)
	1255
	1256	ENDIF
	1257
	1258	END DO
[7992]	1259
	1260	#if defined key_bdy
	1261	! Create a mask grid points in boundary rim
	1262	IF (ld_bound_reject) THEN
	1263	zbdymask(:,:) = 1.0_wp
	1264	DO ji = 1, nb_bdy
	1265	DO jj = 1, idx_bdy(ji)%nblen(1)
	1266	zbdymask(idx_bdy(ji)%nbi(jj,1),idx_bdy(ji)%nbj(jj,1)) = 0.0_wp
	1267	ENDDO
	1268	ENDDO
	1269	ENDIF
	1270
[8035]	1271	CALL obs_int_comm_2d( 2, 2, kprofno, kpi, kpj, igrdi, igrdj, zbdymask, zbmsk )
[7992]	1272	#endif
[2128]	1273
[7992]	1274	CALL obs_int_comm_3d( 2, 2, kprofno, kpi, kpj, kpk, igrdi, igrdj, pmask, zgmsk )
	1275	CALL obs_int_comm_2d( 2, 2, kprofno, kpi, kpj, igrdi, igrdj, plam, zglam )
	1276	CALL obs_int_comm_2d( 2, 2, kprofno, kpi, kpj, igrdi, igrdj, pphi, zgphi )
[12506]	1277	CALL obs_int_comm_2d( 2, 2, kprofno, kpi, kpj, igrdi, igrdj, REAL(mbathy), zbathy )
[7992]	1278	! Need to know the bathy depth for each observation for sco
	1279	CALL obs_int_comm_3d( 2, 2, kprofno, kpi, kpj, kpk, igrdi, igrdj, fsdepw(:,:,:), &
	1280	& zgdepw )
[11461]	1281	CALL obs_int_comm_3d( 2, 2, kprofno, kpi, kpj, kpk, igrdi, igrdj, fsdept(:,:,:), &
	1282	& zgdept )
[2128]	1283
	1284	DO jobs = 1, kprofno
	1285
	1286	! Skip bad profiles
[7992]	1287	IF ( kpobsqc(jobs) >= kqc_cutoff ) CYCLE
[2128]	1288
	1289	! Check if this observation is on a grid point
	1290
	1291	lgridobs = .FALSE.
	1292	iig = -1
	1293	ijg = -1
	1294	DO jj = 1, 2
	1295	DO ji = 1, 2
	1296	IF ( ( ABS( zgphi(ji,jj,jobs) - pobsphi(jobs) ) < 1.0e-6_wp ) &
	1297	& .AND. &
[8222]	1298	& ( ABS( MOD( zglam(ji,jj,jobs) - pobslam(jobs),360.0) ) < 1.0e-6_wp ) &
[2128]	1299	& ) THEN
	1300	lgridobs = .TRUE.
	1301	iig = ji
	1302	ijg = jj
	1303	ENDIF
	1304	END DO
	1305	END DO
	1306
[7992]	1307	! Check if next to land
	1308	IF ( ANY( zgmsk(1:2,1:2,1,jobs) == 0.0_wp ) ) THEN
	1309	ll_next_to_land=.TRUE.
	1310	ELSE
	1311	ll_next_to_land=.FALSE.
	1312	ENDIF
	1313
[2128]	1314	! Reject observations
	1315
	1316	DO jobsp = kpstart(jobs), kpend(jobs)
	1317
[12506]	1318	! Calculate max T and W depths of 2x2 grid
	1319	maxdept=zgdept(1,1,NINT(zbathy(1,1,jobs)),jobs)
	1320	maxdepw=zgdepw(1,1,NINT(zbathy(1,1,jobs))+1,jobs)
	1321	DO jj = 1, 2
	1322	DO ji = 1, 2
	1323	IF ( zgdept(ji,jj,NINT(zbathy(ji,jj,jobs)),jobs) > maxdept ) THEN
	1324	maxdept = zgdept(ji,jj,NINT(zbathy(ji,jj,jobs)),jobs)
	1325	END IF
	1326	IF ( zgdepw(ji,jj,NINT(zbathy(ji,jj,jobs))+1,jobs) > maxdepw ) THEN
	1327	maxdepw = zgdepw(ji,jj,NINT(zbathy(ji,jj,jobs))+1,jobs)
	1328	END IF
	1329	END DO
	1330	END DO
	1331
[2128]	1332	! Flag if the observation falls outside the model spatial domain
	1333	IF ( ( pobslam(jobs) < -180. ) &
	1334	& .OR. ( pobslam(jobs) > 180. ) &
	1335	& .OR. ( pobsphi(jobs) < -90. ) &
	1336	& .OR. ( pobsphi(jobs) > 90. ) &
	1337	& .OR. ( pobsdep(jobsp) < 0.0 ) &
[12506]	1338	& .OR. ( pobsdep(jobsp) >= maxdepw ) ) THEN
[7992]	1339	kobsqc(jobsp) = IBSET(kobsqc(jobsp),11)
[2128]	1340	kosdobs = kosdobs + 1
	1341	CYCLE
	1342	ENDIF
	1343
[7992]	1344	! To check if an observations falls within land there are two cases:
	1345	! 1: z-coordibnates, where the check uses the mask
	1346	! 2: terrain following (eg s-coordinates),
	1347	! where we use the depth of the bottom cell to mask observations
	1348
	1349	IF( (.NOT. lk_vvl) .AND. ( ln_zps .OR. ln_zco ) ) THEN !(CASE 1)
	1350
	1351	! Flag if the observation falls with a model land cell
	1352	IF ( SUM( zgmsk(1:2,1:2,kobsk(jobsp)-1:kobsk(jobsp),jobs) ) &
	1353	& == 0.0_wp ) THEN
	1354	kobsqc(jobsp) = IBSET(kobsqc(jobsp),10)
	1355	klanobs = klanobs + 1
	1356	CYCLE
	1357	ENDIF
	1358
	1359	! Flag if the observation is close to land
	1360	IF ( MINVAL( zgmsk(1:2,1:2,kobsk(jobsp)-1:kobsk(jobsp),jobs) ) == &
	1361	& 0.0_wp) THEN
	1362	knlaobs = knlaobs + 1
	1363	IF (ld_nea) THEN
	1364	kobsqc(jobsp) = IBSET(kobsqc(jobsp),10)
	1365	ENDIF
	1366	ENDIF
	1367
	1368	ELSE ! Case 2
	1369	! Flag if the observation is deeper than the bathymetry
	1370	! Or if it is within the mask
	1371	IF ( ANY( zgdepw(1:2,1:2,kpk,jobs) < pobsdep(jobsp) ) &
	1372	& .OR. &
	1373	& ( SUM( zgmsk(1:2,1:2,kobsk(jobsp)-1:kobsk(jobsp),jobs) ) &
	1374	& == 0.0_wp) ) THEN
	1375	kobsqc(jobsp) = IBSET(kobsqc(jobsp),10)
	1376	klanobs = klanobs + 1
	1377	CYCLE
	1378	ENDIF
	1379
	1380	! Flag if the observation is close to land
	1381	IF ( ll_next_to_land ) THEN
	1382	knlaobs = knlaobs + 1
	1383	IF (ld_nea) THEN
	1384	kobsqc(jobsp) = IBSET(kobsqc(jobsp),10)
	1385	ENDIF
	1386	ENDIF
	1387
[2128]	1388	ENDIF
	1389
	1390	! For observations on the grid reject them if their are at
	1391	! a masked point
	1392
	1393	IF (lgridobs) THEN
	1394	IF (zgmsk(iig,ijg,kobsk(jobsp)-1,jobs) == 0.0_wp ) THEN
[7992]	1395	kobsqc(jobsp) = IBSET(kobsqc(jobs),10)
[2128]	1396	klanobs = klanobs + 1
	1397	CYCLE
	1398	ENDIF
	1399	ENDIF
	1400
[11443]	1401	! Set observation depth equal to that of the first model depth
[11461]	1402	IF ( pobsdep(jobsp) < MINVAL(zgdept(1:2,1:2,1,jobs) ) ) THEN
	1403	pobsdep(jobsp) = MINVAL(zgdept(1:2,1:2,1,jobs))
[11443]	1404	ENDIF
[12506]	1405
	1406	! Set observation depth equal to that of the last wet T-point
	1407	IF ( ( pobsdep(jobsp) > maxdept ) .AND. &
	1408	& ( pobsdep(jobsp) < maxdepw ) ) THEN
	1409	pobsdep(jobsp) = maxdept
	1410	END IF
[11443]	1411
[7992]	1412	#if defined key_bdy
	1413	! Flag if the observation falls close to the boundary rim
	1414	IF (ld_bound_reject) THEN
	1415	IF ( MINVAL( zbmsk(1:2,1:2,jobs) ) == 0.0_wp ) THEN
	1416	kobsqc(jobsp) = IBSET(kobsqc(jobs),8)
	1417	kbdyobs = kbdyobs + 1
	1418	CYCLE
	1419	ENDIF
	1420	! for observations on the grid...
	1421	IF (lgridobs) THEN
	1422	IF (zbmsk(iig,ijg,jobs) == 0.0_wp ) THEN
	1423	kobsqc(jobsp) = IBSET(kobsqc(jobs),8)
	1424	kbdyobs = kbdyobs + 1
	1425	CYCLE
	1426	ENDIF
	1427	ENDIF
	1428	ENDIF
	1429	#endif
	1430
[2128]	1431	END DO
	1432	END DO
	1433
	1434	END SUBROUTINE obs_coo_spc_3d
	1435
[7992]	1436	SUBROUTINE obs_pro_rej( profdata, kqc_cutoff )
[2128]	1437	!!----------------------------------------------------------------------
	1438	!! * ROUTINE obs_pro_rej *
	1439	!!
	1440	!! ** Purpose : Reject all data within a rejected profile
	1441	!!
	1442	!! ** Method :
	1443	!!
	1444	!! ** Action :
	1445	!!
	1446	!! References :
	1447	!!
	1448	!! History :
	1449	!! ! 2007-10 (K. Mogensen) Original code
	1450	!!----------------------------------------------------------------------
	1451	!! * Modules used
	1452	!! * Arguments
	1453	TYPE(obs_prof), INTENT(INOUT) :: profdata ! Profile data
[7992]	1454	INTEGER, INTENT(IN) :: kqc_cutoff ! QC cutoff value
	1455
[2128]	1456	!! * Local declarations
	1457	INTEGER :: jprof
	1458	INTEGER :: jvar
	1459	INTEGER :: jobs
	1460
	1461	! Loop over profiles
	1462
	1463	DO jprof = 1, profdata%nprof
	1464
[7992]	1465	IF ( profdata%nqc(jprof) > kqc_cutoff ) THEN
[2128]	1466
	1467	DO jvar = 1, profdata%nvar
	1468
	1469	DO jobs = profdata%npvsta(jprof,jvar), &
	1470	& profdata%npvend(jprof,jvar)
	1471
	1472	profdata%var(jvar)%nvqc(jobs) = &
[7992]	1473	& IBSET(profdata%var(jvar)%nvqc(jobs),14)
[2128]	1474
	1475	END DO
	1476
	1477	END DO
	1478
	1479	ENDIF
	1480
	1481	END DO
	1482
	1483	END SUBROUTINE obs_pro_rej
	1484
[13393]	1485	SUBROUTINE obs_uv_rej_pro( profdata, knumu, knumv, kqc_cutoff )
[2128]	1486	!!----------------------------------------------------------------------
[13393]	1487	!! * ROUTINE obs_uv_rej_pro *
[2128]	1488	!!
	1489	!! ** Purpose : Reject u if v is rejected and vice versa
	1490	!!
	1491	!! ** Method :
	1492	!!
	1493	!! ** Action :
	1494	!!
	1495	!! References :
	1496	!!
	1497	!! History :
	1498	!! ! 2009-2 (K. Mogensen) Original code
	1499	!!----------------------------------------------------------------------
	1500	!! * Modules used
	1501	!! * Arguments
	1502	TYPE(obs_prof), INTENT(INOUT) :: profdata ! Profile data
	1503	INTEGER, INTENT(INOUT) :: knumu ! Number of u rejected
	1504	INTEGER, INTENT(INOUT) :: knumv ! Number of v rejected
[7992]	1505	INTEGER, INTENT(IN) :: kqc_cutoff ! QC cutoff value
	1506
[2128]	1507	!! * Local declarations
	1508	INTEGER :: jprof
	1509	INTEGER :: jvar
	1510	INTEGER :: jobs
	1511
	1512	! Loop over profiles
	1513
	1514	DO jprof = 1, profdata%nprof
	1515
	1516	IF ( ( profdata%npvsta(jprof,1) /= profdata%npvsta(jprof,2) ) .OR. &
	1517	& ( profdata%npvend(jprof,1) /= profdata%npvend(jprof,2) ) ) THEN
	1518
	1519	CALL ctl_stop('U,V profiles inconsistent in obs_uv_rej')
	1520	RETURN
	1521
	1522	ENDIF
	1523
	1524	DO jobs = profdata%npvsta(jprof,1), profdata%npvend(jprof,1)
	1525
[7992]	1526	IF ( ( profdata%var(1)%nvqc(jobs) > kqc_cutoff ) .AND. &
	1527	& ( profdata%var(2)%nvqc(jobs) <= kqc_cutoff) ) THEN
	1528	profdata%var(2)%nvqc(jobs) = IBSET(profdata%var(1)%nvqc(jobs),15)
[2128]	1529	knumv = knumv + 1
	1530	ENDIF
[7992]	1531	IF ( ( profdata%var(2)%nvqc(jobs) > kqc_cutoff ) .AND. &
	1532	& ( profdata%var(1)%nvqc(jobs) <= kqc_cutoff) ) THEN
	1533	profdata%var(1)%nvqc(jobs) = IBSET(profdata%var(1)%nvqc(jobs),15)
[2128]	1534	knumu = knumu + 1
	1535	ENDIF
	1536
	1537	END DO
	1538
	1539	END DO
	1540
[13393]	1541	END SUBROUTINE obs_uv_rej_pro
[2128]	1542
	1543	END MODULE obs_prep

Note: See TracBrowser for help on using the repository browser.

Download in other formats: