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_grid.F90 in branches/dev_1784_OBS/NEMO/OPA_SRC/OBS – NEMO

Context Navigation

source: branches/dev_1784_OBS/NEMO/OPA_SRC/OBS/obs_grid.F90 @ 2001

Last change on this file since 2001 was 2001, checked in by djlea, 14 years ago
Adding observation operator code
File size: 46.8 KB

Line
1	MODULE obs_grid
2	!!======================================================================
3	!! * MODULE obs_grid *
4	!! Observation diagnostics: Various tools for grid searching etc.
5	!!======================================================================
6	!!----------------------------------------------------------------------
7	!! obs_grid_search : Find i,j on the ORCA grid from lat,lon
8	!! obs_level_search : Find level from depth
9	!! obs_zlevel_search : Find depth level from observed depth
10	!! obs_tlevel_search : Find temperature level from observed temp
11	!! obs_rlevel_search : Find density level from observed density
12	!!----------------------------------------------------------------------
13	!! * Modules used
14	USE par_kind, ONLY : & ! Precision variables
15	& wp
16	USE par_oce, ONLY : & ! Ocean parameters
17	& jpk, &
18	& jpni, &
19	& jpnj, &
20	& jpnij
21	USE dom_oce ! Ocean space and time domain variables
22	USE obs_mpp, ONLY : & ! MPP support routines for observation diagnostics
23	& obs_mpp_find_obs_proc, &
24	& mpp_global_max, &
25	& obs_mpp_max_integer
26	USE phycst, ONLY : & ! Physical constants
27	& rad
28	USE obs_utils, ONLY : & ! Observation operator utility functions
29	& grt_cir_dis, &
30	& chkerr
31	USE in_out_manager ! Printing support
32	USE netcdf
33	USE obs_const, ONLY : &
34	& obfillflt ! Fillvalue
35
36	IMPLICIT NONE
37
38	!! * Routine accessibility
39	PUBLIC obs_grid_setup, & ! Setup grid searching
40	& obs_grid_search, & ! Find i, j on the ORCA grid from lat, lon
41	& obs_grid_deallocate, & ! Deallocate the look up table
42	& obs_level_search ! Find level from depth
43
44	PRIVATE linquad, & ! Determine whether a point lies within a cell
45	& maxdist, & ! Find the maximum distance between 2 pts in a cell
46	& obs_grid_search_bruteforce, & ! Find i, j on the ORCA grid from lat, lon
47	& obs_grid_search_lookup ! Find i, j on the ORCA grid from lat, lon quicker
48
49	!!* Module variables
50
51	!! Default values
52	REAL, PUBLIC :: grid_search_res = 0.5 ! Resolution of grid
53	INTEGER, PRIVATE :: gsearch_nlons_def ! Num of longitudes
54	INTEGER, PRIVATE :: gsearch_nlats_def ! Num of latitudes
55	REAL(wp), PRIVATE :: gsearch_lonmin_def ! Min longitude
56	REAL(wp), PRIVATE :: gsearch_latmin_def ! Min latitude
57	REAL(wp), PRIVATE :: gsearch_dlon_def ! Lon spacing
58	REAL(wp), PRIVATE :: gsearch_dlat_def ! Lat spacing
59	!! Variable versions
60	INTEGER, PRIVATE :: nlons ! Num of longitudes
61	INTEGER, PRIVATE :: nlats ! Num of latitudes
62	REAL(wp), PRIVATE :: lonmin ! Min longitude
63	REAL(wp), PRIVATE :: latmin ! Min latitude
64	REAL(wp), PRIVATE :: dlon ! Lon spacing
65	REAL(wp), PRIVATE :: dlat ! Lat spacing
66
67	INTEGER, PRIVATE :: maxxdiff, maxydiff ! Max diffs between model points
68	INTEGER, PRIVATE :: limxdiff, limydiff
69
70	! Data storage
71	REAL(wp), PRIVATE, DIMENSION(:,:), ALLOCATABLE :: &
72	& lons, &
73	& lats
74	INTEGER, PRIVATE, DIMENSION(:,:), ALLOCATABLE :: &
75	& ixpos, &
76	& iypos, &
77	& iproc
78
79	! Switches
80	LOGICAL, PUBLIC :: &
81	& ln_grid_search_lookup ! Use lookup table to speed up grid search
82	LOGICAL, PUBLIC :: &
83	& ln_grid_global ! Use global distribution of observations
84	CHARACTER(LEN=44), PUBLIC :: &
85	& grid_search_file ! file name head for grid search lookup
86
87	CONTAINS
88
89	SUBROUTINE obs_grid_search( kobsin, plam, pphi, kobsi, kobsj, kproc, &
90	& cdgrid )
91	!!----------------------------------------------------------------------
92	!! * ROUTINE obs_grid_search *
93	!!
94	!! ** Purpose : Search local gridpoints to find the grid box containing
95	!! the observations calls either
96	!! obs_grid_search_bruteforce - the original brute force search
97	!! or
98	!! obs_grid_search_lookup - uses a lookup table to do a fast
99	!!search
100	!!History :
101	!! ! 2007-12 (D. Lea)
102	!!------------------------------------------------------------------------
103
104	!! * Arguments
105	INTEGER :: &
106	& kobsin ! Size of the observation arrays
107	REAL(KIND=wp), DIMENSION(kobsin), INTENT(IN) :: &
108	& plam, & ! Longitude of obsrvations
109	& pphi ! Latitude of observations
110	INTEGER, DIMENSION(kobsin), INTENT(OUT) :: &
111	& kobsi, & ! I-index of observations
112	& kobsj, & ! J-index of observations
113	& kproc ! Processor number of observations
114	CHARACTER(LEN=1) :: &
115	& cdgrid ! Grid to search
116
117	IF(kobsin > 0) THEN
118
119	IF ( ln_grid_search_lookup .AND. ( cdgrid == 'T' ) ) THEN
120	CALL obs_grid_search_lookup( kobsin, plam, pphi, &
121	& kobsi, kobsj, kproc )
122	ELSE
123	IF ( cdgrid == 'T' ) THEN
124	CALL obs_grid_search_bruteforce( jpi, jpj, jpiglo, jpjglo, &
125	& nldi, nlei,nldj, nlej, &
126	& nproc, jpnij, &
127	& glamt, gphit, tmask, &
128	& kobsin, plam, pphi, &
129	& kobsi, kobsj, kproc )
130	ELSEIF ( cdgrid == 'U' ) THEN
131	CALL obs_grid_search_bruteforce( jpi, jpj, jpiglo, jpjglo, &
132	& nldi, nlei,nldj, nlej, &
133	& nproc, jpnij, &
134	& glamu, gphiu, umask, &
135	& kobsin, plam, pphi, &
136	& kobsi, kobsj, kproc )
137	ELSEIF ( cdgrid == 'V' ) THEN
138	CALL obs_grid_search_bruteforce( jpi, jpj, jpiglo, jpjglo, &
139	& nldi, nlei,nldj, nlej, &
140	& nproc, jpnij, &
141	& glamv, gphiv, vmask, &
142	& kobsin, plam, pphi, &
143	& kobsi, kobsj, kproc )
144	ELSEIF ( cdgrid == 'F' ) THEN
145	CALL obs_grid_search_bruteforce( jpi, jpj, jpiglo, jpjglo, &
146	& nldi, nlei,nldj, nlej, &
147	& nproc, jpnij, &
148	& glamf, gphif, fmask, &
149	& kobsin, plam, pphi, &
150	& kobsi, kobsj, kproc )
151	ELSE
152	CALL ctl_stop( 'Grid not supported' )
153	ENDIF
154	ENDIF
155
156	ENDIF
157
158	END SUBROUTINE obs_grid_search
159
160	#include "obs_grid_search_bruteforce.h90"
161
162	SUBROUTINE obs_grid_search_lookup( kobs, plam, pphi, kobsi, kobsj, kproc )
163	!!----------------------------------------------------------------------
164	!! * ROUTINE obs_grid_search_lookup *
165	!!
166	!! ** Purpose : Search local gridpoints to find the grid box containing
167	!! the observations (much faster then obs_grid_search_bruteforce)
168	!!
169	!! ** Method : Call to linquad
170	!!
171	!! ** Action : Return kproc holding the observation and kiobsi,kobsj
172	!! valid on kproc=nproc processor only.
173	!!
174	!! History :
175	!! ! 2007-12 (D. Lea) new routine based on obs_grid_search
176	!!! updated with fixes from new version of obs_grid_search_bruteforce
177	!!! speeded up where points are not near a "difficult" region like an edge
178	!!----------------------------------------------------------------------
179
180	!! * Arguments
181	INTEGER :: &
182	& kobs ! Size of the observation arrays
183	REAL(KIND=wp), DIMENSION(kobs), INTENT(IN) :: &
184	& plam, & ! Longitude of obsrvations
185	& pphi ! Latitude of observations
186	INTEGER, DIMENSION(kobs), INTENT(OUT) :: &
187	& kobsi, & ! I-index of observations
188	& kobsj, & ! J-index of observations
189	& kproc ! Processor number of observations
190
191	!! * Local declarations
192	REAL(KIND=wp), DIMENSION(:), ALLOCATABLE :: &
193	& zplam
194	REAL(wp) :: &
195	& zlammax, &
196	& zlam
197	INTEGER :: &
198	& ji, &
199	& jj, &
200	& jk, &
201	& jo, &
202	& isx, &
203	& isy, &
204	& jimin, &
205	& jimax, &
206	& jjmin, &
207	& jjmax, &
208	& jojimin, &
209	& jojimax, &
210	& jojjmin, &
211	& jojjmax, &
212	& ipx1, &
213	& ipy1, &
214	& ip, &
215	& jp, &
216	& ipx, &
217	& ipy, &
218	& ipmx, &
219	& jlon, &
220	& jlat, &
221	& joffset, &
222	& jostride
223	REAL(KIND=wp), DIMENSION(:,:), ALLOCATABLE :: &
224	& zlamg, &
225	& zphig, &
226	& zmskg, &
227	& zphitmax,&
228	& zphitmin,&
229	& zlamtmax,&
230	& zlamtmin
231	LOGICAL, DIMENSION(:,:), ALLOCATABLE :: &
232	& llinvalidcell
233	REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: &
234	& zlamtm, &
235	& zphitm
236	LOGICAL :: &
237	& llfourflag
238	INTEGER :: &
239	& ifourflagcountt, &
240	& ifourflagcountf
241	INTEGER, DIMENSION(5) :: &
242	& ifourflagcountr
243
244	!-----------------------------------------------------------------------
245	! Define grid for grid search
246	!-----------------------------------------------------------------------
247	IF (ln_grid_global) THEN
248	jlon = jpiglo
249	jlat = jpjglo
250	joffset = nproc
251	jostride = jpnij
252	ELSE
253	jlon = jpi
254	jlat = jpj
255	joffset = 0
256	jostride = 1
257	ENDIF
258	!-----------------------------------------------------------------------
259	! Set up data for grid search
260	!-----------------------------------------------------------------------
261	ALLOCATE( &
262	& zlamg(jlon,jlat), &
263	& zphig(jlon,jlat), &
264	& zmskg(jlon,jlat), &
265	& zphitmax(jlon-1,jlat-1), &
266	& zphitmin(jlon-1,jlat-1), &
267	& zlamtmax(jlon-1,jlat-1), &
268	& zlamtmin(jlon-1,jlat-1), &
269	& llinvalidcell(jlon-1,jlat-1), &
270	& zlamtm(4,jlon-1,jlat-1), &
271	& zphitm(4,jlon-1,jlat-1) &
272	& )
273	!-----------------------------------------------------------------------
274	! Copy data to local arrays
275	!-----------------------------------------------------------------------
276	IF (ln_grid_global) THEN
277	zlamg(:,:) = -1.e+10
278	zphig(:,:) = -1.e+10
279	zmskg(:,:) = -1.e+10
280	! Add various grids here.
281	DO jj = nldj, nlej
282	DO ji = nldi, nlei
283	zlamg(mig(ji),mjg(jj)) = glamt(ji,jj)
284	zphig(mig(ji),mjg(jj)) = gphit(ji,jj)
285	zmskg(mig(ji),mjg(jj)) = tmask(ji,jj,1)
286	END DO
287	END DO
288	CALL mpp_global_max( zlamg )
289	CALL mpp_global_max( zphig )
290	CALL mpp_global_max( zmskg )
291	ELSE
292	! Add various grids here.
293	DO jj = 1, jlat
294	DO ji = 1, jlon
295	zlamg(ji,jj) = glamt(ji,jj)
296	zphig(ji,jj) = gphit(ji,jj)
297	zmskg(ji,jj) = tmask(ji,jj,1)
298	END DO
299	END DO
300	ENDIF
301	!-----------------------------------------------------------------------
302	! Copy longitudes
303	!-----------------------------------------------------------------------
304	ALLOCATE( &
305	& zplam(kobs) &
306	& )
307	DO jo = 1, kobs
308	zplam(jo) = plam(jo)
309	ENDDO
310	!-----------------------------------------------------------------------
311	! Set default values for output
312	!-----------------------------------------------------------------------
313	kproc(:) = -1
314	kobsi(:) = -1
315	kobsj(:) = -1
316	!-----------------------------------------------------------------------
317	! Copy grid positions to temporary arrays and renormalize to 0 to 360.
318	!-----------------------------------------------------------------------
319	DO jj = 1, jlat-1
320	DO ji = 1, jlon-1
321	zlamtm(1,ji,jj) = zlamg(ji ,jj )
322	zphitm(1,ji,jj) = zphig(ji ,jj )
323	zlamtm(2,ji,jj) = zlamg(ji+1,jj )
324	zphitm(2,ji,jj) = zphig(ji+1,jj )
325	zlamtm(3,ji,jj) = zlamg(ji+1,jj+1)
326	zphitm(3,ji,jj) = zphig(ji+1,jj+1)
327	zlamtm(4,ji,jj) = zlamg(ji ,jj+1)
328	zphitm(4,ji,jj) = zphig(ji ,jj+1)
329	END DO
330	END DO
331	WHERE ( zlamtm(:,:,:) < 0.0_wp )
332	zlamtm(:,:,:) = zlamtm(:,:,:) + 360.0_wp
333	END WHERE
334	WHERE ( zlamtm(:,:,:) > 360.0_wp )
335	zlamtm(:,:,:) = zlamtm(:,:,:) - 360.0_wp
336	END WHERE
337	!-----------------------------------------------------------------------
338	! Handle case of the wraparound; beware, not working with orca180
339	!-----------------------------------------------------------------------
340	DO jj = 1, jlat-1
341	DO ji = 1, jlon-1
342	zlammax = MAXVAL( zlamtm(:,ji,jj) )
343	WHERE (zlammax - zlamtm(:, ji, jj) > 180 ) &
344	& zlamtm(:,ji,jj) = zlamtm(:,ji,jj) + 360._wp
345	zphitmax(ji,jj) = MAXVAL(zphitm(:,ji,jj))
346	zphitmin(ji,jj) = MINVAL(zphitm(:,ji,jj))
347	zlamtmax(ji,jj) = MAXVAL(zlamtm(:,ji,jj))
348	zlamtmin(ji,jj) = MINVAL(zlamtm(:,ji,jj))
349	END DO
350	END DO
351	!-----------------------------------------------------------------------
352	! Search for boxes with only land points mark them invalid
353	!-----------------------------------------------------------------------
354	llinvalidcell(:,:) = .FALSE.
355	DO jj = 1, jlat-1
356	DO ji = 1, jlon-1
357	llinvalidcell(ji,jj) = &
358	& zmskg(ji ,jj ) == 0.0_wp .AND. &
359	& zmskg(ji+1,jj ) == 0.0_wp .AND. &
360	& zmskg(ji+1,jj+1) == 0.0_wp .AND. &
361	& zmskg(ji ,jj+1) == 0.0_wp
362	END DO
363	END DO
364
365	if(lwp) WRITE(numout,*) 'obs_grid_search_lookup do coordinate search using lookup table'
366
367	!-----------------------------------------------------------------------
368	! Do coordinate search using lookup table with local searches.
369	! - For land points kproc is set to number of the processor + 1000000
370	! and we continue the search.
371	! - For ocean points kproc is set to the number of the processor
372	! and we stop the search.
373	!-----------------------------------------------------------------------
374	ifourflagcountt = 0
375	ifourflagcountf = 0
376	ifourflagcountr(:) = 0
377
378	!------------------------------------------------------------------------
379	! Master loop for grid search
380	!------------------------------------------------------------------------
381
382	gpkobs: DO jo = 1+joffset, kobs, jostride
383	! Normal case
384	! specify 4 points which surround the lat lon of interest
385	! x i,j+1 x i+1, j+1
386	!
387	!
388	! * lon,lat
389	! x i,j x i+1,j
390
391	! bottom corner point
392	ipx1 = INT( ( zplam(jo) - lonmin ) / dlon + 1.0 )
393	ipy1 = INT( ( pphi (jo) - latmin ) / dlat + 1.0 )
394
395	ipx = ipx1 + 1
396	ipy = ipy1 + 1
397
398	! flag for searching around four points separately
399	! default to false
400	llfourflag = .FALSE.
401
402	! check for point fully outside of region
403	IF ( (ipx1 > nlons) .OR. (ipy1 > nlats) .OR. &
404	& (ipx < 1) .OR. (ipy < 1) ) THEN
405	CYCLE
406	ENDIF
407	! check wrap around
408	IF ( (ipx > nlons) .OR. (ipy > nlats) .OR. &
409	& (ipx1 < 1) .OR. (ipy1 < 1) ) THEN
410	llfourflag=.TRUE.
411	ifourflagcountr(1) = ifourflagcountr(1) + 1
412	ENDIF
413
414	IF (MAXVAL(ixpos(ipx1:ipx,ipy1:ipy)) == -1) CYCLE! cycle if no lookup points found
415
416	jimin = 0
417	jimax = 0
418	jjmin = 0
419	jjmax = 0
420
421	IF (.NOT. llfourflag) THEN
422
423	! calculate points range
424	! define a square region encompassing the four corner points
425	! do I need the -1 points?
426
427	jojimin = MINVAL(ixpos(ipx1:ipx,ipy1:ipy)) - 1
428	jojimax = MAXVAL(ixpos(ipx1:ipx,ipy1:ipy)) + 1
429	jojjmin = MINVAL(iypos(ipx1:ipx,ipy1:ipy)) - 1
430	jojjmax = MAXVAL(iypos(ipx1:ipx,ipy1:ipy)) + 1
431
432	jimin = jojimin - 1
433	jimax = jojimax + 1
434	jjmin = jojjmin - 1
435	jjmax = jojjmax + 1
436
437	IF ( jojimin < 0 .OR. jojjmin < 0) THEN
438	llfourflag = .TRUE.
439	ifourflagcountr(2) = ifourflagcountr(2) + 1
440	ENDIF
441	IF ( jojimax - jojimin > maxxdiff) THEN
442	llfourflag = .TRUE.
443	ifourflagcountr(3) = ifourflagcountr(3) + 1
444	ENDIF
445	IF ( jojjmax - jojjmin > maxydiff) THEN
446	llfourflag = .TRUE.
447	ifourflagcountr(4) = ifourflagcountr(4) + 1
448	ENDIF
449
450	ENDIF
451
452	ipmx = 0
453	IF (llfourflag) ipmx = 1
454
455	IF (llfourflag) THEN
456	ifourflagcountt = ifourflagcountt + 1
457	ELSE
458	ifourflagcountf = ifourflagcountf + 1
459	ENDIF
460
461	gridpointsn : DO ip = 0, ipmx
462	DO jp = 0, ipmx
463
464	IF ( kproc(jo) /= -1 ) EXIT gridpointsn
465
466	ipx = ipx1 + ip
467	ipy = ipy1 + jp
468
469	IF (llfourflag) THEN
470
471	! deal with wrap around
472	IF ( ipx > nlons ) ipx = 1
473	IF ( ipy > nlats ) ipy = 1
474	IF ( ipx < 1 ) ipx = nlons
475	IF ( ipy < 1 ) ipy = nlats
476
477	! get i,j
478	isx = ixpos(ipx,ipy)
479	isy = iypos(ipx,ipy)
480
481	! estimate appropriate search region (use max/min values)
482	jimin = isx - maxxdiff - 1
483	jimax = isx + maxxdiff + 1
484	jjmin = isy - maxydiff - 1
485	jjmax = isy + maxydiff + 1
486
487	ENDIF
488
489	IF ( jimin < 1 ) jimin = 1
490	IF ( jimax > jlon-1 ) jimax = jlon-1
491	IF ( jjmin < 1 ) jjmin = 1
492	IF ( jjmax > jlat-1 ) jjmax = jlat-1
493
494	!---------------------------------------------------------------
495	! Ensure that all observation longtiudes are between 0 and 360
496	!---------------------------------------------------------------
497
498	IF ( zplam(jo) < 0.0_wp ) zplam(jo) = zplam(jo) + 360.0_wp
499	IF ( zplam(jo) > 360.0_wp ) zplam(jo) = zplam(jo) - 360.0_wp
500
501	!---------------------------------------------------------------
502	! Find observations which are on within 1e-6 of a grid point
503	!---------------------------------------------------------------
504
505	gridloop: DO jj = jjmin, jjmax
506	DO ji = jimin, jimax
507	IF ( ABS( zphig(ji,jj) - pphi(jo) ) < 1e-6 ) THEN
508	zlam = zlamg(ji,jj)
509	IF ( zlam < 0.0_wp ) zlam = zlam + 360.0_wp
510	IF ( zlam > 360.0_wp ) zlam = zlam - 360.0_wp
511	IF ( ABS( zlam - zplam(jo) ) < 1e-6 ) THEN
512	IF ( llinvalidcell(ji,jj) ) THEN
513	kproc(jo) = nproc + 1000000
514	kobsi(jo) = ji + 1
515	kobsj(jo) = jj + 1
516	CYCLE
517	ELSE
518	kproc(jo) = nproc
519	kobsi(jo) = ji + 1
520	kobsj(jo) = jj + 1
521	EXIT gridloop
522	END IF
523	ENDIF
524	ENDIF
525	ENDDO
526	ENDDO gridloop
527
528	!---------------------------------------------------------------
529	! Ensure that all observation longtiudes are between -180/180
530	!---------------------------------------------------------------
531
532	IF ( zplam(jo) > 180 ) zplam(jo) = zplam(jo) - 360.0_wp
533
534	IF ( kproc(jo) == -1 ) THEN
535
536	! Normal case
537	gridpoints : DO jj = jjmin, jjmax
538	DO ji = jimin, jimax
539
540
541	IF ( ( zplam(jo) > zlamtmax(ji,jj) ) .OR. &
542	& ( zplam(jo) < zlamtmin(ji,jj) ) ) CYCLE
543
544	IF ( ABS( pphi(jo) ) < 85 ) THEN
545	IF ( ( pphi(jo) > zphitmax(ji,jj) ) .OR. &
546	& ( pphi(jo) < zphitmin(ji,jj) ) ) CYCLE
547	ENDIF
548
549	IF ( linquad( zplam(jo), pphi(jo), &
550	& zlamtm(:,ji,jj), zphitm(:,ji,jj) ) ) THEN
551	IF ( llinvalidcell(ji,jj) ) THEN
552	kproc(jo) = nproc + 1000000
553	kobsi(jo) = ji + 1
554	kobsj(jo) = jj + 1
555	CYCLE
556	ELSE
557	kproc(jo) = nproc
558	kobsi(jo) = ji + 1
559	kobsj(jo) = jj + 1
560	EXIT gridpoints
561	ENDIF
562	ENDIF
563
564	END DO
565	END DO gridpoints
566	ENDIF
567
568	! In case of failure retry for obs. longtiude + 360.
569	IF ( kproc(jo) == -1 ) THEN
570	gridpoints_greenwich : DO jj = jjmin, jjmax
571	DO ji = jimin, jimax
572
573	IF ( ( zplam(jo)+360.0_wp > zlamtmax(ji,jj) ) .OR. &
574	& ( zplam(jo)+360.0_wp < zlamtmin(ji,jj) ) ) CYCLE
575
576	IF ( ABS( pphi(jo) ) < 85 ) THEN
577	IF ( ( pphi(jo) > zphitmax(ji,jj) ) .OR. &
578	& ( pphi(jo) < zphitmin(ji,jj) ) ) CYCLE
579	ENDIF
580
581	IF ( linquad( zplam(jo)+360.0_wp, pphi(jo), &
582	& zlamtm(:,ji,jj), zphitm(:,ji,jj) ) ) THEN
583	IF ( llinvalidcell(ji,jj) ) THEN
584	kproc(jo) = nproc + 1000000
585	kobsi(jo) = ji + 1
586	kobsj(jo) = jj + 1
587	CYCLE
588	ELSE
589	kproc(jo) = nproc
590	kobsi(jo) = ji + 1
591	kobsj(jo) = jj + 1
592	EXIT gridpoints_greenwich
593	ENDIF
594	ENDIF
595
596	END DO
597	END DO gridpoints_greenwich
598
599	ENDIF ! kproc
600
601	END DO
602	END DO gridpointsn
603	END DO gpkobs ! kobs
604
605	!----------------------------------------------------------------------
606	! Synchronize kproc on all processors
607	!----------------------------------------------------------------------
608	IF ( ln_grid_global ) THEN
609	CALL obs_mpp_max_integer( kproc, kobs )
610	CALL obs_mpp_max_integer( kobsi, kobs )
611	CALL obs_mpp_max_integer( kobsj, kobs )
612	ELSE
613	CALL obs_mpp_find_obs_proc( kproc, kobsi, kobsj, kobs )
614	ENDIF
615
616	WHERE( kproc(:) >= 1000000 )
617	kproc(:) = kproc(:) - 1000000
618	END WHERE
619
620	DEALLOCATE( &
621	& zlamg, &
622	& zphig, &
623	& zmskg, &
624	& zphitmax, &
625	& zphitmin, &
626	& zlamtmax, &
627	& zlamtmin, &
628	& llinvalidcell, &
629	& zlamtm, &
630	& zphitm, &
631	& zplam &
632	& )
633
634	END SUBROUTINE obs_grid_search_lookup
635
636
637	SUBROUTINE obs_grid_setup
638	!!----------------------------------------------------------------------
639	!! * ROUTINE obs_grid_setup *
640	!!
641	!! ** Purpose : Setup a lookup table to reduce the searching required
642	!! for converting lat lons to grid point location
643	!! produces or reads in a preexisting file for use in
644	!! obs_grid_search_lookup_local
645	!!
646	!! ** Method : calls obs_grid_search_bruteforce_local with a array
647	!! of lats and lons
648	!!
649	!! History :
650	!! ! 2007-12 (D. Lea) new routine
651	!!----------------------------------------------------------------------
652
653	!! * Local declarations
654	CHARACTER(LEN=15), PARAMETER :: &
655	& cpname = 'obs_grid_setup'
656	CHARACTER(LEN=40) :: cfname
657	INTEGER :: &
658	& ji, &
659	& jj, &
660	& jk, &
661	& jo
662	INTEGER :: &
663	& idfile, idny, idnx, idxpos, idypos, idlat, idlon, fileexist
664	INTEGER, DIMENSION(2) :: incdim
665	CHARACTER(LEN=20) :: datestr=" ",timestr=" "
666	REAL(wp) :: tmpx1, tmpx2, tmpy1, tmpy2
667	REAL(wp) :: meanxdiff, meanydiff
668	REAL(wp) :: meanxdiff1, meanydiff1
669	REAL(wp) :: meanxdiff2, meanydiff2
670	INTEGER :: numx1, numx2, numy1, numy2, df
671	INTEGER :: &
672	& jimin, jimax, jjmin, jjmax
673	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: &
674	& lonsi, &
675	& latsi
676	INTEGER, DIMENSION(:,:), ALLOCATABLE :: &
677	& ixposi, &
678	& iyposi, &
679	& iproci
680	INTEGER, PARAMETER :: histsize=90
681	INTEGER, DIMENSION(histsize) :: &
682	& histx1, histx2, histy1, histy2
683	REAL, DIMENSION(histsize) :: &
684	& fhistx1, fhistx2, fhisty1, fhisty2
685	REAL(wp) :: histtol
686
687	IF (ln_grid_search_lookup) THEN
688
689	WRITE(numout,*) 'Calling obs_grid_setup'
690
691	#ifdef key_datetime_out
692	IF (lwp) THEN
693	CALL DATE_AND_TIME(datestr,timestr)
694	WRITE(numout,*) 'obs_grid_setup begin date_and_time ',datestr,' ',timestr
695	ENDIF
696	#endif
697
698	IF(lwp) WRITE(numout,*)
699	IF(lwp) WRITE(numout,*)'Grid search resolution : ', grid_search_res
700
701	gsearch_nlons_def = NINT( 360.0_wp / grid_search_res )
702	gsearch_nlats_def = NINT( 180.0_wp / grid_search_res )
703	gsearch_lonmin_def = -180.0_wp + 0.5_wp * grid_search_res
704	gsearch_latmin_def = -90.0_wp + 0.5_wp * grid_search_res
705	gsearch_dlon_def = grid_search_res
706	gsearch_dlat_def = grid_search_res
707
708	IF (lwp) THEN
709	WRITE(numout,*)'Grid search gsearch_nlons_def = ',gsearch_nlons_def
710	WRITE(numout,*)'Grid search gsearch_nlats_def = ',gsearch_nlats_def
711	WRITE(numout,*)'Grid search gsearch_lonmin_def = ',gsearch_lonmin_def
712	WRITE(numout,*)'Grid search gsearch_latmin_def = ',gsearch_latmin_def
713	WRITE(numout,*)'Grid search gsearch_dlon_def = ',gsearch_dlon_def
714	WRITE(numout,*)'Grid search gsearch_dlat_def = ',gsearch_dlat_def
715	ENDIF
716
717	IF ( ln_grid_global ) THEN
718	WRITE(cfname, FMT="(A,'_',A)") &
719	& TRIM(grid_search_file), 'global.nc'
720	ELSE
721	WRITE(cfname, FMT="(A,'_',I4.4,'of',I4.4,'by',I4.4,'.nc')") &
722	& TRIM(grid_search_file), nproc, jpni, jpnj
723	ENDIF
724
725	fileexist=nf90_open( TRIM( cfname ), nf90_nowrite, &
726	& idfile )
727
728	IF ( fileexist == nf90_noerr ) THEN
729
730	! read data
731	! initially assume size is as defined (to be fixed)
732
733	WRITE(numout,*) 'Reading: ',cfname
734
735	CALL chkerr( nf90_open( TRIM( cfname ), nf90_nowrite, idfile ), &
736	& cpname, __LINE__ )
737	CALL chkerr( nf90_get_att( idfile, nf90_global, 'maxxdiff', maxxdiff ), &
738	& cpname, __LINE__ )
739	CALL chkerr( nf90_get_att( idfile, nf90_global, 'maxydiff', maxydiff ), &
740	& cpname, __LINE__ )
741	CALL chkerr( nf90_get_att( idfile, nf90_global, 'dlon', dlon ), &
742	& cpname, __LINE__ )
743	CALL chkerr( nf90_get_att( idfile, nf90_global, 'dlat', dlat ), &
744	& cpname, __LINE__ )
745	CALL chkerr( nf90_get_att( idfile, nf90_global, 'lonmin', lonmin ), &
746	& cpname, __LINE__ )
747	CALL chkerr( nf90_get_att( idfile, nf90_global, 'latmin', latmin ), &
748	& cpname, __LINE__ )
749
750	CALL chkerr( nf90_inq_dimid(idfile, 'nx' , idnx), &
751	& cpname, __LINE__ )
752	CALL chkerr( nf90_inquire_dimension( idfile, idnx, len = nlons ), &
753	& cpname, __LINE__ )
754	CALL chkerr( nf90_inq_dimid(idfile, 'ny' , idny), &
755	& cpname, __LINE__ )
756	CALL chkerr( nf90_inquire_dimension( idfile, idny, len = nlats ), &
757	& cpname, __LINE__ )
758
759	ALLOCATE( &
760	& lons(nlons,nlats), &
761	& lats(nlons,nlats), &
762	& ixpos(nlons,nlats), &
763	& iypos(nlons,nlats), &
764	& iproc(nlons,nlats) &
765	& )
766
767	CALL chkerr( nf90_inq_varid( idfile, 'XPOS', idxpos ), &
768	& cpname, __LINE__ )
769	CALL chkerr( nf90_get_var ( idfile, idxpos, ixpos), &
770	& cpname, __LINE__ )
771	CALL chkerr( nf90_inq_varid( idfile, 'YPOS', idypos ), &
772	& cpname, __LINE__ )
773	CALL chkerr( nf90_get_var ( idfile, idypos, iypos), &
774	& cpname, __LINE__ )
775
776	CALL chkerr( nf90_close( idfile ), cpname, __LINE__ )
777
778	! setup arrays
779
780	DO ji = 1, nlons
781	DO jj = 1, nlats
782	lons(ji,jj) = lonmin + (ji-1) * dlon
783	lats(ji,jj) = latmin + (jj-1) * dlat
784	ENDDO
785	ENDDO
786
787	! if we are not reading the file we need to create it
788	! create new obs grid search lookup file
789
790	ELSE
791
792	! call obs_grid_search
793
794	IF (lwp) THEN
795	WRITE(numout,*) 'creating: ',cfname
796	WRITE(numout,) 'calling obs_grid_search: ',nlonsnlats
797	ENDIF
798
799	! set parameters from default values
800	nlons = gsearch_nlons_def
801	nlats = gsearch_nlats_def
802	lonmin = gsearch_lonmin_def
803	latmin = gsearch_latmin_def
804	dlon = gsearch_dlon_def
805	dlat = gsearch_dlat_def
806
807	! setup arrays
808
809	ALLOCATE( &
810	& lonsi(nlons,nlats), &
811	& latsi(nlons,nlats), &
812	& ixposi(nlons,nlats), &
813	& iyposi(nlons,nlats), &
814	& iproci(nlons,nlats) &
815	& )
816
817	DO ji = 1, nlons
818	DO jj = 1, nlats
819	lonsi(ji,jj) = lonmin + (ji-1) * dlon
820	latsi(ji,jj) = latmin + (jj-1) * dlat
821	ENDDO
822	ENDDO
823
824	CALL obs_grid_search_bruteforce( jpi, jpj, jpiglo, jpjglo, &
825	& nldi, nlei,nldj, nlej, &
826	& nproc, jpnij, &
827	& glamt, gphit, tmask, &
828	& nlons*nlats, lonsi, latsi, &
829	& ixposi, iyposi, iproci )
830
831	! minimise file size by removing regions with no data from xypos file
832	! should be able to just use xpos (ypos will have the same areas of missing data)
833
834	jimin=1
835	jimax=nlons
836	jjmin=1
837	jjmax=nlats
838
839	minlon_xpos: DO ji= 1, nlons
840	IF (COUNT(ixposi(ji,:) >= 0) > 0) THEN
841	jimin=ji
842	EXIT minlon_xpos
843	ENDIF
844	ENDDO minlon_xpos
845
846	maxlon_xpos: DO ji= nlons, 1, -1
847	IF (COUNT(ixposi(ji,:) >= 0) > 0) THEN
848	jimax=ji
849	EXIT maxlon_xpos
850	ENDIF
851	ENDDO maxlon_xpos
852
853	minlat_xpos: DO jj= 1, nlats
854	IF (COUNT(ixposi(:,jj) >= 0) > 0) THEN
855	jjmin=jj
856	EXIT minlat_xpos
857	ENDIF
858	ENDDO minlat_xpos
859
860	maxlat_xpos: DO jj= nlats, 1, -1
861	IF (COUNT(ixposi(:,jj) >= 0) > 0) THEN
862	jjmax=jj
863	EXIT maxlat_xpos
864	ENDIF
865	ENDDO maxlat_xpos
866
867	lonmin = lonsi(jimin,jjmin)
868	latmin = latsi(jimin,jjmin)
869	nlons = jimax-jimin+1
870	nlats = jjmax-jjmin+1
871
872	! construct new arrays
873
874	ALLOCATE( &
875	& lons(nlons,nlats), &
876	& lats(nlons,nlats), &
877	& ixpos(nlons,nlats), &
878	& iypos(nlons,nlats), &
879	& iproc(nlons,nlats) &
880	& )
881
882	lons(:,:) = lonsi(jimin:jimax,jjmin:jjmax)
883	lats(:,:) = latsi(jimin:jimax,jjmin:jjmax)
884	ixpos(:,:) = ixposi(jimin:jimax,jjmin:jjmax)
885	iypos(:,:) = iyposi(jimin:jimax,jjmin:jjmax)
886	iproc(:,:) = iproci(jimin:jimax,jjmin:jjmax)
887
888	DEALLOCATE(lonsi,latsi,ixposi,iyposi,iproci)
889
890	! calculate (estimate) maxxdiff, maxydiff
891	! this is to help define the search area for obs_grid_search_lookup
892
893	maxxdiff = 1
894	maxydiff = 1
895
896	tmpx1 = 0
897	tmpx2 = 0
898	tmpy1 = 0
899	tmpy2 = 0
900
901	numx1 = 0
902	numx2 = 0
903	numy1 = 0
904	numy2 = 0
905
906	! calculate the mean absolute xdiff and ydiff
907	! also calculate a histogram
908	! note the reason why looking for xdiff and ydiff in both directions
909	! is to allow for rotated grids
910
911	DO ji = 1, nlons-1
912	DO jj = 1, nlats-1
913	IF ( ixpos(ji,jj) > 0 .AND. iypos(ji,jj) > 0 ) THEN
914	IF ( ixpos(ji+1,jj) > 0 ) THEN
915	df = ABS( ixpos(ji+1,jj) - ixpos(ji,jj) )
916	tmpx1 = tmpx1+df
917	numx1 = numx1+1
918	IF ( df < histsize ) histx1(df+1) = histx1(df+1) + 1
919	ENDIF
920	IF ( ixpos(ji,jj+1) > 0 ) THEN
921	df = ABS( ixpos(ji,jj+1) - ixpos(ji,jj) )
922	tmpx2 = tmpx2 + df
923	numx2 = numx2 + 1
924	IF ( df < histsize ) histx2(df+1) = histx2(df+1) + 1
925	ENDIF
926	IF (iypos(ji+1,jj) > 0) THEN
927	df = ABS( iypos(ji+1,jj) - iypos(ji,jj) )
928	tmpy1 = tmpy1 + df
929	numy1 = numy1 + 1
930	IF ( df < histsize ) histy1(df+1) = histy1(df+1) + 1
931	ENDIF
932	IF ( iypos(ji,jj+1) > 0 ) THEN
933	df = ABS( iypos(ji,jj+1) - iypos(ji,jj) )
934	tmpy2 = tmpy2 + df
935	numy2 = numy2 + 1
936	IF ( df < histsize ) histy2(df+1) = histy2(df+1) + 1
937	ENDIF
938	ENDIF
939	ENDDO
940	ENDDO
941
942	IF (lwp) THEN
943	WRITE(numout,*) 'histograms'
944	WRITE(numout,*) '0 1 2 3 4 5 6 7 8 9 10 ...'
945	WRITE(numout,*) 'histx1'
946	WRITE(numout,*) histx1
947	WRITE(numout,*) 'histx2'
948	WRITE(numout,*) histx2
949	WRITE(numout,*) 'histy1'
950	WRITE(numout,*) histy1
951	WRITE(numout,*) 'histy2'
952	WRITE(numout,*) histy2
953	ENDIF
954
955	meanxdiff1 = tmpx1 / numx1
956	meanydiff1 = tmpy1 / numy1
957	meanxdiff2 = tmpx2 / numx2
958	meanydiff2 = tmpy2 / numy2
959
960	meanxdiff = MAXVAL((/ meanxdiff1, meanxdiff2 /))
961	meanydiff = MAXVAL((/ meanydiff1, meanydiff2 /))
962
963	IF (lwp) THEN
964	WRITE(numout,*) tmpx1, tmpx2, tmpy1, tmpy2
965	WRITE(numout,*) numx1, numx2, numy1, numy2
966	WRITE(numout,*) 'meanxdiff: ',meanxdiff, meanxdiff1, meanxdiff2
967	WRITE(numout,*) 'meanydiff: ',meanydiff, meanydiff1, meanydiff2
968	ENDIF
969
970	tmpx1 = 0
971	tmpx2 = 0
972	tmpy1 = 0
973	tmpy2 = 0
974
975	numx1 = 0
976	numx2 = 0
977	numy1 = 0
978	numy2 = 0
979
980	histx1(:) = 0
981	histx2(:) = 0
982	histy1(:) = 0
983	histy2(:) = 0
984
985	limxdiff = meanxdiff * 4! limit the difference to avoid picking up wraparound
986	limydiff = meanydiff * 4
987
988	DO ji = 1, nlons-1
989	DO jj = 1, nlats-1
990	IF ( ixpos(ji,jj) > 0 .AND. iypos(ji,jj) > 0 ) THEN
991
992	IF ( ixpos(ji+1,jj) > 0 ) THEN
993	df = ABS( ixpos(ji+1,jj)-ixpos(ji,jj) )
994	tmpx1 = df
995	IF ( df < limxdiff ) numx1 = numx1+1
996	IF ( df < histsize ) histx1(df+1) = histx1(df+1) + 1
997	ENDIF
998	IF ( ixpos(ji,jj+1) > 0 ) THEN
999	df = ABS( ixpos(ji,jj+1) - ixpos(ji,jj) )
1000	tmpx2 = df
1001	IF ( df < limxdiff ) numx2 = numx2 + 1
1002	IF ( df < histsize ) histx2(df+1) = histx2(df+1) + 1
1003	ENDIF
1004	IF (iypos(ji+1,jj) > 0) THEN
1005	df = ABS( iypos(ji+1,jj) - iypos(ji,jj) )
1006	tmpy1 = df
1007	IF ( df < limxdiff ) numy1 = numy1 + 1
1008	IF ( df < histsize ) histy1(df+1) = histy1(df+1) + 1
1009	ENDIF
1010	IF (iypos(ji,jj+1) > 0) THEN
1011	df = ABS( iypos(ji,jj+1) - iypos(ji,jj) )
1012	tmpy2 = df
1013	IF ( df < limxdiff ) numy2 = numy2+1
1014	IF ( df < histsize ) histy2(df+1) = histy2(df+1)+1
1015	ENDIF
1016
1017	IF ( maxxdiff < tmpx1 .AND. tmpx1 < limxdiff ) &
1018	& maxxdiff = tmpx1
1019	IF ( maxxdiff < tmpx2 .AND. tmpx2 < limxdiff ) &
1020	& maxxdiff = tmpx2
1021	IF ( maxydiff < tmpy1 .AND. tmpy1 < limydiff ) &
1022	& maxydiff = tmpy1
1023	IF ( maxydiff < tmpy2 .AND. tmpy2 < limydiff ) &
1024	& maxydiff = tmpy2
1025
1026	ENDIF
1027	END DO
1028	END DO
1029
1030	! cumulative histograms
1031
1032	DO ji = 1, histsize - 1
1033	histx1(ji+1) = histx1(ji+1) + histx1(ji)
1034	histx2(ji+1) = histx2(ji+1) + histx2(ji)
1035	histy1(ji+1) = histy1(ji+1) + histy1(ji)
1036	histy2(ji+1) = histy2(ji+1) + histy2(ji)
1037	ENDDO
1038
1039	fhistx1(:) = histx1(:) * 1.0 / numx1
1040	fhistx2(:) = histx2(:) * 1.0 / numx2
1041	fhisty1(:) = histy1(:) * 1.0 / numy1
1042	fhisty2(:) = histy2(:) * 1.0 / numy2
1043
1044	! output new histograms
1045
1046	IF (lwp) THEN
1047	WRITE(numout,*) 'cumulative histograms'
1048	WRITE(numout,*) '0 1 2 3 4 5 6 7 8 9 10 ...'
1049	WRITE(numout,*) 'fhistx1'
1050	WRITE(numout,*) fhistx1
1051	WRITE(numout,*) 'fhistx2'
1052	WRITE(numout,*) fhistx2
1053	WRITE(numout,*) 'fhisty1'
1054	WRITE(numout,*) fhisty1
1055	WRITE(numout,*) 'fhisty2'
1056	WRITE(numout,*) fhisty2
1057	ENDIF
1058
1059	! calculate maxxdiff and maxydiff based on cumulative histograms
1060	! where > 0.999 of points are
1061
1062	! maxval just converts 1x1 vector return from maxloc to a scalar
1063
1064	histtol = 0.999
1065	tmpx1 = MAXVAL( MAXLOC( fhistx1(:), mask = ( fhistx1(:) <= histtol ) ) )
1066	tmpx2 = MAXVAL( MAXLOC( fhistx2(:), mask = ( fhistx2(:) <= histtol ) ) )
1067	tmpy1 = MAXVAL( MAXLOC( fhisty1(:), mask = ( fhisty1(:) <= histtol ) ) )
1068	tmpy2 = MAXVAL( MAXLOC( fhisty2(:), mask = ( fhisty2(:) <= histtol ) ) )
1069
1070	maxxdiff = MAXVAL( (/ tmpx1, tmpx2 /) ) + 1
1071	maxydiff = MAXVAL( (/ tmpy1, tmpy2 /) ) + 1
1072
1073	! Write out data
1074
1075	IF ( ( .NOT. ln_grid_global ) .OR. &
1076	& ( ( ln_grid_global ) .AND. ( nproc==0 ) ) ) THEN
1077
1078	CALL chkerr( nf90_create (TRIM(cfname), nf90_clobber, idfile), &
1079	& cpname, __LINE__ )
1080	CALL chkerr( nf90_put_att( idfile, nf90_global, 'title', &
1081	& 'Mapping file from lon/lat to model grid point' ),&
1082	& cpname,__LINE__ )
1083	CALL chkerr( nf90_put_att( idfile, nf90_global, 'maxxdiff', &
1084	& maxxdiff ), &
1085	& cpname,__LINE__ )
1086	CALL chkerr( nf90_put_att( idfile, nf90_global, 'maxydiff', &
1087	& maxydiff ), &
1088	& cpname,__LINE__ )
1089	CALL chkerr( nf90_put_att( idfile, nf90_global, 'dlon', dlon ),&
1090	& cpname,__LINE__ )
1091	CALL chkerr( nf90_put_att( idfile, nf90_global, 'dlat', dlat ),&
1092	& cpname,__LINE__ )
1093	CALL chkerr( nf90_put_att( idfile, nf90_global, 'lonmin', &
1094	& lonmin ), &
1095	& cpname,__LINE__ )
1096	CALL chkerr( nf90_put_att( idfile, nf90_global, 'latmin', &
1097	& latmin ), &
1098	& cpname,__LINE__ )
1099
1100	CALL chkerr( nf90_def_dim(idfile, 'nx' , nlons, idnx), &
1101	& cpname,__LINE__ )
1102	CALL chkerr( nf90_def_dim(idfile, 'ny' , nlats, idny), &
1103	& cpname,__LINE__ )
1104
1105	incdim(1) = idnx
1106	incdim(2) = idny
1107
1108	CALL chkerr( nf90_def_var( idfile, 'LON', nf90_float, incdim, &
1109	& idlon ), &
1110	& cpname, __LINE__ )
1111	CALL chkerr( nf90_put_att( idfile, idlon, 'long_name', &
1112	& 'longitude' ), &
1113	& cpname, __LINE__ )
1114
1115	CALL chkerr( nf90_def_var( idfile, 'LAT', nf90_float, incdim, &
1116	& idlat ), &
1117	& cpname, __LINE__ )
1118	CALL chkerr( nf90_put_att( idfile, idlat, 'long_name', &
1119	& 'latitude' ), &
1120	& cpname, __LINE__ )
1121
1122	CALL chkerr( nf90_def_var( idfile, 'XPOS', nf90_int, incdim, &
1123	& idxpos ), &
1124	& cpname, __LINE__ )
1125	CALL chkerr( nf90_put_att( idfile, idxpos, 'long_name', &
1126	& 'x position' ), &
1127	& cpname, __LINE__ )
1128	CALL chkerr( nf90_put_att( idfile, idxpos, '_FillValue', -1 ), &
1129	& cpname, __LINE__ )
1130
1131	CALL chkerr( nf90_def_var( idfile, 'YPOS', nf90_int, incdim, &
1132	& idypos ), &
1133	& cpname, __LINE__ )
1134	CALL chkerr( nf90_put_att( idfile, idypos, 'long_name', &
1135	& 'y position' ), &
1136	& cpname, __LINE__ )
1137	CALL chkerr( nf90_put_att( idfile, idypos, '_FillValue', -1 ), &
1138	& cpname, __LINE__ )
1139
1140	CALL chkerr( nf90_enddef( idfile ), cpname, __LINE__ )
1141
1142	CALL chkerr( nf90_put_var( idfile, idlon, lons), &
1143	& cpname, __LINE__ )
1144	CALL chkerr( nf90_put_var( idfile, idlat, lats), &
1145	& cpname, __LINE__ )
1146	CALL chkerr( nf90_put_var( idfile, idxpos, ixpos), &
1147	& cpname, __LINE__ )
1148	CALL chkerr( nf90_put_var( idfile, idypos, iypos), &
1149	& cpname, __LINE__ )
1150
1151	CALL chkerr( nf90_close( idfile ), cpname, __LINE__ )
1152
1153	! should also output max i, max j spacing for use in
1154	! obs_grid_search_lookup
1155
1156	END IF
1157
1158	ENDIF
1159
1160	ENDIF
1161
1162	END SUBROUTINE obs_grid_setup
1163
1164	SUBROUTINE obs_grid_deallocate( )
1165	!!----------------------------------------------------------------------
1166	!! * ROUTINE obs_grid_setup *
1167	!!
1168	!! ** Purpose : Deallocate arrays setup by obs_grid_setup
1169	!!
1170	!! History :
1171	!! ! 2007-12 (D. Lea) new routine
1172	!!-----------------------------------------------------------------------
1173
1174	IF (ln_grid_search_lookup) THEN
1175	DEALLOCATE( lons, lats, ixpos, iypos, iproc )
1176	ENDIF
1177
1178	END SUBROUTINE obs_grid_deallocate
1179
1180	#include "obs_level_search.h90"
1181
1182	#include "linquad.h90"
1183
1184	#include "maxdist.h90"
1185
1186	#include "find_obs_proc.h90"
1187
1188	END MODULE obs_grid
1189

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

Download in other formats: