Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Normal
Revision Log

obs_read_prof.F90 @ 11863

Last change on this file since 11863 was 11863, checked in by kingr, 4 years ago
Merged recent changes to dev_r5518_obs_oper_update
File size: 31.4 KB

Rev	Line
[2128]	1	MODULE obs_read_prof
	2	!!======================================================================
	3	!! * MODULE obs_read_prof *
	4	!! Observation diagnostics: Read the T and S profile observations
	5	!!======================================================================
	6
	7	!!----------------------------------------------------------------------
	8	!! obs_rea_pro_dri : Driver for reading profile obs
	9	!!----------------------------------------------------------------------
	10
	11	!! * Modules used
	12	USE par_kind ! Precision variables
	13	USE par_oce ! Ocean parameters
	14	USE in_out_manager ! I/O manager
	15	USE dom_oce ! Ocean space and time domain variables
	16	USE obs_mpp ! MPP support routines for observation diagnostics
	17	USE julian ! Julian date routines
	18	USE obs_utils ! Observation operator utility functions
	19	USE obs_prep ! Prepare observation arrays
	20	USE obs_grid ! Grid search
	21	USE obs_sort ! Sorting observation arrays
	22	USE obs_profiles_def ! Profile definitions
	23	USE obs_conv ! Various conversion routines
	24	USE obs_types ! Observation type definitions
	25	USE netcdf ! NetCDF library
	26	USE obs_oper ! Observation operators
[2715]	27	USE lib_mpp ! For ctl_warn/stop
[7992]	28	USE obs_fbm ! Feedback routines
[2128]	29
	30	IMPLICIT NONE
	31
	32	!! * Routine accessibility
	33	PRIVATE
	34
[7992]	35	PUBLIC obs_rea_prof ! Read the profile observations
[2128]	36
[2287]	37	!!----------------------------------------------------------------------
	38	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
	39	!! $Id$
	40	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
	41	!!----------------------------------------------------------------------
	42
[2128]	43	CONTAINS
[7992]	44
	45	SUBROUTINE obs_rea_prof( profdata, knumfiles, cdfilenames, &
	46	& kvars, kextr, kstp, ddobsini, ddobsend, &
[9306]	47	& ldvar, ldignmis, ldsatt, &
[11863]	48	& ldmod, ldclim, cdvars, kdailyavtypes )
[2128]	49	!!---------------------------------------------------------------------
	50	!!
[7992]	51	!! * ROUTINE obs_rea_prof *
[2128]	52	!!
	53	!! ** Purpose : Read from file the profile observations
	54	!!
[7992]	55	!! ** Method : Read feedback data in and transform to NEMO internal
	56	!! profile data structure
[2128]	57	!!
	58	!! ** Action :
	59	!!
	60	!! References :
	61	!!
	62	!! History :
	63	!! ! : 2009-09 (K. Mogensen) : New merged version of old routines
[7992]	64	!! ! : 2015-08 (M. Martin) : Merged profile and velocity routines
[2128]	65	!!----------------------------------------------------------------------
[7992]	66
[2128]	67	!! * Arguments
[7992]	68	TYPE(obs_prof), INTENT(OUT) :: &
	69	& profdata ! Profile data to be read
	70	INTEGER, INTENT(IN) :: knumfiles ! Number of files to read
[2128]	71	CHARACTER(LEN=128), INTENT(IN) :: &
[7992]	72	& cdfilenames(knumfiles) ! File names to read in
[2128]	73	INTEGER, INTENT(IN) :: kvars ! Number of variables in profdata
[7992]	74	INTEGER, INTENT(IN) :: kextr ! Number of extra fields for each var
	75	INTEGER, INTENT(IN) :: kstp ! Ocean time-step index
[9306]	76	LOGICAL, DIMENSION(kvars), INTENT(IN) :: ldvar ! Observed variables switches
[7992]	77	LOGICAL, INTENT(IN) :: ldignmis ! Ignore missing files
	78	LOGICAL, INTENT(IN) :: ldsatt ! Compute salinity at all temperature points
	79	LOGICAL, INTENT(IN) :: ldmod ! Initialize model from input data
[11468]	80	LOGICAL, INTENT(IN) :: ldclim ! Set flag to show climatology will be output
[7992]	81	REAL(dp), INTENT(IN) :: ddobsini ! Obs. ini time in YYYYMMDD.HHMMSS
	82	REAL(dp), INTENT(IN) :: ddobsend ! Obs. end time in YYYYMMDD.HHMMSS
[11863]	83	CHARACTER(len=8), DIMENSION(kvars), INTENT(IN) :: cdvars
[2128]	84	INTEGER, DIMENSION(imaxavtypes), OPTIONAL :: &
[7992]	85	& kdailyavtypes ! Types of daily average observations
[2128]	86
	87	!! * Local declarations
[7992]	88	CHARACTER(LEN=15), PARAMETER :: cpname='obs_rea_prof'
	89	CHARACTER(len=8) :: clrefdate
[11863]	90	CHARACTER(len=8), DIMENSION(:), ALLOCATABLE :: clvarsin
[2128]	91	INTEGER :: jvar
	92	INTEGER :: ji
	93	INTEGER :: jj
	94	INTEGER :: jk
	95	INTEGER :: ij
	96	INTEGER :: iflag
	97	INTEGER :: inobf
	98	INTEGER :: i_file_id
	99	INTEGER :: inowin
	100	INTEGER :: iyea
	101	INTEGER :: imon
	102	INTEGER :: iday
	103	INTEGER :: ihou
	104	INTEGER :: imin
	105	INTEGER :: isec
[7992]	106	INTEGER :: iprof
	107	INTEGER :: iproftot
[9306]	108	INTEGER, DIMENSION(kvars) :: ivart0
	109	INTEGER, DIMENSION(kvars) :: ivart
[7992]	110	INTEGER :: ip3dt
	111	INTEGER :: ios
	112	INTEGER :: ioserrcount
[9306]	113	INTEGER, DIMENSION(kvars) :: ivartmpp
[7992]	114	INTEGER :: ip3dtmpp
	115	INTEGER :: itype
[2128]	116	INTEGER, DIMENSION(knumfiles) :: &
	117	& irefdate
[9306]	118	INTEGER, DIMENSION(ntyp1770+1,kvars) :: &
	119	& itypvar, &
	120	& itypvarmpp
	121	INTEGER, DIMENSION(:,:), ALLOCATABLE :: &
	122	& iobsi, &
	123	& iobsj, &
	124	& iproc
[2128]	125	INTEGER, DIMENSION(:), ALLOCATABLE :: &
	126	& iindx, &
	127	& ifileidx, &
	128	& iprofidx
	129	INTEGER, DIMENSION(imaxavtypes) :: &
	130	& idailyavtypes
[7992]	131	INTEGER, DIMENSION(kvars) :: &
	132	& iv3dt
[2128]	133	REAL(wp), DIMENSION(:), ALLOCATABLE :: &
	134	& zphi, &
	135	& zlam
[7992]	136	REAL(wp), DIMENSION(:), ALLOCATABLE :: &
[2128]	137	& zdat
[7992]	138	REAL(wp), DIMENSION(knumfiles) :: &
	139	& djulini, &
	140	& djulend
[2128]	141	LOGICAL :: llvalprof
[7992]	142	LOGICAL :: lldavtimset
[9306]	143	LOGICAL :: llcycle
[2128]	144	TYPE(obfbdata), POINTER, DIMENSION(:) :: &
	145	& inpfiles
[7992]	146
[2128]	147	! Local initialization
	148	iprof = 0
[9306]	149	ivart0(:) = 0
[2128]	150	ip3dt = 0
	151
	152	! Daily average types
[7992]	153	lldavtimset = .FALSE.
[2128]	154	IF ( PRESENT(kdailyavtypes) ) THEN
	155	idailyavtypes(:) = kdailyavtypes(:)
[7992]	156	IF ( ANY (idailyavtypes(:) /= -1) ) lldavtimset = .TRUE.
[2128]	157	ELSE
	158	idailyavtypes(:) = -1
	159	ENDIF
	160
	161	!-----------------------------------------------------------------------
[7992]	162	! Count the number of files needed and allocate the obfbdata type
[2128]	163	!-----------------------------------------------------------------------
	164
	165	inobf = knumfiles
[7992]	166
[2128]	167	ALLOCATE( inpfiles(inobf) )
	168
	169	prof_files : DO jj = 1, inobf
[7992]	170
[2128]	171	!---------------------------------------------------------------------
	172	! Prints
	173	!---------------------------------------------------------------------
	174	IF(lwp) THEN
	175	WRITE(numout,*)
	176	WRITE(numout,*) ' obs_rea_pro_dri : Reading from file = ', &
[7992]	177	& TRIM( TRIM( cdfilenames(jj) ) )
[2128]	178	WRITE(numout,*) ' ~~~~~~~~~~~~~~~'
	179	WRITE(numout,*)
	180	ENDIF
	181
	182	!---------------------------------------------------------------------
	183	! Initialization: Open file and get dimensions only
	184	!---------------------------------------------------------------------
[7992]	185
	186	iflag = nf90_open( TRIM( cdfilenames(jj) ), nf90_nowrite, &
[2128]	187	& i_file_id )
[7992]	188
[2128]	189	IF ( iflag /= nf90_noerr ) THEN
	190
	191	IF ( ldignmis ) THEN
	192	inpfiles(jj)%nobs = 0
[7992]	193	CALL ctl_warn( 'File ' // TRIM( cdfilenames(jj) ) // &
[2128]	194	& ' not found' )
	195	ELSE
[7992]	196	CALL ctl_stop( 'File ' // TRIM( cdfilenames(jj) ) // &
[2128]	197	& ' not found' )
	198	ENDIF
	199
	200	ELSE
[7992]	201
[2128]	202	!------------------------------------------------------------------
[7992]	203	! Close the file since it is opened in read_obfbdata
[2128]	204	!------------------------------------------------------------------
[7992]	205
[2128]	206	iflag = nf90_close( i_file_id )
	207
	208	!------------------------------------------------------------------
	209	! Read the profile file into inpfiles
	210	!------------------------------------------------------------------
[7992]	211	CALL init_obfbdata( inpfiles(jj) )
	212	CALL read_obfbdata( TRIM( cdfilenames(jj) ), inpfiles(jj), &
	213	& ldgrid = .TRUE. )
	214
[9306]	215	IF ( inpfiles(jj)%nvar /= kvars ) THEN
[7992]	216	CALL ctl_stop( 'Feedback format error: ', &
[9306]	217	& ' unexpected number of vars in profile file' )
[7992]	218	ENDIF
	219
	220	IF ( ldmod .AND. ( inpfiles(jj)%nadd == 0 ) ) THEN
	221	CALL ctl_stop( 'Model not in input data' )
	222	ENDIF
	223
	224	IF ( jj == 1 ) THEN
[11863]	225	ALLOCATE( clvarsin( inpfiles(jj)%nvar ) )
[7992]	226	DO ji = 1, inpfiles(jj)%nvar
[11863]	227	clvarsin(ji) = inpfiles(jj)%cname(ji)
	228	IF ( clvarsin(ji) /= cdvars(ji) ) THEN
	229	CALL ctl_stop( 'Feedback file variables do not match', &
	230	& ' expected variable names for this type' )
	231	ENDIF
[7992]	232	END DO
[2128]	233	ELSE
[7992]	234	DO ji = 1, inpfiles(jj)%nvar
[11863]	235	IF ( inpfiles(jj)%cname(ji) /= clvarsin(ji) ) THEN
[7992]	236	CALL ctl_stop( 'Feedback file variables not consistent', &
	237	& ' with previous files for this type' )
	238	ENDIF
	239	END DO
[2128]	240	ENDIF
[7992]	241
[2128]	242	!------------------------------------------------------------------
	243	! Change longitude (-180,180)
	244	!------------------------------------------------------------------
	245
	246	DO ji = 1, inpfiles(jj)%nobs
	247
	248	IF ( inpfiles(jj)%plam(ji) < -180. ) &
	249	& inpfiles(jj)%plam(ji) = inpfiles(jj)%plam(ji) + 360.
	250
	251	IF ( inpfiles(jj)%plam(ji) > 180. ) &
	252	& inpfiles(jj)%plam(ji) = inpfiles(jj)%plam(ji) - 360.
	253
	254	END DO
	255
	256	!------------------------------------------------------------------
	257	! Calculate the date (change eventually)
	258	!------------------------------------------------------------------
[7992]	259	clrefdate=inpfiles(jj)%cdjuldref(1:8)
	260	READ(clrefdate,'(I8)') irefdate(jj)
	261
[2128]	262	CALL ddatetoymdhms( ddobsini, iyea, imon, iday, ihou, imin, isec )
	263	CALL greg2jul( isec, imin, ihou, iday, imon, iyea, djulini(jj), &
	264	& krefdate = irefdate(jj) )
	265	CALL ddatetoymdhms( ddobsend, iyea, imon, iday, ihou, imin, isec )
	266	CALL greg2jul( isec, imin, ihou, iday, imon, iyea, djulend(jj), &
	267	& krefdate = irefdate(jj) )
	268
[4990]	269	ioserrcount=0
[7992]	270	IF ( lldavtimset ) THEN
	271
	272	IF ( ANY ( idailyavtypes(:) /= -1 ) .AND. lwp) THEN
	273	WRITE(numout,*)' Resetting time of daily averaged', &
	274	& ' observations to the end of the day'
	275	ENDIF
	276
[2128]	277	DO ji = 1, inpfiles(jj)%nobs
[4990]	278	READ( inpfiles(jj)%cdtyp(ji), '(I4)', IOSTAT = ios, ERR = 900 ) itype
	279	900 IF ( ios /= 0 ) THEN
[7992]	280	! Set type to zero if there is a problem in the string conversion
	281	itype = 0
[4990]	282	ENDIF
[7992]	283
	284	IF ( ANY ( idailyavtypes(:) == itype ) ) THEN
	285	! for daily averaged data force the time
	286	! to be the last time-step of the day, but still within the day.
	287	IF ( inpfiles(jj)%ptim(ji) >= 0. ) THEN
	288	inpfiles(jj)%ptim(ji) = &
	289	& INT(inpfiles(jj)%ptim(ji)) + 0.9999
	290	ELSE
	291	inpfiles(jj)%ptim(ji) = &
	292	& INT(inpfiles(jj)%ptim(ji)) - 0.0001
	293	ENDIF
[2128]	294	ENDIF
[7992]	295
[2128]	296	END DO
[7992]	297
[2128]	298	ENDIF
[7992]	299
[2128]	300	IF ( inpfiles(jj)%nobs > 0 ) THEN
[7992]	301	inpfiles(jj)%iproc(:,:) = -1
	302	inpfiles(jj)%iobsi(:,:) = -1
	303	inpfiles(jj)%iobsj(:,:) = -1
[2128]	304	ENDIF
	305	inowin = 0
	306	DO ji = 1, inpfiles(jj)%nobs
[7992]	307	IF ( BTEST(inpfiles(jj)%ioqc(ji),2 ) ) CYCLE
[9306]	308	llcycle = .TRUE.
	309	DO jvar = 1, kvars
	310	IF ( .NOT. ( BTEST(inpfiles(jj)%ivqc(ji,jvar),2) ) ) THEN
	311	llcycle = .FALSE.
	312	EXIT
	313	ENDIF
	314	END DO
	315	IF ( llcycle ) CYCLE
[2128]	316	IF ( ( inpfiles(jj)%ptim(ji) > djulini(jj) ) .AND. &
	317	& ( inpfiles(jj)%ptim(ji) <= djulend(jj) ) ) THEN
	318	inowin = inowin + 1
	319	ENDIF
	320	END DO
	321	ALLOCATE( zlam(inowin) )
	322	ALLOCATE( zphi(inowin) )
[9306]	323	ALLOCATE( iobsi(inowin,kvars) )
	324	ALLOCATE( iobsj(inowin,kvars) )
	325	ALLOCATE( iproc(inowin,kvars) )
[2128]	326	inowin = 0
	327	DO ji = 1, inpfiles(jj)%nobs
[7992]	328	IF ( BTEST(inpfiles(jj)%ioqc(ji),2 ) ) CYCLE
[9306]	329	llcycle = .TRUE.
	330	DO jvar = 1, kvars
	331	IF ( .NOT. ( BTEST(inpfiles(jj)%ivqc(ji,jvar),2) ) ) THEN
	332	llcycle = .FALSE.
	333	EXIT
	334	ENDIF
	335	END DO
	336	IF ( llcycle ) CYCLE
[2128]	337	IF ( ( inpfiles(jj)%ptim(ji) > djulini(jj) ) .AND. &
	338	& ( inpfiles(jj)%ptim(ji) <= djulend(jj) ) ) THEN
	339	inowin = inowin + 1
	340	zlam(inowin) = inpfiles(jj)%plam(ji)
	341	zphi(inowin) = inpfiles(jj)%pphi(ji)
	342	ENDIF
	343	END DO
	344
[9306]	345	! Assume anything other than velocity is on T grid
	346	IF ( TRIM( inpfiles(jj)%cname(1) ) == 'UVEL' ) THEN
	347	CALL obs_grid_search( inowin, zlam, zphi, iobsi(:,1), iobsj(:,1), &
	348	& iproc(:,1), 'U' )
	349	CALL obs_grid_search( inowin, zlam, zphi, iobsi(:,2), iobsj(:,2), &
	350	& iproc(:,2), 'V' )
	351	ELSE
	352	CALL obs_grid_search( inowin, zlam, zphi, iobsi(:,1), iobsj(:,1), &
	353	& iproc(:,1), 'T' )
	354	IF ( kvars > 1 ) THEN
	355	DO jvar = 2, kvars
	356	iobsi(:,jvar) = iobsi(:,1)
	357	iobsj(:,jvar) = iobsj(:,1)
	358	iproc(:,jvar) = iproc(:,1)
	359	END DO
	360	ENDIF
[7992]	361	ENDIF
[2128]	362
	363	inowin = 0
	364	DO ji = 1, inpfiles(jj)%nobs
[7992]	365	IF ( BTEST(inpfiles(jj)%ioqc(ji),2 ) ) CYCLE
[9306]	366	llcycle = .TRUE.
	367	DO jvar = 1, kvars
	368	IF ( .NOT. ( BTEST(inpfiles(jj)%ivqc(ji,jvar),2) ) ) THEN
	369	llcycle = .FALSE.
	370	EXIT
	371	ENDIF
	372	END DO
	373	IF ( llcycle ) CYCLE
[2128]	374	IF ( ( inpfiles(jj)%ptim(ji) > djulini(jj) ) .AND. &
	375	& ( inpfiles(jj)%ptim(ji) <= djulend(jj) ) ) THEN
	376	inowin = inowin + 1
[9306]	377	DO jvar = 1, kvars
	378	inpfiles(jj)%iproc(ji,jvar) = iproc(inowin,jvar)
	379	inpfiles(jj)%iobsi(ji,jvar) = iobsi(inowin,jvar)
	380	inpfiles(jj)%iobsj(ji,jvar) = iobsj(inowin,jvar)
	381	END DO
	382	IF ( kvars > 1 ) THEN
	383	DO jvar = 2, kvars
	384	IF ( inpfiles(jj)%iproc(ji,jvar) /= &
	385	& inpfiles(jj)%iproc(ji,1) ) THEN
	386	CALL ctl_stop( 'Error in obs_read_prof:', &
	387	& 'observation on different processors for different vars')
	388	ENDIF
	389	END DO
[7992]	390	ENDIF
[2128]	391	ENDIF
	392	END DO
[9306]	393	DEALLOCATE( zlam, zphi, iobsi, iobsj, iproc )
[2128]	394
	395	DO ji = 1, inpfiles(jj)%nobs
[7992]	396	IF ( BTEST(inpfiles(jj)%ioqc(ji),2 ) ) CYCLE
[9306]	397	llcycle = .TRUE.
	398	DO jvar = 1, kvars
	399	IF ( .NOT. ( BTEST(inpfiles(jj)%ivqc(ji,jvar),2) ) ) THEN
	400	llcycle = .FALSE.
	401	EXIT
	402	ENDIF
	403	END DO
	404	IF ( llcycle ) CYCLE
[2128]	405	IF ( ( inpfiles(jj)%ptim(ji) > djulini(jj) ) .AND. &
	406	& ( inpfiles(jj)%ptim(ji) <= djulend(jj) ) ) THEN
	407	IF ( nproc == 0 ) THEN
	408	IF ( inpfiles(jj)%iproc(ji,1) > nproc ) CYCLE
	409	ELSE
	410	IF ( inpfiles(jj)%iproc(ji,1) /= nproc ) CYCLE
	411	ENDIF
	412	llvalprof = .FALSE.
[9306]	413	DO jvar = 1, kvars
	414	IF ( ldvar(jvar) ) THEN
	415	DO ij = 1,inpfiles(jj)%nlev
	416	IF ( inpfiles(jj)%pdep(ij,ji) >= 6000. ) &
	417	& CYCLE
	418	IF ( .NOT. BTEST(inpfiles(jj)%ivlqc(ij,ji,jvar),2) .AND. &
	419	& .NOT. BTEST(inpfiles(jj)%idqc(ij,ji),2) ) THEN
	420	ivart0(jvar) = ivart0(jvar) + 1
	421	ENDIF
	422	END DO
	423	ENDIF
	424	END DO
	425	DO ij = 1,inpfiles(jj)%nlev
[2128]	426	IF ( inpfiles(jj)%pdep(ij,ji) >= 6000. ) &
	427	& CYCLE
[9306]	428	DO jvar = 1, kvars
	429	IF ( ( .NOT. BTEST(inpfiles(jj)%ivlqc(ij,ji,jvar),2) .AND. &
	430	& .NOT. BTEST(inpfiles(jj)%idqc(ij,ji),2) .AND. &
	431	& ldvar(jvar) ) ) THEN
	432	ip3dt = ip3dt + 1
	433	llvalprof = .TRUE.
	434	EXIT
	435	ENDIF
	436	END DO
	437	END DO
[2128]	438
	439	IF ( llvalprof ) iprof = iprof + 1
	440
	441	ENDIF
	442	END DO
	443
	444	ENDIF
	445
	446	END DO prof_files
[7992]	447
[2128]	448	!-----------------------------------------------------------------------
	449	! Get the time ordered indices of the input data
	450	!-----------------------------------------------------------------------
	451
	452	!---------------------------------------------------------------------
	453	! Loop over input data files to count total number of profiles
	454	!---------------------------------------------------------------------
	455	iproftot = 0
	456	DO jj = 1, inobf
	457	DO ji = 1, inpfiles(jj)%nobs
[7992]	458	IF ( BTEST(inpfiles(jj)%ioqc(ji),2 ) ) CYCLE
[9306]	459	llcycle = .TRUE.
	460	DO jvar = 1, kvars
	461	IF ( .NOT. ( BTEST(inpfiles(jj)%ivqc(ji,jvar),2) ) ) THEN
	462	llcycle = .FALSE.
	463	EXIT
	464	ENDIF
	465	END DO
	466	IF ( llcycle ) CYCLE
[2128]	467	IF ( ( inpfiles(jj)%ptim(ji) > djulini(jj) ) .AND. &
	468	& ( inpfiles(jj)%ptim(ji) <= djulend(jj) ) ) THEN
	469	iproftot = iproftot + 1
	470	ENDIF
	471	END DO
	472	END DO
	473
	474	ALLOCATE( iindx(iproftot), ifileidx(iproftot), &
	475	& iprofidx(iproftot), zdat(iproftot) )
	476	jk = 0
	477	DO jj = 1, inobf
	478	DO ji = 1, inpfiles(jj)%nobs
[7992]	479	IF ( BTEST(inpfiles(jj)%ioqc(ji),2 ) ) CYCLE
[9306]	480	llcycle = .TRUE.
	481	DO jvar = 1, kvars
	482	IF ( .NOT. ( BTEST(inpfiles(jj)%ivqc(ji,jvar),2) ) ) THEN
	483	llcycle = .FALSE.
	484	EXIT
	485	ENDIF
	486	END DO
	487	IF ( llcycle ) CYCLE
[2128]	488	IF ( ( inpfiles(jj)%ptim(ji) > djulini(jj) ) .AND. &
	489	& ( inpfiles(jj)%ptim(ji) <= djulend(jj) ) ) THEN
	490	jk = jk + 1
	491	ifileidx(jk) = jj
	492	iprofidx(jk) = ji
	493	zdat(jk) = inpfiles(jj)%ptim(ji)
	494	ENDIF
	495	END DO
	496	END DO
	497	CALL sort_dp_indx( iproftot, &
	498	& zdat, &
	499	& iindx )
[7992]	500
[2128]	501	iv3dt(:) = -1
	502	IF (ldsatt) THEN
[9306]	503	iv3dt(:) = ip3dt
[2128]	504	ELSE
[9306]	505	iv3dt(:) = ivart0(:)
[2128]	506	ENDIF
	507	CALL obs_prof_alloc( profdata, kvars, kextr, iprof, iv3dt, &
[11468]	508	& kstp, jpi, jpj, jpk, ldclim )
[7992]	509
[2128]	510	! * Read obs/positions, QC, all variable and assign to profdata
	511
	512	profdata%nprof = 0
	513	profdata%nvprot(:) = 0
[11863]	514	profdata%cvars(:) = clvarsin(:)
[2128]	515	iprof = 0
	516
	517	ip3dt = 0
[9306]	518	ivart(:) = 0
	519	itypvar (:,:) = 0
	520	itypvarmpp(:,:) = 0
[7992]	521
	522	ioserrcount = 0
[2128]	523	DO jk = 1, iproftot
[7992]	524
[2128]	525	jj = ifileidx(iindx(jk))
	526	ji = iprofidx(iindx(jk))
	527
[9306]	528	IF ( BTEST(inpfiles(jj)%ioqc(ji),2 ) ) CYCLE
	529	llcycle = .TRUE.
	530	DO jvar = 1, kvars
	531	IF ( .NOT. ( BTEST(inpfiles(jj)%ivqc(ji,jvar),2) ) ) THEN
	532	llcycle = .FALSE.
	533	EXIT
	534	ENDIF
	535	END DO
	536	IF ( llcycle ) CYCLE
[2128]	537
	538	IF ( ( inpfiles(jj)%ptim(ji) > djulini(jj) ) .AND. &
	539	& ( inpfiles(jj)%ptim(ji) <= djulend(jj) ) ) THEN
[7992]	540
[2128]	541	IF ( nproc == 0 ) THEN
	542	IF ( inpfiles(jj)%iproc(ji,1) > nproc ) CYCLE
	543	ELSE
	544	IF ( inpfiles(jj)%iproc(ji,1) /= nproc ) CYCLE
	545	ENDIF
[7992]	546
[2128]	547	llvalprof = .FALSE.
	548
	549	IF ( inpfiles(jj)%ioqc(ji) > 2 ) CYCLE
	550
[7992]	551	IF ( BTEST(inpfiles(jj)%ioqc(ji),2 ) ) CYCLE
[9306]	552	llcycle = .TRUE.
	553	DO jvar = 1, kvars
	554	IF ( .NOT. ( BTEST(inpfiles(jj)%ivqc(ji,jvar),2) ) ) THEN
	555	llcycle = .FALSE.
	556	EXIT
	557	ENDIF
	558	END DO
	559	IF ( llcycle ) CYCLE
[2128]	560
	561	loop_prof : DO ij = 1, inpfiles(jj)%nlev
[7992]	562
[2128]	563	IF ( inpfiles(jj)%pdep(ij,ji) >= 6000. ) &
	564	& CYCLE
[7992]	565
[9306]	566	DO jvar = 1, kvars
	567	IF ( .NOT. BTEST(inpfiles(jj)%ivlqc(ij,ji,jvar),2) .AND. &
	568	& .NOT. BTEST(inpfiles(jj)%idqc(ij,ji),2) ) THEN
[7992]	569
[9306]	570	llvalprof = .TRUE.
	571	EXIT loop_prof
[7992]	572
[9306]	573	ENDIF
	574	END DO
[7992]	575
[2128]	576	END DO loop_prof
[7992]	577
[2128]	578	! Set profile information
[7992]	579
[2128]	580	IF ( llvalprof ) THEN
[7992]	581
[2128]	582	iprof = iprof + 1
	583
	584	CALL jul2greg( isec, &
	585	& imin, &
	586	& ihou, &
	587	& iday, &
	588	& imon, &
	589	& iyea, &
	590	& inpfiles(jj)%ptim(ji), &
	591	& irefdate(jj) )
	592
	593
	594	! Profile time coordinates
	595	profdata%nyea(iprof) = iyea
	596	profdata%nmon(iprof) = imon
	597	profdata%nday(iprof) = iday
	598	profdata%nhou(iprof) = ihou
	599	profdata%nmin(iprof) = imin
[7992]	600
[2128]	601	! Profile space coordinates
	602	profdata%rlam(iprof) = inpfiles(jj)%plam(ji)
	603	profdata%rphi(iprof) = inpfiles(jj)%pphi(ji)
	604
	605	! Coordinate search parameters
[9306]	606	DO jvar = 1, kvars
	607	profdata%mi (iprof,jvar) = inpfiles(jj)%iobsi(ji,jvar)
	608	profdata%mj (iprof,jvar) = inpfiles(jj)%iobsj(ji,jvar)
	609	END DO
[7992]	610
[2128]	611	! Profile WMO number
	612	profdata%cwmo(iprof) = inpfiles(jj)%cdwmo(ji)
[7992]	613
[2128]	614	! Instrument type
[4990]	615	READ( inpfiles(jj)%cdtyp(ji), '(I4)', IOSTAT = ios, ERR = 901 ) itype
	616	901 IF ( ios /= 0 ) THEN
	617	IF (ioserrcount == 0) CALL ctl_warn ( 'Problem converting an instrument type to integer. Setting type to zero' )
	618	ioserrcount = ioserrcount + 1
	619	itype = 0
	620	ENDIF
[7992]	621
[2128]	622	profdata%ntyp(iprof) = itype
[7992]	623
[2128]	624	! QC stuff
	625
	626	profdata%nqc(iprof) = inpfiles(jj)%ioqc(ji)
	627	profdata%nqcf(:,iprof) = inpfiles(jj)%ioqcf(:,ji)
	628	profdata%ipqc(iprof) = inpfiles(jj)%ipqc(ji)
	629	profdata%ipqcf(:,iprof) = inpfiles(jj)%ipqcf(:,ji)
	630	profdata%itqc(iprof) = inpfiles(jj)%itqc(ji)
	631	profdata%itqcf(:,iprof) = inpfiles(jj)%itqcf(:,ji)
	632	profdata%ivqc(iprof,:) = inpfiles(jj)%ivqc(ji,:)
	633	profdata%ivqcf(:,iprof,:) = inpfiles(jj)%ivqcf(:,ji,:)
	634
	635	! Bookkeeping data to match profiles
	636	profdata%npidx(iprof) = iprof
	637	profdata%npfil(iprof) = iindx(jk)
	638
	639	! Observation QC flag (whole profile)
	640	profdata%nqc(iprof) = 0 !TODO
	641
[7992]	642	loop_p : DO ij = 1, inpfiles(jj)%nlev
	643
[2128]	644	IF ( inpfiles(jj)%pdep(ij,ji) >= 6000. ) &
	645	& CYCLE
	646
	647	IF (ldsatt) THEN
	648
[9306]	649	DO jvar = 1, kvars
	650	IF ( ( .NOT. BTEST(inpfiles(jj)%ivlqc(ij,ji,jvar),2) .AND. &
	651	& .NOT. BTEST(inpfiles(jj)%idqc(ij,ji),2) .AND. &
	652	& ldvar(jvar) ) ) THEN
	653	ip3dt = ip3dt + 1
	654	EXIT
	655	ELSE IF ( jvar == kvars ) THEN
	656	CYCLE loop_p
	657	ENDIF
	658	END DO
[7992]	659
[2128]	660	ENDIF
	661
[9306]	662	DO jvar = 1, kvars
	663
	664	IF ( ( .NOT. BTEST(inpfiles(jj)%ivlqc(ij,ji,jvar),2) .AND. &
	665	& .NOT. BTEST(inpfiles(jj)%idqc(ij,ji),2) .AND. &
	666	& ldvar(jvar) ) .OR. ldsatt ) THEN
[7992]	667
[9306]	668	IF (ldsatt) THEN
[2128]	669
[9306]	670	ivart(jvar) = ip3dt
[2128]	671
[9306]	672	ELSE
[2128]	673
[9306]	674	ivart(jvar) = ivart(jvar) + 1
[7992]	675
[9306]	676	ENDIF
[2128]	677
[9306]	678	! Depth of jvar observation
	679	profdata%var(jvar)%vdep(ivart(jvar)) = &
	680	& inpfiles(jj)%pdep(ij,ji)
[7992]	681
[9306]	682	! Depth of jvar observation QC
	683	profdata%var(jvar)%idqc(ivart(jvar)) = &
	684	& inpfiles(jj)%idqc(ij,ji)
[7992]	685
[9306]	686	! Depth of jvar observation QC flags
	687	profdata%var(jvar)%idqcf(:,ivart(jvar)) = &
	688	& inpfiles(jj)%idqcf(:,ij,ji)
[7992]	689
[9306]	690	! Profile index
	691	profdata%var(jvar)%nvpidx(ivart(jvar)) = iprof
[7992]	692
[9306]	693	! Vertical index in original profile
	694	profdata%var(jvar)%nvlidx(ivart(jvar)) = ij
[2128]	695
[9306]	696	! Profile jvar value
	697	IF ( .NOT. BTEST(inpfiles(jj)%ivlqc(ij,ji,jvar),2) .AND. &
	698	& .NOT. BTEST(inpfiles(jj)%idqc(ij,ji),2) ) THEN
	699	profdata%var(jvar)%vobs(ivart(jvar)) = &
	700	& inpfiles(jj)%pob(ij,ji,jvar)
	701	IF ( ldmod ) THEN
	702	profdata%var(jvar)%vmod(ivart(jvar)) = &
	703	& inpfiles(jj)%padd(ij,ji,1,jvar)
	704	ENDIF
[11468]	705	IF ( profdata%lclim ) THEN
	706	profdata%var(jvar)%vclm(ivart(jvar)) = fbrmdi
	707	ENDIF
[9306]	708	! Count number of profile var1 data as function of type
	709	itypvar( profdata%ntyp(iprof) + 1, jvar ) = &
	710	& itypvar( profdata%ntyp(iprof) + 1, jvar ) + 1
	711	ELSE
	712	profdata%var(jvar)%vobs(ivart(jvar)) = fbrmdi
[2128]	713	ENDIF
	714
[9306]	715	! Profile jvar qc
	716	profdata%var(jvar)%nvqc(ivart(jvar)) = &
	717	& inpfiles(jj)%ivlqc(ij,ji,jvar)
[2128]	718
[9306]	719	! Profile jvar qc flags
	720	profdata%var(jvar)%nvqcf(:,ivart(jvar)) = &
	721	& inpfiles(jj)%ivlqcf(:,ij,ji,jvar)
[2128]	722
[9306]	723	! Profile insitu T value
	724	IF ( TRIM( inpfiles(jj)%cname(jvar) ) == 'POTM' ) THEN
	725	profdata%var(jvar)%vext(ivart(jvar),1) = &
	726	& inpfiles(jj)%pext(ij,ji,1)
	727	ENDIF
[7992]	728
[2128]	729	ENDIF
[9306]	730
	731	END DO
[2128]	732
	733	END DO loop_p
	734
	735	ENDIF
	736
	737	ENDIF
	738
	739	END DO
	740
	741	!-----------------------------------------------------------------------
	742	! Sum up over processors
	743	!-----------------------------------------------------------------------
[7992]	744
[9306]	745	DO jvar = 1, kvars
	746	CALL obs_mpp_sum_integer ( ivart0(jvar), ivartmpp(jvar) )
	747	END DO
[7992]	748	CALL obs_mpp_sum_integer ( ip3dt, ip3dtmpp )
	749
[9306]	750	DO jvar = 1, kvars
	751	CALL obs_mpp_sum_integers( itypvar(:,jvar), itypvarmpp(:,jvar), ntyp1770 + 1 )
	752	END DO
[7992]	753
[2128]	754	!-----------------------------------------------------------------------
	755	! Output number of observations.
	756	!-----------------------------------------------------------------------
	757	IF(lwp) THEN
	758	WRITE(numout,*)
[7992]	759	WRITE(numout,'(A)') ' Profile data'
[2128]	760	WRITE(numout,'(1X,A)') '------------'
	761	WRITE(numout,*)
[9306]	762	DO jvar = 1, kvars
	763	WRITE(numout,'(1X,A)') 'Profile data, '//TRIM( profdata%cvars(jvar) )
	764	WRITE(numout,'(1X,A)') '------------------------'
	765	DO ji = 0, ntyp1770
	766	IF ( itypvarmpp(ji+1,jvar) > 0 ) THEN
	767	WRITE(numout,'(1X,A3,1X,A48,A3,I8)') ctypshort(ji), &
	768	& cwmonam1770(ji)(1:52),' = ', &
	769	& itypvarmpp(ji+1,jvar)
	770	ENDIF
	771	END DO
	772	WRITE(numout,'(1X,A)') &
	773	& '---------------------------------------------------------------'
	774	WRITE(numout,'(1X,A55,I8)') &
	775	& 'Total profile data for variable '//TRIM( profdata%cvars(jvar) )// &
	776	& ' = ', ivartmpp(jvar)
	777	WRITE(numout,'(1X,A)') &
	778	& '---------------------------------------------------------------'
	779	WRITE(numout,*)
[2128]	780	END DO
	781	ENDIF
[7992]	782
[2128]	783	IF (ldsatt) THEN
[9306]	784	profdata%nvprot(:) = ip3dt
	785	profdata%nvprotmpp(:) = ip3dtmpp
[2128]	786	ELSE
[9306]	787	DO jvar = 1, kvars
	788	profdata%nvprot(jvar) = ivart(jvar)
	789	profdata%nvprotmpp(jvar) = ivartmpp(jvar)
	790	END DO
[2128]	791	ENDIF
	792	profdata%nprof = iprof
	793
	794	!-----------------------------------------------------------------------
	795	! Model level search
	796	!-----------------------------------------------------------------------
[9306]	797	DO jvar = 1, kvars
	798	IF ( ldvar(jvar) ) THEN
	799	CALL obs_level_search( jpk, gdept_1d, &
	800	& profdata%nvprot(jvar), profdata%var(jvar)%vdep, &
	801	& profdata%var(jvar)%mvk )
	802	ENDIF
	803	END DO
[7992]	804
[2128]	805	!-----------------------------------------------------------------------
	806	! Set model equivalent to 99999
	807	!-----------------------------------------------------------------------
	808	IF ( .NOT. ldmod ) THEN
	809	DO jvar = 1, kvars
	810	profdata%var(jvar)%vmod(:) = fbrmdi
	811	END DO
	812	ENDIF
	813	!-----------------------------------------------------------------------
	814	! Deallocate temporary data
	815	!-----------------------------------------------------------------------
[11863]	816	DEALLOCATE( ifileidx, iprofidx, zdat, clvarsin )
[2128]	817
	818	!-----------------------------------------------------------------------
	819	! Deallocate input data
	820	!-----------------------------------------------------------------------
	821	DO jj = 1, inobf
	822	CALL dealloc_obfbdata( inpfiles(jj) )
	823	END DO
	824	DEALLOCATE( inpfiles )
	825
[7992]	826	END SUBROUTINE obs_rea_prof
[2128]	827
	828	END MODULE obs_read_prof

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

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

Context Navigation

source: branches/UKMO/dev_r5518_obs_oper_update_SlaAssim/NEMOGCM/NEMO/OPA_SRC/OBS/obs_read_prof.F90 @ 11863

Download in other formats: