Context Navigation

obs_read_prof.F90 @ 2358

Last change on this file since 2358 was 2287, checked in by smasson, 14 years ago
update licence of all NEMO files...
Property svn:keywords set to `Id`
File size: 31.9 KB

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

Download in other formats: