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.

grid.f90 in utils/tools/SIREN/src – NEMO

Context Navigation

source: utils/tools/SIREN/src/grid.f90 @ 12080

Last change on this file since 12080 was 12080, checked in by jpaul, 4 years ago
update nemo trunk
File size: 226.9 KB

Line
1	!----------------------------------------------------------------------
2	! NEMO system team, System and Interface for oceanic RElocable Nesting
3	!----------------------------------------------------------------------
4	!
5	! DESCRIPTION:
6	!> @brief This module is grid manager.
7	!>
8	!> @details
9	!> to get NEMO pivot point index:<br/>
10	!> @code
11	!> il_pivot=grid_get_pivot(td_file)
12	!> @endcode
13	!> - il_pivot is NEMO pivot point index F(0), T(1)
14	!> - td_file is mpp structure
15	!>
16	!> to get NEMO periodicity index:<br/>
17	!> @code
18	!> il_perio=grid_get_perio(td_file)
19	!> @endcode
20	!> - il_perio is NEMO periodicity index (0,1,2,3,4,5,6)
21	!> - td_file is mpp structure
22	!>
23	!> to check domain validity:<br/>
24	!> @code
25	!> CALL grid_check_dom(td_coord, id_imin, id_imax, id_jmin, id_jmax)
26	!> @endcode
27	!> - td_coord is coordinates mpp structure
28	!> - id_imin is i-direction lower left point indice
29	!> - id_imax is i-direction upper right point indice
30	!> - id_jmin is j-direction lower left point indice
31	!> - id_jmax is j-direction upper right point indice
32	!>
33	!> to get closest coarse grid indices of fine grid domain:<br/>
34	!> @code
35	!> il_index(:,:)=grid_get_coarse_index(td_coord0, td_coord1,
36	!> [id_rho,] [cd_point])
37	!> @endcode
38	!> or
39	!> @code
40	!> il_index(:,:)=grid_get_coarse_index(td_lon0, td_lat0, td_coord1,
41	!> [id_rho,] [cd_point])
42	!> @endcode
43	!> or
44	!> @code
45	!> il_index(:,:)=grid_get_coarse_index(td_coord0, td_lon1, td_lat1,
46	!> [id_rho,] [cd_point])
47	!> @endcode
48	!> or
49	!> @code
50	!> il_index(:,:)=grid_get_coarse_index(td_lon0, td_lat0, td_lon1, td_lat1,
51	!> [id_rho,] [cd_point])
52	!> @endcode
53	!> - il_index(:,:) is coarse grid indices (/ (/ imin0, imax0 /),
54	!> (/ jmin0, jmax0 /) /)
55	!> - td_coord0 is coarse grid coordinate mpp structure
56	!> - td_coord1 is fine grid coordinate mpp structure
57	!> - td_lon0 is coarse grid longitude variable structure
58	!> - td_lat0 is coarse grid latitude variable structure
59	!> - td_lon1 is fine grid longitude variable structure
60	!> - td_lat1 is fine grid latitude variable structure
61	!> - id_rho is array of refinment factor (default 1)
62	!> - cd_point is Arakawa grid point (default 'T')
63	!>
64	!> to know if grid is global:<br/>
65	!> @code
66	!> ll_global=grid_is_global(td_lon, td_lat)
67	!> @endcode
68	!> - td_lon is longitude variable structure
69	!> - td_lat is latitude variable structure
70	!>
71	!> to know if grid contains north fold:<br/>
72	!> @code
73	!> ll_north=grid_is_north_fold(td_lat)
74	!> @endcode
75	!> - td_lat is latitude variable structure
76	!>
77	!> to get coarse grid indices of the closest point from one fine grid
78	!> point:<br/>
79	!> @code
80	!> il_index(:)=grid_get_closest(dd_lon0(:,:), dd_lat0(:,:), dd_lon1, dd_lat1
81	!> [,dd_fill] [,cd_pos])
82	!> @endcode
83	!> - il_index(:) is coarse grid indices (/ i0, j0 /)
84	!> - dd_lon0 is coarse grid array of longitude value (real(8))
85	!> - dd_lat0 is coarse grid array of latitude value (real(8))
86	!> - dd_lon1 is fine grid longitude value (real(8))
87	!> - dd_lat1 is fine grid latitude value (real(8))
88	!> - dd_fill
89	!> - cd_pos
90	!>
91	!> to compute distance between a point A and grid points:<br/>
92	!> @code
93	!> il_dist(:,:)=grid_distance(dd_lon, dd_lat, dd_lonA, dd_latA)
94	!> @endcode
95	!> - il_dist(:,:) is array of distance between point A and grid points
96	!> - dd_lon is array of longitude value (real(8))
97	!> - dd_lat is array of longitude value (real(8))
98	!> - dd_lonA is longitude of point A (real(8))
99	!> - dd_latA is latitude of point A (real(8))
100	!>
101	!> to get offset between fine grid and coarse grid:<br/>
102	!> @code
103	!> il_offset(:,:)=grid_get_fine_offset(td_coord0,
104	!> id_imin0, id_jmin0, id_imax0, id_jmax0,
105	!> td_coord1
106	!> [,id_rho] [,cd_point])
107	!> @endcode
108	!> or
109	!> @code
110	!> il_offset(:,:)=grid_get_fine_offset(dd_lon0, dd_lat0,
111	!> id_imin0, id_jmin0,id_imax0, id_jmax0,
112	!> td_coord1
113	!> [,id_rho] [,cd_point])
114	!> @endcode
115	!> or
116	!> @code
117	!> il_offset(:,:)=grid_get_fine_offset(td_coord0,
118	!> id_imin0, id_jmin0, id_imax0, id_jmax0,
119	!> dd_lon1, dd_lat1
120	!> [,id_rho] [,cd_point])
121	!> @endcode
122	!> or
123	!> @code
124	!> il_offset(:,:)=grid_get_fine_offset(dd_lon0, dd_lat0,
125	!> id_imin0, id_jmin0, id_imax0, id_jmax0,
126	!> dd_lon1, dd_lat1
127	!> [,id_rho] [,cd_point])
128	!> @endcode
129	!> - il_offset(:,:) is offset array
130	!> (/ (/ i_offset_left, i_offset_right /), (/ j_offset_lower, j_offset_upper /) /)
131	!> - td_coord0 is coarse grid coordinate mpp structure
132	!> - dd_lon0 is coarse grid longitude array (real(8))
133	!> - dd_lat0 is coarse grid latitude array (real(8))
134	!> - id_imin0 is coarse grid lower left corner i-indice of fine grid
135	!> domain
136	!> - id_jmin0 is coarse grid lower left corner j-indice of fine grid
137	!> domain
138	!> - id_imax0 is coarse grid upper right corner i-indice of fine grid
139	!> domain
140	!> - id_jmax0 is coarse grid upper right corner j-indice of fine grid
141	!> domain
142	!> - td_coord1 is fine grid coordinate mpp structure
143	!> - dd_lon1 is fine grid longitude array (real(8))
144	!> - dd_lat1 is fine grid latitude array (real(8))
145	!> - id_rho is array of refinment factor (default 1)
146	!> - cd_point is Arakawa grid point (default 'T')
147	!>
148	!> to check fine and coarse grid coincidence:<br/>
149	!> @code
150	!> CALL grid_check_coincidence(td_coord0, td_coord1,
151	!> id_imin0, id_imax0, id_jmin0, id_jmax0
152	!> ,id_rho)
153	!> @endcode
154	!> - td_coord0 is coarse grid coordinate mpp structure
155	!> - td_coord1 is fine grid coordinate mpp structure
156	!> - id_imin0 is coarse grid lower left corner i-indice of fine grid
157	!> domain
158	!> - id_imax0 is coarse grid upper right corner i-indice of fine grid
159	!> domain
160	!> - id_jmin0 is coarse grid lower left corner j-indice of fine grid
161	!> domain
162	!> - id_jmax0 is coarse grid upper right corner j-indice of fine grid
163	!> domain
164	!> - id_rho is array of refinement factor
165	!>
166	!> to add ghost cell at boundaries:<br/>
167	!> @code
168	!> CALL grid_add_ghost(td_var, id_ghost)
169	!> @endcode
170	!> - td_var is array of variable structure
171	!> - id_ghost is 2D array of ghost cell factor
172	!>
173	!> to delete ghost cell at boundaries:<br/>
174	!> @code
175	!> CALL grid_del_ghost(td_var, id_ghost)
176	!> @endcode
177	!> - td_var is array of variable structure
178	!> - id_ghost is 2D array of ghost cell factor
179	!>
180	!> to get ghost cell factor (use or not):<br/>
181	!> @code
182	!> il_factor(:)= grid_get_ghost( td_var )
183	!> @endcode
184	!> or
185	!> @code
186	!> il_factor(:)= grid_get_ghost( td_mpp )
187	!> @endcode
188	!> - il_factor(:) is array of ghost cell factor (0 or 1)
189	!> - td_var is variable structure
190	!> - td_mpp is mpp sturcture
191	!>
192	!> to compute closed sea domain:<br/>
193	!> @code
194	!> il_mask(:,:)=grid_split_domain(td_var, [id_level])
195	!> @endcode
196	!> - il_mask(:,:) is domain mask
197	!> - td_var is variable strucutre
198	!> - id_level is level to be used [optional]
199	!>
200	!> to fill small closed sea with _FillValue:<br/>
201	!> @code
202	!> CALL grid_fill_small_dom(td_var, id_mask, [id_minsize])
203	!> @endcode
204	!> - td_var is variable structure
205	!> - id_mask is domain mask (from grid_split_domain)
206	!> - id_minsize is minimum size of sea to be kept [optional]
207	!>
208	!> @author
209	!> J.Paul
210	!>
211	!> @date November, 2013 - Initial Version
212	!> @date September, 2014
213	!> - add header
214	!> @date October, 2014
215	!> - use mpp file structure instead of file
216	!> @date February, 2015
217	!> - add function grid_fill_small_msk to fill small domain inside bigger one
218	!> @date February, 2016
219	!> - improve way to check coincidence (bug fix)
220	!> - manage grid cases for T,U,V or F point, with even or odd refinment (bug fix)
221	!> @date April, 2016
222	!> - add function to get closest grid point using coarse grid coordinates strucutre
223	!> @date May, 2019
224	!> - define as module variable im_max_overlap
225	!>
226	!> @note Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
227	!----------------------------------------------------------------------
228	MODULE grid
229
230	USE netcdf
231	USE kind ! F90 kind parameter
232	USE fct ! basic usefull function
233	USE global ! global parameter
234	USE phycst ! physical constant
235	USE logger ! log file manager
236	USE file ! file manager
237	USE att ! attribute manager
238	USE var ! variable manager
239	USE dim ! dimension manager
240	USE iom ! I/O manager
241	USE mpp ! MPP manager
242	USE dom ! domain manager
243	USE iom_mpp ! MPP I/O manager
244	USE iom_dom ! DOM I/O manager
245
246	IMPLICIT NONE
247	! NOTE_avoid_public_variables_if_possible
248
249	! type and variable
250	INTEGER(i4), PARAMETER :: im_max_overlap = 5
251
252	! function and subroutine
253	PUBLIC :: grid_get_info !< get information about mpp global domain (pivot, perio, ew)
254	PUBLIC :: grid_get_pivot !< get NEMO pivot point index
255	PUBLIC :: grid_get_perio !< get NEMO periodicity index
256	PUBLIC :: grid_get_ew_overlap !< get East West overlap
257	PUBLIC :: grid_check_dom !< check domain validity
258	PUBLIC :: grid_get_coarse_index !< get closest coarse grid indices of fine grid domain.
259	PUBLIC :: grid_is_global !< check if grid is global or not
260	PUBLIC :: grid_is_north_fold
261	PUBLIC :: grid_get_closest !< return closest coarse grid point from another point
262	PUBLIC :: grid_distance !< compute grid distance to a point
263	PUBLIC :: grid_get_fine_offset !< get fine grid offset
264	PUBLIC :: grid_check_coincidence !< check fine and coarse grid coincidence
265	PUBLIC :: grid_add_ghost !< add ghost cell at boundaries.
266	PUBLIC :: grid_del_ghost !< delete ghost cell at boundaries.
267	PUBLIC :: grid_get_ghost !< return ghost cell factor
268	PUBLIC :: grid_split_domain !< compute closed sea domain
269	PUBLIC :: grid_fill_small_dom !< fill small closed sea with fill value
270	PUBLIC :: grid_fill_small_msk !< fill small domain inside bigger one
271
272	! get closest coarse grid indices of fine grid domain
273	PRIVATE :: grid__get_coarse_index_ff ! - using coarse and fine grid coordinates files
274	PRIVATE :: grid__get_coarse_index_cf ! - using coarse grid array of lon,lat and fine grid coordinates files
275	PRIVATE :: grid__get_coarse_index_fc ! - using coarse grid coordinates files, and fine grid array of lon,lat
276	PRIVATE :: grid__get_coarse_index_cc ! - using coarse and fine grid array of lon,lat
277
278	! return closest coarse grid point from another point
279	PRIVATE :: grid__get_closest_str ! - using coarse grid coordinates strucutre
280	PRIVATE :: grid__get_closest_arr ! - using coarse grid array of lon,lat
281
282	! get offset between fine and coarse grid
283	PRIVATE :: grid__get_fine_offset_ff ! - using coarse and fine grid coordinates files
284	PRIVATE :: grid__get_fine_offset_cf ! - using coarse grid array of lon,lat and fine grid coordinates files
285	PRIVATE :: grid__get_fine_offset_fc ! - using coarse grid coordinates files, and fine grid array of lon,lat
286	PRIVATE :: grid__get_fine_offset_cc ! - using coarse and fine grid array of lon,lat
287
288	! get information about global domain (pivot, perio, ew)
289	PRIVATE :: grid__get_info_mpp ! - using mpp files structure
290	PRIVATE :: grid__get_info_file ! - using files structure
291
292	! get NEMO pivot point index
293	PRIVATE :: grid__get_pivot_mpp ! - using mpp files structure
294	PRIVATE :: grid__get_pivot_file ! - using files structure
295	PRIVATE :: grid__get_pivot_var ! - using variable structure
296	PRIVATE :: grid__get_pivot_varT ! compute NEMO pivot point index for variable on grid T
297	PRIVATE :: grid__get_pivot_varU ! compute NEMO pivot point index for variable on grid U
298	PRIVATE :: grid__get_pivot_varV ! compute NEMO pivot point index for variable on grid V
299	PRIVATE :: grid__get_pivot_varF ! compute NEMO pivot point index for variable on grid F
300
301	! get NEMO periodicity index
302	PRIVATE :: grid__get_perio_mpp ! - using mpp files structure
303	PRIVATE :: grid__get_perio_file ! - using files structure
304	PRIVATE :: grid__get_perio_var ! - using variable structure
305
306	! get East West overlap
307	PRIVATE :: grid__get_ew_overlap_mpp ! - using mpp files structure
308	PRIVATE :: grid__get_ew_overlap_file ! - using files structure
309	PRIVATE :: grid__get_ew_overlap_var ! - using longitude variable structure
310
311	! return ghost cell factor
312	PRIVATE :: grid__get_ghost_mpp ! - using mpp files structure
313	PRIVATE :: grid__get_ghost_var ! - using array of lon,lat
314	PRIVATE :: grid__check_corner ! check that fine grid is inside coarse grid
315	PRIVATE :: grid__check_lat ! check that fine grid latitude are inside coarse grid latitude
316
317	INTERFACE grid_get_info
318	MODULE PROCEDURE grid__get_info_mpp
319	MODULE PROCEDURE grid__get_info_file
320	END INTERFACE grid_get_info
321
322	INTERFACE grid_get_pivot
323	MODULE PROCEDURE grid__get_pivot_mpp
324	MODULE PROCEDURE grid__get_pivot_file
325	MODULE PROCEDURE grid__get_pivot_var
326	END INTERFACE grid_get_pivot
327
328	INTERFACE grid_get_perio
329	MODULE PROCEDURE grid__get_perio_mpp
330	MODULE PROCEDURE grid__get_perio_file
331	MODULE PROCEDURE grid__get_perio_var
332	END INTERFACE grid_get_perio
333
334	INTERFACE grid_get_ew_overlap
335	MODULE PROCEDURE grid__get_ew_overlap_mpp
336	MODULE PROCEDURE grid__get_ew_overlap_file
337	MODULE PROCEDURE grid__get_ew_overlap_var
338	END INTERFACE grid_get_ew_overlap
339
340	INTERFACE grid_get_ghost
341	MODULE PROCEDURE grid__get_ghost_var
342	MODULE PROCEDURE grid__get_ghost_mpp
343	END INTERFACE grid_get_ghost
344
345	INTERFACE grid_get_closest
346	MODULE PROCEDURE grid__get_closest_str
347	MODULE PROCEDURE grid__get_closest_arr
348	END INTERFACE grid_get_closest
349
350	INTERFACE grid_get_coarse_index
351	MODULE PROCEDURE grid__get_coarse_index_ff
352	MODULE PROCEDURE grid__get_coarse_index_cf
353	MODULE PROCEDURE grid__get_coarse_index_fc
354	MODULE PROCEDURE grid__get_coarse_index_cc
355	END INTERFACE grid_get_coarse_index
356
357	INTERFACE grid_get_fine_offset
358	MODULE PROCEDURE grid__get_fine_offset_ff
359	MODULE PROCEDURE grid__get_fine_offset_fc
360	MODULE PROCEDURE grid__get_fine_offset_cf
361	MODULE PROCEDURE grid__get_fine_offset_cc
362	END INTERFACE grid_get_fine_offset
363
364	CONTAINS
365	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
366	SUBROUTINE grid__get_info_file(td_file)
367	!-------------------------------------------------------------------
368	!> @brief This subroutine get information about global domain, given file
369	!> strucutre.
370	!>
371	!> @details
372	!> open edge files then:
373	!> - compute NEMO pivot point
374	!> - compute NEMO periodicity
375	!> - compute East West overlap
376	!>
377	!> @note need all processor files to be there
378	!> @author J.Paul
379	!> @date October, 2014 - Initial Version
380	!>
381	!> @param[inout] td_file file structure
382	!-------------------------------------------------------------------
383
384	IMPLICIT NONE
385
386	! Argument
387	TYPE(TFILE), INTENT(INOUT) :: td_file
388
389	! local variable
390	INTEGER(i4) :: il_ew
391	INTEGER(i4) :: il_pivot
392	INTEGER(i4) :: il_perio
393	INTEGER(i4) :: il_attid
394
395	TYPE(TATT) :: tl_att
396
397	TYPE(TFILE) :: tl_file
398
399	! loop indices
400	INTEGER(i4) :: ji
401	!----------------------------------------------------------------
402	! initialise
403	il_pivot=-1
404	il_perio=-1
405	il_ew =-1
406
407	! copy structure
408	tl_file=file_copy(td_file)
409
410	! open file to be used
411	CALL iom_open(tl_file)
412
413	IF( td_file%i_perio >= 0 .AND. td_file%i_perio <= 6 )THEN
414	il_perio=td_file%i_perio
415	ELSE
416	! look for attribute in file
417	il_attid=att_get_index(tl_file%t_att(:),'periodicity')
418	IF( il_attid /= 0 )THEN
419	il_perio=INT(tl_file%t_att(il_attid)%d_value(1),i4)
420	ENDIF
421	ENDIF
422
423	IF( td_file%i_ew >= 0 )THEN
424	il_ew=td_file%i_ew
425	ELSE
426	! look for attribute in file
427	il_attid=att_get_index(tl_file%t_att(:),'ew_overlap')
428	IF( il_attid /= 0 )THEN
429	il_ew=INT(tl_file%t_att(il_attid)%d_value(1),i4)
430	ENDIF
431	ENDIF
432
433	SELECT CASE(il_perio)
434	CASE(3,4)
435	il_pivot=0
436	CASE(5,6)
437	il_pivot=1
438	CASE(0,1,2)
439	il_pivot=1
440	END SELECT
441
442	IF( il_pivot < 0 .OR. il_pivot > 1 )THEN
443	! get pivot
444	il_pivot=grid_get_pivot(tl_file)
445	ENDIF
446
447	IF( il_perio < 0 .OR. il_perio > 6 )THEN
448	! get periodicity
449	il_perio=grid_get_perio(tl_file, il_pivot)
450	ENDIF
451
452	IF( il_ew < 0 )THEN
453	! get periodicity
454	il_ew=grid_get_ew_overlap(tl_file)
455	ENDIF
456
457	! close
458	CALL iom_close(tl_file)
459
460	!save in file structure
461	td_file%i_ew=il_ew
462	td_file%i_pivot=il_pivot
463	td_file%i_perio=il_perio
464
465	! save in variable of file structure
466	tl_att=att_init("ew_overlap",il_ew)
467	DO ji=1,td_file%i_nvar
468	IF( td_file%t_var(ji)%t_dim(jp_I)%l_use )THEN
469	CALL var_move_att(td_file%t_var(ji),tl_att)
470	ENDIF
471	ENDDO
472
473	! clean
474	CALL file_clean(tl_file)
475	CALL att_clean(tl_att)
476
477	IF( td_file%i_perio == -1 )THEN
478	CALL logger_fatal("GRID GET INFO: can not read or compute "//&
479	& "domain periodicity from file "//TRIM(td_file%c_name)//"."//&
480	& " you have to inform periodicity in namelist.")
481	ENDIF
482
483	END SUBROUTINE grid__get_info_file
484	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
485	SUBROUTINE grid__get_info_mpp(td_mpp)
486	!-------------------------------------------------------------------
487	!> @brief This subroutine get information about global domain, given mpp
488	!> structure.
489	!>
490	!> @details
491	!> open edge files then:
492	!> - compute NEMO pivot point
493	!> - compute NEMO periodicity
494	!> - compute East West overlap
495	!>
496	!> @note need all processor files
497	!> @author J.Paul
498	!> @date October, 2014 - Initial Version
499	!>
500	!> @param[in] td_mpp mpp structure
501	!-------------------------------------------------------------------
502
503	IMPLICIT NONE
504
505	! Argument
506	TYPE(TMPP) , INTENT(INOUT) :: td_mpp
507
508	! local variable
509	INTEGER(i4) :: il_ew
510	INTEGER(i4) :: il_pivot
511	INTEGER(i4) :: il_perio
512	INTEGER(i4) :: il_attid
513
514	TYPE(TATT) :: tl_att
515
516	TYPE(TMPP) :: tl_mpp
517
518	! loop indices
519	INTEGER(i4) :: ji
520	INTEGER(i4) :: jj
521	!----------------------------------------------------------------
522	! initialise
523	il_pivot=-1
524	il_perio=-1
525	il_ew =-1
526
527	CALL logger_info("GRID GET INFO: look for "//TRIM(td_mpp%c_name))
528	! copy structure
529	tl_mpp=mpp_copy(td_mpp)
530	! select edge files
531	CALL mpp_get_contour(tl_mpp)
532	! open mpp file to be used
533	CALL iom_mpp_open(tl_mpp)
534
535	IF( td_mpp%i_perio >= 0 .AND. td_mpp%i_perio <= 6 )THEN
536	il_perio=td_mpp%i_perio
537	ELSE
538	! look for attribute in mpp files
539	il_attid=att_get_index(tl_mpp%t_proc(1)%t_att(:),'periodicity')
540	IF( il_attid /= 0 )THEN
541	il_perio=INT(tl_mpp%t_proc(1)%t_att(il_attid)%d_value(1),i4)
542	ENDIF
543	ENDIF
544
545	IF( td_mpp%i_ew >= 0 )THEN
546	il_ew=td_mpp%i_ew
547	ELSE
548	! look for attribute in mpp files
549	il_attid=att_get_index(tl_mpp%t_proc(1)%t_att(:),'ew_overlap')
550	IF( il_attid /= 0 )THEN
551	il_ew=INT(tl_mpp%t_proc(1)%t_att(il_attid)%d_value(1),i4)
552	ENDIF
553	ENDIF
554
555	CALL logger_info("GRID GET INFO: perio "//TRIM(fct_str(il_perio)))
556
557	SELECT CASE(il_perio)
558	CASE(3,4)
559	il_pivot=1
560	CASE(5,6)
561	il_pivot=0
562	CASE(0,1,2)
563	il_pivot=1
564	END SELECT
565
566	IF( il_pivot < 0 .OR. il_pivot > 1 )THEN
567	! get pivot
568	CALL logger_info("GRID GET INFO: look for pivot ")
569	il_pivot=grid_get_pivot(tl_mpp)
570	ENDIF
571
572	IF( il_perio < 0 .OR. il_perio > 6 )THEN
573	! get periodicity
574	CALL logger_info("GRID GET INFO: look for perio ")
575	il_perio=grid_get_perio(tl_mpp, il_pivot)
576	ENDIF
577
578	IF( il_ew < 0 )THEN
579	! get periodicity
580	CALL logger_info("GRID GET INFO: look for overlap ")
581	il_ew=grid_get_ew_overlap(tl_mpp)
582	ENDIF
583
584	! close
585	CALL iom_mpp_close(tl_mpp)
586
587	!save in mpp structure
588	td_mpp%i_ew=il_ew
589	td_mpp%i_pivot=il_pivot
590	td_mpp%i_perio=il_perio
591
592	! save in variable of mpp structure
593	IF( ASSOCIATED(td_mpp%t_proc) )THEN
594	tl_att=att_init("ew_overlap",il_ew)
595	DO jj=1,td_mpp%i_nproc
596	DO ji=1,td_mpp%t_proc(jj)%i_nvar
597	IF( td_mpp%t_proc(jj)%t_var(ji)%t_dim(jp_I)%l_use )THEN
598	CALL var_move_att(td_mpp%t_proc(jj)%t_var(ji),tl_att)
599	ENDIF
600	ENDDO
601	ENDDO
602	ENDIF
603
604	! clean
605	CALL mpp_clean(tl_mpp)
606	CALL att_clean(tl_att)
607
608	IF( td_mpp%i_perio == -1 )THEN
609	CALL logger_fatal("GRID GET INFO: can not read or compute "//&
610	& "domain periodicity from mpp "//TRIM(td_mpp%c_name)//"."//&
611	& " you have to inform periodicity in namelist.")
612	ENDIF
613
614	END SUBROUTINE grid__get_info_mpp
615	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
616	FUNCTION grid__get_pivot_var(td_var) &
617	& RESULT (if_pivot)
618	!-------------------------------------------------------------------
619	!> @brief
620	!> This function compute NEMO pivot point index of the input variable.
621	!> - F-point : 0
622	!> - T-point : 1
623	!>
624	!> @details
625	!> check north points of latitude grid (indices jpj to jpj-3) depending on which grid point
626	!> (T,F,U,V) variable is defined
627	!>
628	!> @note variable must be at least 2D variable, and should not be coordinate
629	!> variable (i.e lon, lat)
630	!>
631	!> @warning
632	!> - do not work with ORCA2 grid (T-point)
633	!>
634	!> @author J.Paul
635	!> @date November, 2013 - Initial version
636	!> @date September, 2014
637	!> - add dummy loop in case variable not over right point.
638	!> @date October, 2014
639	!> - work on variable structure instead of file structure
640	!>
641	!> @param[in] td_lat latitude variable structure
642	!> @param[in] td_var variable structure
643	!> @return pivot point index
644	!-------------------------------------------------------------------
645
646	IMPLICIT NONE
647
648	! Argument
649	TYPE(TVAR), INTENT(IN) :: td_var
650
651	! function
652	INTEGER(i4) :: if_pivot
653
654	! local variable
655	INTEGER(i4), DIMENSION(ip_maxdim) :: il_dim
656
657	REAL(dp) , DIMENSION(:,:,:,:), ALLOCATABLE :: dl_value
658
659	! loop indices
660	INTEGER(i4) :: jj
661	!----------------------------------------------------------------
662	! intitalise
663	if_pivot=-1
664
665	IF( .NOT. ASSOCIATED(td_var%d_value) .OR. &
666	&.NOT. ALL(td_var%t_dim(1:2)%l_use) )THEN
667	CALL logger_error("GRID GET PIVOT: can not compute pivot point"//&
668	& " with variable "//TRIM(td_var%c_name)//"."//&
669	& " no value associated or missing dimension.")
670	ELSE
671	il_dim(:)=td_var%t_dim(:)%i_len
672
673	ALLOCATE(dl_value(il_dim(1),4,1,1))
674	! extract value
675	dl_value(:,:,:,:)=td_var%d_value( 1:il_dim(1), &
676	& il_dim(2)-3:il_dim(2),&
677	& 1:1, &
678	& 1:1 )
679
680	SELECT CASE(TRIM(td_var%c_point))
681	CASE('T')
682	if_pivot=grid__get_pivot_varT(dl_value)
683	CASE('U')
684	if_pivot=grid__get_pivot_varU(dl_value)
685	CASE('V')
686	if_pivot=grid__get_pivot_varV(dl_value)
687	CASE('F')
688	if_pivot=grid__get_pivot_varF(dl_value)
689	END SELECT
690
691	! dummy loop in case variable not over right point
692	! (ex: nav_lon over U-point)
693	IF( if_pivot == -1 )THEN
694
695	! no pivot point found
696	CALL logger_warn("GRID GET PIVOT: something wrong "//&
697	& "when computing pivot point with variable "//&
698	& TRIM(td_var%c_name))
699
700	DO jj=1,ip_npoint
701	SELECT CASE(TRIM(cp_grid_point(jj)))
702	CASE('T')
703	CALL logger_debug("GRID GET PIVOT: check variable on point T")
704	if_pivot=grid__get_pivot_varT(dl_value)
705	CASE('U')
706	CALL logger_debug("GRID GET PIVOT: check variable on point U")
707	if_pivot=grid__get_pivot_varU(dl_value)
708	CASE('V')
709	CALL logger_debug("GRID GET PIVOT: check variable on point V")
710	if_pivot=grid__get_pivot_varV(dl_value)
711	CASE('F')
712	CALL logger_debug("GRID GET PIVOT: check variable on point F")
713	if_pivot=grid__get_pivot_varF(dl_value)
714	END SELECT
715
716	IF( if_pivot /= -1 )THEN
717	CALL logger_info("GRID GET PIVOT: variable "//&
718	& TRIM(td_var%c_name)//" seems to be on grid point "//&
719	& TRIM(cp_grid_point(jj)) )
720	EXIT
721	ENDIF
722
723	ENDDO
724	ENDIF
725
726	IF( if_pivot == -1 )THEN
727	CALL logger_warn("GRID GET PIVOT: not able to found pivot point. "//&
728	& "Force to use pivot point T.")
729	if_pivot = 1
730	ENDIF
731
732	! clean
733	DEALLOCATE(dl_value)
734
735	ENDIF
736
737	END FUNCTION grid__get_pivot_var
738	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
739	FUNCTION grid__get_pivot_varT(dd_value) &
740	& RESULT (if_pivot)
741	!-------------------------------------------------------------------
742	!> @brief
743	!> This function compute NEMO pivot point index for variable on grid T.
744	!>
745	!> @details
746	!> - F-point : 0
747	!> - T-point : 1
748	!>
749	!> @note array of value must be only the top border of the domain.
750	!>
751	!> @author J.Paul
752	!> @date October, 2014 - Initial version
753	!>
754	!> @param[in] dd_value array of value
755	!> @return pivot point index
756	!-------------------------------------------------------------------
757
758	IMPLICIT NONE
759
760	! Argument
761	REAL(dp), DIMENSION(:,:,:,:), INTENT(IN) :: dd_value
762
763	! function
764	INTEGER(i4) :: if_pivot
765
766	! local variable
767	INTEGER(i4) :: il_midT
768	INTEGER(i4) :: il_midF
769
770	INTEGER(i4) :: it1
771	INTEGER(i4) :: it2
772	INTEGER(i4) :: jt1
773	INTEGER(i4) :: jt2
774
775	INTEGER(i4) :: if1
776	INTEGER(i4) :: if2
777	INTEGER(i4) :: jf1
778	INTEGER(i4) :: jf2
779
780	INTEGER(i4), DIMENSION(ip_maxdim) :: il_dim
781
782	LOGICAL :: ll_check
783
784	! loop indices
785	INTEGER(i4) :: ji
786	!----------------------------------------------------------------
787	! intitalise
788	if_pivot=-1
789
790	il_dim(:)=SHAPE(dd_value(:,:,:,:))
791
792	! T-point pivot !case of ORCA2, ORCA025, ORCA12 grid
793	jt1=4 ; jt2=2
794	il_midT=il_dim(1)/2+1
795
796	! F-point pivot !case of ORCA05 grid
797	jf1=4 ; jf2=3
798	il_midF=il_dim(1)/2
799
800	! check T-point pivot
801	DO ji=2,il_midT
802	ll_check=.TRUE.
803	it1=ji
804	it2=il_dim(1)-(ji-2)
805	IF( dd_value(it1,jt1,1,1) /= dd_value(it2,jt2,1,1) )THEN
806	ll_check=.FALSE.
807	EXIT
808	ENDIF
809	ENDDO
810
811	IF( ll_check )THEN
812	CALL logger_info("GRID GET PIVOT: T-pivot")
813	if_pivot=1
814	ELSE
815
816	! check F-point pivot
817	DO ji=1,il_midF
818	ll_check=.TRUE.
819	if1=ji
820	if2=il_dim(1)-(ji-1)
821	IF( dd_value(if1,jf1,1,1) /= dd_value(if2,jf2,1,1) )THEN
822	ll_check=.FALSE.
823	EXIT
824	ENDIF
825	ENDDO
826
827	IF( ll_check )THEN
828	CALL logger_info("GRID GET PIVOT: F-pivot")
829	if_pivot=0
830	ENDIF
831
832	ENDIF
833
834	END FUNCTION grid__get_pivot_varT
835	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
836	FUNCTION grid__get_pivot_varU(dd_value) &
837	& RESULT (if_pivot)
838	!-------------------------------------------------------------------
839	!> @brief
840	!> This function compute NEMO pivot point index for variable on grid U.
841	!>
842	!> @details
843	!> - F-point : 0
844	!> - T-point : 1
845	!>
846	!> @note array of value must be only the top border of the domain.
847	!>
848	!> @author J.Paul
849	!> @date October, 2014 - Initial version
850	!>
851	!> @param[in] dd_value array of value
852	!> @return pivot point index
853	!-------------------------------------------------------------------
854
855	IMPLICIT NONE
856
857	! Argument
858	REAL(dp), DIMENSION(:,:,:,:), INTENT(IN) :: dd_value
859
860	! function
861	INTEGER(i4) :: if_pivot
862
863	! local variable
864	INTEGER(i4) :: il_midT
865	INTEGER(i4) :: il_midF
866
867	INTEGER(i4) :: it1
868	INTEGER(i4) :: it2
869	INTEGER(i4) :: jt1
870	INTEGER(i4) :: jt2
871
872	INTEGER(i4) :: if1
873	INTEGER(i4) :: if2
874	INTEGER(i4) :: jf1
875	INTEGER(i4) :: jf2
876
877	INTEGER(i4), DIMENSION(ip_maxdim) :: il_dim
878
879	LOGICAL :: ll_check
880
881	! loop indices
882	INTEGER(i4) :: ji
883	!----------------------------------------------------------------
884	! intitalise
885	if_pivot=-1
886
887	il_dim(:)=SHAPE(dd_value(:,:,:,:))
888
889	! T-point pivot !case of ORCA2, ORCA025, ORCA12 grid
890	jt1=4 ; jt2=2
891	il_midT=il_dim(1)/2+1
892
893	! F-point pivot !case of ORCA05 grid
894	jf1=4 ; jf2=3
895	il_midF=il_dim(1)/2
896
897	! check T-point pivot
898	DO ji=1,il_midT
899	ll_check=.TRUE.
900	it1=ji
901	it2=il_dim(1)-(ji-2)
902	IF( dd_value(it1,jt1,1,1) /= dd_value(it2-1,jt2,1,1) )THEN
903	ll_check=.FALSE.
904	EXIT
905	ENDIF
906	ENDDO
907
908	IF( ll_check )THEN
909	CALL logger_info("GRID GET PIVOT: T-pivot")
910	if_pivot=1
911	ELSE
912
913	! check F-point pivot
914	DO ji=1,il_midF
915	ll_check=.TRUE.
916	if1=ji
917	if2=il_dim(1)-(ji-1)
918	IF( dd_value(if1,jf1,1,1) /= dd_value(if2-1,jf2,1,1) )THEN
919	ll_check=.FALSE.
920	EXIT
921	ENDIF
922	ENDDO
923
924	IF( ll_check )THEN
925	CALL logger_info("GRID GET PIVOT: F-pivot")
926	if_pivot=0
927	ENDIF
928
929	ENDIF
930
931	END FUNCTION grid__get_pivot_varU
932	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
933	FUNCTION grid__get_pivot_varV(dd_value) &
934	& RESULT (if_pivot)
935	!-------------------------------------------------------------------
936	!> @brief
937	!> This function compute NEMO pivot point index for variable on grid V.
938	!>
939	!> @details
940	!> - F-point : 0
941	!> - T-point : 1
942	!>
943	!> @note array of value must be only the top border of the domain.
944	!>
945	!> @author J.Paul
946	!> @date October, 2014 - Initial version
947	!>
948	!> @param[in] dd_value array of value
949	!> @return pivot point index
950	!-------------------------------------------------------------------
951
952	IMPLICIT NONE
953
954	! Argument
955	REAL(dp), DIMENSION(:,:,:,:), INTENT(IN) :: dd_value
956
957	! function
958	INTEGER(i4) :: if_pivot
959
960	! local variable
961	INTEGER(i4) :: il_midT
962	INTEGER(i4) :: il_midF
963
964	INTEGER(i4) :: it1
965	INTEGER(i4) :: it2
966	INTEGER(i4) :: jt1
967	INTEGER(i4) :: jt2
968
969	INTEGER(i4) :: if1
970	INTEGER(i4) :: if2
971	INTEGER(i4) :: jf1
972	INTEGER(i4) :: jf2
973
974	INTEGER(i4), DIMENSION(ip_maxdim) :: il_dim
975
976	LOGICAL :: ll_check
977
978	! loop indices
979	INTEGER(i4) :: ji
980	!----------------------------------------------------------------
981	! intitalise
982	if_pivot=-1
983
984	il_dim(:)=SHAPE(dd_value(:,:,:,:))
985
986	! T-point pivot !case of ORCA2, ORCA025, ORCA12 grid
987	jt1=4 ; jt2=2
988	il_midT=il_dim(1)/2+1
989
990	! F-point pivot !case of ORCA05 grid
991	jf1=4 ; jf2=3
992	il_midF=il_dim(1)/2
993
994	! check T-point pivot
995	DO ji=2,il_midT
996	ll_check=.TRUE.
997	it1=ji
998	it2=il_dim(1)-(ji-2)
999	IF( dd_value(it1,jt1,1,1) /= dd_value(it2,jt2-1,1,1) )THEN
1000	ll_check=.FALSE.
1001	EXIT
1002	ENDIF
1003	ENDDO
1004
1005	IF( ll_check )THEN
1006	CALL logger_info("GRID GET PIVOT: T-pivot")
1007	if_pivot=1
1008	ELSE
1009
1010	! check F-point pivot
1011	DO ji=1,il_midF
1012	ll_check=.TRUE.
1013	if1=ji
1014	if2=il_dim(1)-(ji-1)
1015	IF( dd_value(if1,jf1,1,1) /= dd_value(if2,jf2-1,1,1) )THEN
1016	ll_check=.FALSE.
1017	EXIT
1018	ENDIF
1019	ENDDO
1020
1021	IF( ll_check )THEN
1022	CALL logger_info("GRID GET PIVOT: F-pivot")
1023	if_pivot=0
1024	ENDIF
1025
1026	ENDIF
1027
1028	END FUNCTION grid__get_pivot_varV
1029	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1030	FUNCTION grid__get_pivot_varF(dd_value) &
1031	& RESULT (if_pivot)
1032	!-------------------------------------------------------------------
1033	!> @brief
1034	!> This function compute NEMO pivot point index for variable on grid F.
1035	!>
1036	!> @details
1037	!> - F-point : 0
1038	!> - T-point : 1
1039	!>
1040	!> @note array of value must be only the top border of the domain.
1041	!>
1042	!> @author J.Paul
1043	!> @date October, 2014 - Initial version
1044	!>
1045	!> @param[in] dd_value array of value
1046	!> @return pivot point index
1047	!-------------------------------------------------------------------
1048
1049	IMPLICIT NONE
1050
1051	! Argument
1052	REAL(dp), DIMENSION(:,:,:,:), INTENT(IN) :: dd_value
1053
1054	! function
1055	INTEGER(i4) :: if_pivot
1056
1057	! local variable
1058	INTEGER(i4) :: il_midT
1059	INTEGER(i4) :: il_midF
1060
1061	INTEGER(i4) :: it1
1062	INTEGER(i4) :: it2
1063	INTEGER(i4) :: jt1
1064	INTEGER(i4) :: jt2
1065
1066	INTEGER(i4) :: if1
1067	INTEGER(i4) :: if2
1068	INTEGER(i4) :: jf1
1069	INTEGER(i4) :: jf2
1070
1071	INTEGER(i4), DIMENSION(ip_maxdim) :: il_dim
1072
1073	LOGICAL :: ll_check
1074
1075	! loop indices
1076	INTEGER(i4) :: ji
1077	!----------------------------------------------------------------
1078	! intitalise
1079	if_pivot=-1
1080
1081	il_dim(:)=SHAPE(dd_value(:,:,:,:))
1082
1083	! T-point pivot !case of ORCA2, ORCA025, ORCA12 grid
1084	jt1=4 ; jt2=2
1085	il_midT=il_dim(1)/2+1
1086
1087	! F-point pivot !case of ORCA05 grid
1088	jf1=4 ; jf2=3
1089	il_midF=il_dim(1)/2
1090
1091	! check T-point pivot
1092	DO ji=1,il_midT
1093	ll_check=.TRUE.
1094	it1=ji
1095	it2=il_dim(1)-(ji-2)
1096	IF( dd_value(it1,jt1,1,1) /= dd_value(it2-1,jt2-1,1,1) )THEN
1097	ll_check=.FALSE.
1098	EXIT
1099	ENDIF
1100	ENDDO
1101
1102	IF( ll_check )THEN
1103	CALL logger_info("GRID GET PIVOT: T-pivot")
1104	if_pivot=1
1105	ELSE
1106
1107	! check F-point pivot
1108	DO ji=1,il_midF
1109	ll_check=.TRUE.
1110	if1=ji
1111	if2=il_dim(1)-(ji-1)
1112	IF( dd_value(if1,jf1,1,1) /= dd_value(if2-1,jf2-1,1,1) )THEN
1113	ll_check=.FALSE.
1114	EXIT
1115	ENDIF
1116	ENDDO
1117
1118	IF( ll_check )THEN
1119	CALL logger_info("GRID GET PIVOT: F-pivot")
1120	if_pivot=0
1121	ENDIF
1122
1123	ENDIF
1124
1125	END FUNCTION grid__get_pivot_varF
1126	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1127	FUNCTION grid__get_pivot_file(td_file) &
1128	& RESULT (if_pivot)
1129	!-------------------------------------------------------------------
1130	!> @brief
1131	!> This function compute NEMO pivot point index from input file variable.
1132	!> - F-point : 0
1133	!> - T-point : 1
1134	!>
1135	!> @details
1136	!> check north points symmetry of a 2D variable (indices jpj to jpj-3), depending on which grid point
1137	!> (T,F,U,V) variable is defined.
1138	!>
1139	!> @warning
1140	!> - do not work with ORCA2 grid (T-point)
1141	!>
1142	!> @author J.Paul
1143	!> @date Ocotber, 2014 - Initial version
1144	!> @date August, 2017
1145	!> - if can't find latitude variable, assume there is a north fold
1146	!> - do not use latitude variable to get pivot (to avoid mistake with regular
1147	!> grid)
1148	!>
1149	!> @param[in] td_file file structure
1150	!> @return pivot point index
1151	!-------------------------------------------------------------------
1152
1153	IMPLICIT NONE
1154
1155	! Argument
1156	TYPE(TFILE), INTENT(IN) :: td_file
1157
1158	! function
1159	INTEGER(i4) :: if_pivot
1160
1161	! local variable
1162	INTEGER(i4) :: il_varid
1163	INTEGER(i4), DIMENSION(ip_maxdim) :: il_dim
1164
1165	LOGICAL :: ll_north
1166
1167	TYPE(TVAR) :: tl_var
1168	TYPE(TVAR) :: tl_lat
1169
1170	! loop indices
1171	INTEGER(i4) :: ji
1172	!----------------------------------------------------------------
1173	! intitalise
1174	if_pivot=-1
1175
1176	! look for north fold
1177	il_varid=var_get_index(td_file%t_var(:), 'latitude')
1178	IF( il_varid == 0 )THEN
1179	CALL logger_error("GRID GET PIVOT: no variable with name "//&
1180	& "or standard name latitude in file structure "//&
1181	& TRIM(td_file%c_name)//". Assume there is north fold and "//&
1182	& "try to get pivot point")
1183
1184	ll_north=.TRUE.
1185	ELSE
1186	IF( ASSOCIATED(td_file%t_var(il_varid)%d_value) )THEN
1187	tl_lat=var_copy(td_file%t_var(il_varid))
1188	ELSE
1189	tl_lat=iom_read_var(td_file, 'latitude')
1190	ENDIF
1191
1192	ll_north=grid_is_north_fold(tl_lat)
1193	! clean
1194	CALL var_clean(tl_lat)
1195	ENDIF
1196
1197	IF( ll_north )THEN
1198	! look for suitable variable
1199	DO ji=1,td_file%i_nvar
1200	IF( .NOT. ALL(td_file%t_var(ji)%t_dim(1:2)%l_use) ) CYCLE
1201
1202	IF( ASSOCIATED(td_file%t_var(ji)%d_value) )THEN
1203	tl_var=var_copy(td_file%t_var(ji))
1204	ELSE
1205	il_dim(:)=td_file%t_var(ji)%t_dim(:)%i_len
1206	tl_var=iom_read_var(td_file, &
1207	& td_file%t_var(ji)%c_name, &
1208	& id_start=(/1,il_dim(2)-3,1,1/), &
1209	& id_count=(/il_dim(1),4,1,1/) )
1210	ENDIF
1211	ENDDO
1212
1213	IF( ASSOCIATED(tl_var%d_value) )THEN
1214
1215	if_pivot=grid_get_pivot(tl_var)
1216
1217	ENDIF
1218
1219	! clean
1220	CALL var_clean(tl_var)
1221	ELSE
1222	CALL logger_warn("GRID GET PIVOT: no north fold. force to use T-PIVOT")
1223	if_pivot=1
1224	ENDIF
1225
1226	END FUNCTION grid__get_pivot_file
1227	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1228	FUNCTION grid__get_pivot_mpp(td_mpp) &
1229	& RESULT (if_pivot)
1230	!-------------------------------------------------------------------
1231	!> @brief
1232	!> This function compute NEMO pivot point index from input mpp variable.
1233	!> - F-point : 0
1234	!> - T-point : 1
1235	!>
1236	!> @details
1237	!> check north points symmetry of a 2D variable (indices jpj to jpj-3), depending
1238	!> on which grid point (T,F,U,V) variable is defined.
1239	!>
1240	!> @warning
1241	!> - do not work with ORCA2 grid (T-point)
1242	!>
1243	!> @author J.Paul
1244	!> @date October, 2014 - Initial version
1245	!> @date August, 2017
1246	!> - if can't find latitude variable, assume there is a north fold
1247	!> - do not use latitude variable to get pivot (to avoid mistake with regular
1248	!> grid)
1249	!>
1250	!> @param[in] td_mpp mpp file structure
1251	!> @return pivot point index
1252	!-------------------------------------------------------------------
1253
1254	IMPLICIT NONE
1255
1256	! Argument
1257	TYPE(TMPP), INTENT(IN) :: td_mpp
1258
1259	! function
1260	INTEGER(i4) :: if_pivot
1261
1262	! local variable
1263	INTEGER(i4) :: il_varid
1264	INTEGER(i4), DIMENSION(ip_maxdim) :: il_dim
1265
1266	LOGICAL :: ll_north
1267
1268	TYPE(TVAR) :: tl_var
1269	TYPE(TVAR) :: tl_lat
1270
1271	! loop indices
1272	INTEGER(i4) :: ji
1273	!----------------------------------------------------------------
1274	! intitalise
1275	if_pivot=-1
1276
1277	! look for north fold
1278	il_varid=var_get_index(td_mpp%t_proc(1)%t_var(:), 'latitude')
1279	IF( il_varid == 0 )THEN
1280	CALL logger_error("GRID GET PIVOT: no variable with name "//&
1281	& "or standard name latitude in mpp structure "//&
1282	& TRIM(td_mpp%c_name)//". Assume there is north fold and "//&
1283	& "try to get pivot point")
1284
1285	ll_north=.TRUE.
1286	ELSE
1287	IF( ASSOCIATED(td_mpp%t_proc(1)%t_var(il_varid)%d_value) )THEN
1288	!
1289	tl_lat=mpp_recombine_var(td_mpp, 'latitude')
1290	ELSE
1291	tl_lat=iom_mpp_read_var(td_mpp, 'latitude')
1292	ENDIF
1293
1294	ll_north=grid_is_north_fold(tl_lat)
1295	! clean
1296	CALL var_clean(tl_lat)
1297	ENDIF
1298
1299	IF( ll_north )THEN
1300
1301	! look for suitable variable
1302	DO ji=1,td_mpp%t_proc(1)%i_nvar
1303	IF(.NOT. ALL(td_mpp%t_proc(1)%t_var(ji)%t_dim(1:2)%l_use)) CYCLE
1304
1305	IF( ASSOCIATED(td_mpp%t_proc(1)%t_var(ji)%d_value) )THEN
1306	CALL logger_debug("GRID GET PIVOT: mpp_recombine_var"//&
1307	& TRIM(td_mpp%t_proc(1)%t_var(ji)%c_name))
1308	tl_var=mpp_recombine_var(td_mpp, &
1309	& TRIM(td_mpp%t_proc(1)%t_var(ji)%c_name))
1310	ELSE
1311	CALL logger_debug("GRID GET PIVOT: iom_mpp_read_var "//&
1312	& TRIM(td_mpp%t_proc(1)%t_var(ji)%c_name))
1313	il_dim(:)=td_mpp%t_dim(:)%i_len
1314
1315	! read variable
1316	tl_var=iom_mpp_read_var(td_mpp, &
1317	& td_mpp%t_proc(1)%t_var(ji)%c_name, &
1318	& id_start=(/1,il_dim(2)-3,1,1/), &
1319	& id_count=(/il_dim(1),4,1,1/) )
1320	ENDIF
1321	EXIT
1322	ENDDO
1323
1324	IF( ASSOCIATED(tl_var%d_value) )THEN
1325
1326	if_pivot=grid_get_pivot(tl_var)
1327
1328	ELSE
1329	CALL logger_warn("GRID GET PIVOT: force to use T-PIVOT")
1330	if_pivot=1
1331	ENDIF
1332
1333	! clean
1334	CALL var_clean(tl_var)
1335	ELSE
1336	CALL logger_warn("GRID GET PIVOT: no north fold. force to use T-PIVOT")
1337	if_pivot=1
1338	ENDIF
1339
1340	END FUNCTION grid__get_pivot_mpp
1341	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1342	FUNCTION grid__get_perio_var(td_var, id_pivot) &
1343	& RESULT (if_perio)
1344	!-------------------------------------------------------------------
1345	!> @brief
1346	!> This subroutine search NEMO periodicity index given variable structure and
1347	!> pivot point index.
1348	!> @details
1349	!> The variable must be on T point.
1350	!>
1351	!> 0: closed boundaries
1352	!> 1: cyclic east-west boundary
1353	!> 2: symmetric boundary condition across the equator
1354	!> 3: North fold boundary (with a T-point pivot)
1355	!> 4: North fold boundary (with a T-point pivot) and cyclic east-west boundary
1356	!> 5: North fold boundary (with a F-point pivot)
1357	!> 6: North fold boundary (with a F-point pivot) and cyclic east-west boundary
1358	!>
1359	!> @warning pivot point should have been computed before run this script. see grid_get_pivot.
1360	!>
1361	!> @author J.Paul
1362	!> @date November, 2013 - Initial version
1363	!> @date October, 2014
1364	!> - work on variable structure instead of file structure
1365	!>
1366	!> @param[in] td_var variable structure
1367	!> @param[in] id_pivot pivot point index
1368	!-------------------------------------------------------------------
1369
1370	IMPLICIT NONE
1371
1372	! Argument
1373	TYPE(TVAR) , INTENT(IN) :: td_var
1374	INTEGER(i4), INTENT(IN) :: id_pivot
1375
1376	! function
1377	INTEGER(i4) :: if_perio
1378
1379	! local variable
1380	INTEGER(i4), DIMENSION(ip_maxdim) :: il_dim
1381
1382	! loop indices
1383	!----------------------------------------------------------------
1384	! intitalise
1385	if_perio=-1
1386
1387	IF( id_pivot < 0 .OR. id_pivot > 1 )THEN
1388	CALL logger_error("GRID GET PERIO: invalid pivot point index. "//&
1389	& "you should use grid_get_pivot to compute it")
1390	ENDIF
1391
1392	IF( .NOT. ASSOCIATED(td_var%d_value) .OR. &
1393	& .NOT. ALL(td_var%t_dim(1:2)%l_use) )THEN
1394	CALL logger_error("GRID GET PERIO: can not compute periodicity"//&
1395	& " with variable "//TRIM(td_var%c_name)//"."//&
1396	& " no value associated or missing dimension.")
1397	ELSE
1398
1399	il_dim(:)=td_var%t_dim(:)%i_len
1400
1401	CALL logger_debug("GRID GET PERIO: use variable "//TRIM(td_var%c_name))
1402	CALL logger_debug("GRID GET PERIO: fill value "//TRIM(fct_str(td_var%d_fill)))
1403	CALL logger_debug("GRID GET PERIO: first value "//TRIM(fct_str(td_var%d_value(1,1,1,1))))
1404
1405	IF(ALL(td_var%d_value( 1 , : ,1,1)/=td_var%d_fill).AND.&
1406	& ALL(td_var%d_value(il_dim(1), : ,1,1)/=td_var%d_fill).AND.&
1407	& ALL(td_var%d_value( : , 1 ,1,1)/=td_var%d_fill).AND.&
1408	& ALL(td_var%d_value( : ,il_dim(2),1,1)/=td_var%d_fill))THEN
1409	! no boundary closed
1410	CALL logger_error("GRID GET PERIO: can't determined periodicity. "//&
1411	& "there is no boundary closed for variable "//&
1412	& TRIM(td_var%c_name) )
1413	! check pivot
1414	SELECT CASE(id_pivot)
1415	CASE(0)
1416	! F pivot
1417	CALL logger_warn("GRID GET PERIO: assume domain is global")
1418	if_perio=6
1419	CASE(1)
1420	! T pivot
1421	CALL logger_warn("GRID GET PERIO: assume domain is global")
1422	if_perio=4
1423	END SELECT
1424	ELSE
1425	if_perio=-1
1426	! check periodicity
1427	IF(ANY(td_var%d_value( 1 ,:,1,1)/=td_var%d_fill).OR.&
1428	& ANY(td_var%d_value(il_dim(1),:,1,1)/=td_var%d_fill))THEN
1429	! East-West cyclic (1,4,6)
1430
1431	IF( ANY(td_var%d_value(:, 1, 1, 1) /= td_var%d_fill) )THEN
1432	! South boundary not closed
1433
1434	CALL logger_debug("GRID GET PERIO: East_West cyclic")
1435	CALL logger_debug("GRID GET PERIO: South boundary not closed")
1436	CALL logger_error("GRID GET PERIO: should have been an "//&
1437	& "impossible case")
1438
1439	ELSE
1440	! South boundary closed (1,4,6)
1441	CALL logger_info("GRID GET PERIO: South boundary closed")
1442
1443	IF(ANY(td_var%d_value(:,il_dim(2),1,1)/=td_var%d_fill) )THEN
1444	! North boundary not closed (4,6)
1445	CALL logger_info("GRID GET PERIO: North boundary not closed")
1446	! check pivot
1447	SELECT CASE(id_pivot)
1448	CASE(0)
1449	! F pivot
1450	if_perio=6
1451	CASE(1)
1452	! T pivot
1453	if_perio=4
1454	CASE DEFAULT
1455	CALL logger_error("GRID GET PERIO: invalid pivot ")
1456	END SELECT
1457	ELSE
1458	! North boundary closed
1459	CALL logger_info("GRID GET PERIO: North boundary closed")
1460	if_perio=1 ! North and South boundaries closed
1461	ENDIF
1462
1463	ENDIF
1464
1465	ELSE
1466	! East-West boundaries closed (0,2,3,5)
1467	CALL logger_info("GRID GET PERIO: East West boundaries closed")
1468
1469	IF( ANY(td_var%d_value(:, 1, 1, 1) /= td_var%d_fill) )THEN
1470	! South boundary not closed (2)
1471	CALL logger_info("GRID GET PERIO: South boundary not closed")
1472
1473	IF(ANY(td_var%d_value(:,il_dim(2),1,1)/=td_var%d_fill))THEN
1474	! North boundary not closed
1475	CALL logger_debug("GRID GET PERIO: East West boundaries "//&
1476	& "closed")
1477	CALL logger_debug("GRID GET PERIO: South boundary not closed")
1478	CALL logger_debug("GRID GET PERIO: North boundary not closed")
1479	CALL logger_error("GRID GET PERIO: should have been "//&
1480	& "an impossible case")
1481	ELSE
1482	! North boundary closed
1483	if_perio=2 ! East-West and North boundaries closed
1484	ENDIF
1485
1486	ELSE
1487	! South boundary closed (0,3,5)
1488	CALL logger_info("GRID GET PERIO: South boundary closed")
1489
1490	IF(ANY(td_var%d_value(:,il_dim(2),1,1)/=td_var%d_fill))THEN
1491	! North boundary not closed (3,5)
1492	CALL logger_info("GRID GET PERIO: North boundary not closed")
1493	! check pivot
1494	SELECT CASE(id_pivot)
1495	CASE(0)
1496	! F pivot
1497	if_perio=5
1498	CASE(1)
1499	! T pivot
1500	if_perio=3
1501	CASE DEFAULT
1502	CALL logger_error("GRID GET PERIO: invalid pivot")
1503	END SELECT
1504	ELSE
1505	! North boundary closed
1506	CALL logger_info("GRID GET PERIO: North boundary closed")
1507	if_perio=0 ! all boundary closed
1508	ENDIF
1509
1510	ENDIF
1511
1512	ENDIF
1513
1514	ENDIF
1515
1516	ENDIF
1517
1518	END FUNCTION grid__get_perio_var
1519	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1520	FUNCTION grid__get_perio_file(td_file, id_pivot) &
1521	& RESULT (if_perio)
1522	!-------------------------------------------------------------------
1523	!> @brief
1524	!> This subroutine search NEMO periodicity index given file structure, and
1525	!> optionaly pivot point index.
1526	!> @details
1527	!> The variable used must be on T point.
1528	!>
1529	!> 0: closed boundaries
1530	!> 1: cyclic east-west boundary
1531	!> 2: symmetric boundary condition across the equator
1532	!> 3: North fold boundary (with a F-point pivot)
1533	!> 4: North fold boundary (with a F-point pivot) and cyclic east-west boundary
1534	!> 5: North fold boundary (with a T-point pivot)
1535	!> 6: North fold boundary (with a T-point pivot) and cyclic east-west boundary
1536	!>
1537	!> @warning pivot point should have been computed before run this script. see grid_get_pivot.
1538	!>
1539	!> @author J.Paul
1540	!> @date October, 2014 - Initial version
1541	!> @date August, 2017
1542	!> - read only grid boundaries to handle huge file
1543	!>
1544	!> @param[in] td_file file structure
1545	!> @param[in] id_pivot pivot point index
1546	!> @return NEMO periodicity index
1547	!-------------------------------------------------------------------
1548
1549	IMPLICIT NONE
1550
1551	! Argument
1552	TYPE(TFILE), INTENT(IN) :: td_file
1553	INTEGER(i4), INTENT(IN), OPTIONAL :: id_pivot
1554
1555	! function
1556	INTEGER(i4) :: if_perio
1557
1558	! local variable
1559	INTEGER(i4) :: il_idx
1560	INTEGER(i4) :: il_pivot
1561	INTEGER(i4), DIMENSION(ip_maxdim) :: il_dim
1562
1563	TYPE(TVAR) :: tl_var
1564	TYPE(TVAR) :: tl_tmp
1565
1566	! loop indices
1567	INTEGER(i4) :: ji
1568	!----------------------------------------------------------------
1569	!initialise
1570	if_perio=-1
1571
1572	IF(PRESENT(id_pivot) )THEN
1573	il_pivot=id_pivot
1574	ELSE
1575	il_pivot=grid_get_pivot(td_file)
1576	ENDIF
1577
1578	IF( il_pivot < 0 .OR. il_pivot > 1 )THEN
1579	CALL logger_error("GRID GET PERIO: invalid pivot point index. "//&
1580	& "you should use grid_get_pivot to compute it")
1581	ENDIF
1582
1583	! look for suitable variable
1584	il_idx=0
1585	DO ji=1,td_file%i_nvar
1586	IF( .NOT. ALL(td_file%t_var(ji)%t_dim(1:2)%l_use) ) CYCLE
1587	SELECT CASE(TRIM(fct_lower(td_file%t_var(ji)%c_stdname)) )
1588	CASE('longitude','latitude')
1589	CASE DEFAULT
1590	il_idx=ji
1591	EXIT
1592	END SELECT
1593	ENDDO
1594
1595	IF( il_idx==0 )THEN
1596
1597	CALL logger_error("GRID GET PERIO: no suitable variable to compute "//&
1598	& " periodicity in file "//TRIM(td_file%c_name))
1599
1600	ELSE
1601
1602	DO ji=1,ip_maxdim
1603	IF( td_file%t_var(il_idx)%t_dim(ji)%l_use )THEN
1604	il_dim(ji)= td_file%t_var(il_idx)%t_dim(ji)%i_len
1605	ELSE
1606	il_dim(ji)=1
1607	ENDIF
1608	ENDDO
1609
1610	! read variable (full array)
1611	!tl_var=iom_read_var(td_file, &
1612	!& td_file%t_var(il_idx)%c_name, &
1613	!& id_start=(/1,1,1,1/), &
1614	!& id_count=(/il_dim(1),il_dim(2),1,1/) )
1615
1616	! read variable (only usefull part)
1617	tl_tmp=iom_read_var(td_file, &
1618	& td_file%t_var(il_idx)%c_name, &
1619	& id_start=(/1,1,1,1/), &
1620	& id_count=(/il_dim(1),1,1,1/) )
1621
1622	! copy variable struct here, to get change done inside read_var too.
1623	tl_var=var_copy(tl_tmp,ld_value=.false.)
1624	! force dimension to be full domain dimension
1625	! (instead of proc dimension)
1626	tl_var%t_dim(:)%i_len=il_dim(:)
1627	ALLOCATE(tl_var%d_value(il_dim(jp_I), &
1628	& il_dim(jp_J), &
1629	& il_dim(jp_K), &
1630	& il_dim(jp_L)))
1631
1632	tl_var%d_value(:,1,1,1)=tl_tmp%d_value(:,1,1,1)
1633	! clean
1634	CALL var_clean(tl_tmp)
1635
1636	! read variable (only usefull part)
1637	tl_tmp=iom_read_var(td_file, &
1638	& td_file%t_var(il_idx)%c_name, &
1639	& id_start=(/1,il_dim(2),1,1/), &
1640	& id_count=(/il_dim(1),1,1,1/) )
1641
1642	tl_var%d_value(:,il_dim(2),1,1)=tl_tmp%d_value(:,1,1,1)
1643	! clean
1644	CALL var_clean(tl_tmp)
1645
1646	! read variable (only usefull part)
1647	tl_tmp=iom_read_var(td_file, &
1648	& td_file%t_var(il_idx)%c_name, &
1649	& id_start=(/1,1,1,1/), &
1650	& id_count=(/1,il_dim(2),1,1/) )
1651
1652	tl_var%d_value(1,:,1,1)=tl_tmp%d_value(1,:,1,1)
1653	! clean
1654	CALL var_clean(tl_tmp)
1655
1656	! read variable (only usefull part)
1657	tl_tmp=iom_read_var(td_file, &
1658	& td_file%t_var(il_idx)%c_name, &
1659	& id_start=(/il_dim(1),1,1,1/), &
1660	& id_count=(/1,il_dim(2),1,1/) )
1661
1662	tl_var%d_value(il_dim(1),:,1,1)=tl_tmp%d_value(1,:,1,1)
1663	! clean
1664	CALL var_clean(tl_tmp)
1665
1666	if_perio=grid_get_perio(tl_var,il_pivot)
1667
1668	! clean
1669	CALL var_clean(tl_var)
1670
1671	ENDIF
1672
1673	END FUNCTION grid__get_perio_file
1674	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1675	FUNCTION grid__get_perio_mpp(td_mpp, id_pivot) &
1676	& RESULT (if_perio)
1677	!-------------------------------------------------------------------
1678	!> @brief
1679	!> This subroutine search NEMO periodicity given mpp structure and optionaly
1680	!> pivot point index.
1681	!> @details
1682	!> The variable used must be on T point.
1683	!>
1684	!> 0: closed boundaries
1685	!> 1: cyclic east-west boundary
1686	!> 2: symmetric boundary condition across the equator
1687	!> 3: North fold boundary (with a T-point pivot)
1688	!> 4: North fold boundary (with a T-point pivot) and cyclic east-west boundary
1689	!> 5: North fold boundary (with a F-point pivot)
1690	!> 6: North fold boundary (with a F-point pivot) and cyclic east-west boundary
1691	!>
1692	!> @warning pivot point should have been computed before run this script. see grid_get_pivot.
1693	!>
1694	!> @author J.Paul
1695	!> @date October, 2014 - Initial version
1696	!> @date August, 2017
1697	!> - read only grid boundaries to handle huge file
1698	!> @date January, 2019
1699	!> - do not use silicalim, or silicamax to get pivot point
1700	!>
1701	!> @todo
1702	!> do not check silicalim, or silicamax
1703	!>
1704	!> @param[in] td_mpp mpp file structure
1705	!> @param[in] id_pivot pivot point index
1706	!> @return NEMO periodicity index
1707	!-------------------------------------------------------------------
1708
1709	IMPLICIT NONE
1710
1711	! Argument
1712	TYPE(TMPP) , INTENT(IN) :: td_mpp
1713	INTEGER(i4), INTENT(IN), OPTIONAL :: id_pivot
1714
1715	! function
1716	INTEGER(i4) :: if_perio
1717
1718	! local variable
1719	INTEGER(i4) :: il_idx
1720	INTEGER(i4) :: il_pivot
1721	INTEGER(i4), DIMENSION(ip_maxdim) :: il_dim
1722
1723	TYPE(TVAR) :: tl_var
1724	TYPE(TVAR) :: tl_tmp
1725
1726	! loop indices
1727	INTEGER(i4) :: ji
1728	!----------------------------------------------------------------
1729	! initialise
1730	if_perio=-1
1731
1732	IF(PRESENT(id_pivot) )THEN
1733	il_pivot=id_pivot
1734	ELSE
1735	il_pivot=grid_get_pivot(td_mpp)
1736	ENDIF
1737
1738	IF( il_pivot < 0 .OR. il_pivot > 1 )THEN
1739	CALL logger_error("GRID GET PERIO: invalid pivot point index. "//&
1740	& "you should use grid_get_pivot to compute it")
1741	ENDIF
1742
1743	! look for suitable variable
1744	il_idx=0
1745	DO ji=1,td_mpp%t_proc(1)%i_nvar
1746	IF( .NOT. ALL(td_mpp%t_proc(1)%t_var(ji)%t_dim(1:2)%l_use) ) CYCLE
1747	SELECT CASE(TRIM(fct_lower(td_mpp%t_proc(1)%t_var(ji)%c_stdname)) )
1748	CASE('longitude','latitude')
1749	CASE DEFAULT
1750	SELECT CASE(TRIM(fct_lower(td_mpp%t_proc(1)%t_var(ji)%c_name)))
1751	CASE('silicalim','silicamax')
1752	CASE DEFAULT
1753	il_idx=ji
1754	EXIT
1755	END SELECT
1756	END SELECT
1757	ENDDO
1758
1759	IF( il_idx==0 )THEN
1760
1761	CALL logger_error("GRID GET PERIO: no suitable variable to compute "//&
1762	& " periodicity in file "//TRIM(td_mpp%c_name))
1763
1764	ELSE
1765
1766	! full domain dimension
1767	DO ji=1,ip_maxdim
1768	IF( td_mpp%t_proc(1)%t_var(il_idx)%t_dim(ji)%l_use )THEN
1769	il_dim(ji)=td_mpp%t_dim(ji)%i_len
1770	ELSE
1771	il_dim(ji)=1
1772	ENDIF
1773	ENDDO
1774
1775	! read variable (full array)
1776	!tl_var=iom_mpp_read_var(td_mpp, td_mpp%t_proc(1)%t_var(il_idx)%c_name)
1777
1778	! read variable (only usefull part)
1779	tl_tmp=iom_mpp_read_var(td_mpp, &
1780	& td_mpp%t_proc(1)%t_var(il_idx)%c_name, &
1781	& id_start=(/1,1,1,1/), &
1782	& id_count=(/il_dim(1),1,1,1/) )
1783
1784	! copy variable struct here, to get change done inside read_var too.
1785	tl_var=var_copy(tl_tmp,ld_value=.false.)
1786	! force dimension to be full domain dimension
1787	! (instead of proc dimension)
1788	tl_var%t_dim(:)%i_len=il_dim(:)
1789	ALLOCATE(tl_var%d_value(il_dim(jp_I), &
1790	& il_dim(jp_J), &
1791	& il_dim(jp_K), &
1792	& il_dim(jp_L)))
1793
1794	tl_var%d_value(:,1,1,1)=tl_tmp%d_value(:,1,1,1)
1795	! clean
1796	CALL var_clean(tl_tmp)
1797
1798	! read variable (only usefull part)
1799	tl_tmp=iom_mpp_read_var(td_mpp, &
1800	& td_mpp%t_proc(1)%t_var(il_idx)%c_name, &
1801	& id_start=(/1,il_dim(2),1,1/), &
1802	& id_count=(/il_dim(1),1,1,1/) )
1803
1804	tl_var%d_value(:,il_dim(2),1,1)=tl_tmp%d_value(:,1,1,1)
1805	! clean
1806	CALL var_clean(tl_tmp)
1807
1808	! read variable (only usefull part)
1809	tl_tmp=iom_mpp_read_var(td_mpp, &
1810	& td_mpp%t_proc(1)%t_var(il_idx)%c_name, &
1811	& id_start=(/1,1,1,1/), &
1812	& id_count=(/1,il_dim(2),1,1/) )
1813
1814	tl_var%d_value(1,:,1,1)=tl_tmp%d_value(1,:,1,1)
1815	! clean
1816	CALL var_clean(tl_tmp)
1817
1818	! read variable (only usefull part)
1819	tl_tmp=iom_mpp_read_var(td_mpp, &
1820	& td_mpp%t_proc(1)%t_var(il_idx)%c_name, &
1821	& id_start=(/il_dim(1),1,1,1/), &
1822	& id_count=(/1,il_dim(2),1,1/) )
1823
1824	tl_var%d_value(il_dim(1),:,1,1)=tl_tmp%d_value(1,:,1,1)
1825	! clean
1826	CALL var_clean(tl_tmp)
1827
1828	if_perio=grid_get_perio(tl_var, il_pivot)
1829
1830	! clean
1831	CALL var_clean(tl_var)
1832	ENDIF
1833
1834	END FUNCTION grid__get_perio_mpp
1835	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1836	FUNCTION grid__get_ew_overlap_var(td_var) &
1837	& RESULT (if_overlap)
1838	!-------------------------------------------------------------------
1839	!> @brief This function get East-West overlap.
1840	!>
1841	!> @details
1842	!> If no East-West wrap return -1,
1843	!> else return the size of the ovarlap band.
1844	!> East-West overlap is computed comparing longitude value of the
1845	!> South part of the domain, to avoid north fold boundary.
1846	!>
1847	!>
1848	!> @author J.Paul
1849	!> @date November, 2013 - Initial Version
1850	!> @date October, 2014
1851	!> - work on mpp file structure instead of file structure
1852	!> @date October, 2016
1853	!> - check longitude as longname
1854	!>
1855	!> @param[in] td_lon longitude variable structure
1856	!> @return East West overlap
1857	!-------------------------------------------------------------------
1858
1859	IMPLICIT NONE
1860
1861	! Argument
1862	TYPE(TVAR), INTENT(INOUT) :: td_var
1863
1864	! function
1865	INTEGER(i4) :: if_overlap
1866
1867	! local variable
1868	REAL(dp), DIMENSION(:,:), ALLOCATABLE :: dl_value
1869	REAL(dp), DIMENSION(:) , ALLOCATABLE :: dl_vare
1870	REAL(dp), DIMENSION(:) , ALLOCATABLE :: dl_varw
1871
1872	REAL(dp) :: dl_delta
1873	REAL(dp) :: dl_varmax
1874	REAL(dp) :: dl_varmin
1875
1876	INTEGER(i4) :: il_east
1877	INTEGER(i4) :: il_west
1878	INTEGER(i4) :: il_jmin
1879	INTEGER(i4) :: il_jmax
1880
1881	! loop indices
1882	INTEGER(i4) :: ji
1883	!----------------------------------------------------------------
1884	! initialise
1885	if_overlap=-1
1886
1887	IF( ASSOCIATED(td_var%d_value) )THEN
1888
1889	IF( td_var%t_dim(1)%i_len > 1 )THEN
1890	il_west=1
1891	il_east=td_var%t_dim(1)%i_len
1892
1893	ALLOCATE( dl_value(td_var%t_dim(1)%i_len, &
1894	& td_var%t_dim(2)%i_len) )
1895
1896	dl_value(:,:)=td_var%d_value(:,:,1,1)
1897
1898	! we do not use jmax as dimension length due to north fold boundary
1899	IF( td_var%t_dim(2)%i_len > 1 )THEN
1900	il_jmin=1+ip_ghost
1901	il_jmax=(td_var%t_dim(2)%i_len-ip_ghost)/2
1902	ELSE
1903	il_jmin=1
1904	il_jmax=1
1905	ENDIF
1906
1907	ALLOCATE( dl_vare(il_jmax-il_jmin+1) )
1908	ALLOCATE( dl_varw(il_jmax-il_jmin+1) )
1909
1910	dl_vare(:)=dl_value(il_east,il_jmin:il_jmax)
1911	dl_varw(:)=dl_value(il_west,il_jmin:il_jmax)
1912
1913	IF( .NOT.( ALL(dl_vare(:)==td_var%d_fill) .AND. &
1914	& ALL(dl_varw(:)==td_var%d_fill) ) )THEN
1915
1916	IF( TRIM(td_var%c_stdname) == 'longitude' .OR. &
1917	& SCAN( TRIM(td_var%c_longname), 'longitude') == 0 )THEN
1918	WHERE( dl_value(:,:) > 180._dp .AND. &
1919	& dl_value(:,:) /= td_var%d_fill )
1920	dl_value(:,:)=360.-dl_value(:,:)
1921	END WHERE
1922
1923	dl_varmax=MAXVAL(dl_value(:,il_jmin:il_jmax))
1924	dl_varmin=MINVAL(dl_value(:,il_jmin:il_jmax))
1925
1926	dl_delta=(dl_varmax-dl_varmin)/td_var%t_dim(1)%i_len
1927
1928	IF( ALL(ABS(dl_vare(:)) - ABS(dl_varw(:)) == dl_delta) )THEN
1929	if_overlap=0
1930	ENDIF
1931	ENDIF
1932
1933	IF( if_overlap == -1 )THEN
1934	DO ji=0,im_max_overlap
1935
1936	IF( il_east-ji == il_west )THEN
1937	! case of small domain
1938	EXIT
1939	ELSE
1940	dl_vare(:)=dl_value(il_east-ji,il_jmin:il_jmax)
1941
1942	IF( ALL( dl_varw(:) == dl_vare(:) ) )THEN
1943	if_overlap=ji+1
1944	EXIT
1945	ENDIF
1946	ENDIF
1947
1948	ENDDO
1949	ENDIF
1950	ENDIF
1951
1952	ENDIF
1953	ELSE
1954	CALL logger_error("GRID GET EW OVERLAP: input variable standard name"//&
1955	& TRIM(td_var%c_stdname)//" can not be used to compute East West "//&
1956	& "overalp. no value associated. ")
1957	ENDIF
1958
1959	END FUNCTION grid__get_ew_overlap_var
1960	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1961	FUNCTION grid__get_ew_overlap_file(td_file) &
1962	& RESULT (if_overlap)
1963	!-------------------------------------------------------------------
1964	!> @brief This function get East-West overlap.
1965	!>
1966	!> @details
1967	!> If no East-West wrap return -1,
1968	!> else return the size of the ovarlap band.
1969	!> East-West overlap is computed comparing longitude value of the
1970	!> South part of the domain, to avoid north fold boundary.
1971	!>
1972	!> @author J.Paul
1973	!> @date October, 2014 - Initial Version
1974	!> @date October, 2016
1975	!> - check varid for longitude_T
1976	!> @date August, 2017
1977	!> - read only grid boundaries to handle huge file
1978	!>
1979	!> @param[in] td_file file structure
1980	!> @return East West overlap
1981	!-------------------------------------------------------------------
1982
1983	IMPLICIT NONE
1984
1985	! Argument
1986	TYPE(TFILE), INTENT(INOUT) :: td_file
1987
1988	! function
1989	INTEGER(i4) :: if_overlap
1990
1991	! local variable
1992	INTEGER(i4) :: il_idx
1993	INTEGER(i4), DIMENSION(ip_maxdim) :: il_dim
1994
1995	TYPE(TVAR) :: tl_var
1996	TYPE(TVAR) :: tl_tmp
1997
1998	! loop indices
1999	INTEGER(i4) :: ji
2000	INTEGER(i4) :: jj
2001	INTEGER(i4) :: i1
2002	INTEGER(i4) :: i2
2003	INTEGER(i4) :: j1
2004	INTEGER(i4) :: j2
2005	INTEGER(i4) :: ic
2006	INTEGER(i4) :: jc
2007	!----------------------------------------------------------------
2008
2009	! look for suitable variable
2010	il_idx=var_get_index(td_file%t_var(:), 'longitude', 'longitude_T')
2011	IF( il_idx == 0 )THEN
2012	DO jj=1,td_file%i_nvar
2013	IF( ALL(td_file%t_var(jj)%t_dim(1:2)%l_use) )THEN
2014	il_idx=jj
2015	EXIT
2016	ENDIF
2017	ENDDO
2018	ENDIF
2019
2020	IF( il_idx==0 )THEN
2021
2022	CALL logger_error("GRID GET EW OVERLAP: no suitable variable to compute "//&
2023	& " east west overlap in file "//TRIM(td_file%c_name))
2024
2025	ELSE
2026	! read variable (full array)
2027	!tl_var=iom_read_var(td_file, td_file%t_var(il_idx)%c_name)
2028
2029	! full domain dimension
2030	DO ji=1,ip_maxdim
2031	IF( td_file%t_var(il_idx)%t_dim(ji)%l_use )THEN
2032	il_dim(ji)= td_file%t_var(il_idx)%t_dim(ji)%i_len
2033	ELSE
2034	il_dim(ji)=1
2035	ENDIF
2036	ENDDO
2037
2038	! read variable (only usefull part)
2039	i1=1 ; j1=1
2040	i2=im_max_overlap ; j2=il_dim(jp_J)
2041	ic=i2-i1+1 ; jc=j2-j1+1
2042	tl_tmp=iom_read_var(td_file, &
2043	& td_file%t_var(il_idx)%c_name, &
2044	& id_start=(/i1,j1,1,1/), &
2045	& id_count=(/ic,jc,1,1/) )
2046
2047	! copy variable struct here, to get change done inside read_var too.
2048	tl_var=var_copy(tl_tmp,ld_value=.false.)
2049	! force dimension to be full domain dimension
2050	! (instead of proc dimension)
2051	tl_var%t_dim(:)%i_len=il_dim(:)
2052	ALLOCATE(tl_var%d_value(il_dim(jp_I), &
2053	& il_dim(jp_J), &
2054	& il_dim(jp_K), &
2055	& il_dim(jp_L)))
2056	! init array
2057	tl_var%d_value(:,:,:,:)=tl_var%d_fill
2058
2059	tl_var%d_value(i1:i2,:,1,1)=tl_tmp%d_value(:,:,1,1)
2060	! clean
2061	CALL var_clean(tl_tmp)
2062
2063	! read variable (only usefull part)
2064	i1=il_dim(jp_I)-im_max_overlap ; j1=1
2065	i2=il_dim(jp_I) ; j2=il_dim(jp_J)
2066	ic=i2-i1+1 ; jc=j2-j1+1
2067	tl_tmp=iom_read_var(td_file, &
2068	& td_file%t_var(il_idx)%c_name, &
2069	& id_start=(/i1,j1,1,1/), &
2070	& id_count=(/ic,jc,1,1/) )
2071
2072	tl_var%d_value(i1:i2,:,1,1)=tl_tmp%d_value(:,:,1,1)
2073	! clean
2074	CALL var_clean(tl_tmp)
2075
2076	if_overlap=grid_get_ew_overlap(tl_var)
2077
2078	! clean
2079	CALL var_clean(tl_var)
2080
2081	ENDIF
2082
2083	END FUNCTION grid__get_ew_overlap_file
2084	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
2085	FUNCTION grid__get_ew_overlap_mpp(td_mpp) &
2086	& RESULT (if_overlap)
2087	!-------------------------------------------------------------------
2088	!> @brief This function get East-West overlap.
2089	!>
2090	!> @details
2091	!> If no East-West wrap return -1,
2092	!> else return the size of the ovarlap band.
2093	!> East-West overlap is computed comparing longitude value of the
2094	!> South part of the domain, to avoid north fold boundary.
2095	!>
2096	!>
2097	!> @author J.Paul
2098	!> @date November, 2013 - Initial Version
2099	!> @date October, 2014
2100	!> - work on mpp file structure instead of file structure
2101	!> @date October, 2016
2102	!> - check varid for longitude_T
2103	!> @date August, 2017
2104	!> - read only grid boundaries to handle huge file
2105	!>
2106	!> @param[in] td_mpp mpp structure
2107	!> @return East West overlap
2108	!-------------------------------------------------------------------
2109
2110	IMPLICIT NONE
2111
2112	! Argument
2113	TYPE(TMPP), INTENT(INOUT) :: td_mpp
2114
2115	! function
2116	INTEGER(i4) :: if_overlap
2117
2118	! local variable
2119	INTEGER(i4) :: il_idx
2120	INTEGER(i4) :: il_ew
2121	INTEGER(i4), DIMENSION(ip_maxdim) :: il_dim
2122
2123	TYPE(TVAR) :: tl_var
2124	TYPE(TVAR) :: tl_tmp
2125
2126	! loop indices
2127	INTEGER(i4) :: ji
2128	INTEGER(i4) :: jj
2129	INTEGER(i4) :: i1
2130	INTEGER(i4) :: i2
2131	INTEGER(i4) :: j1
2132	INTEGER(i4) :: j2
2133	INTEGER(i4) :: ic
2134	INTEGER(i4) :: jc
2135	!----------------------------------------------------------------
2136
2137	! initialise
2138	if_overlap=td_mpp%i_ew
2139
2140	il_idx=var_get_index(td_mpp%t_proc(1)%t_var(:), 'longitude', 'longitude_T')
2141	IF( il_idx == 0 )THEN
2142	DO jj=1,td_mpp%t_proc(1)%i_nvar
2143	IF( ALL(td_mpp%t_proc(1)%t_var(jj)%t_dim(1:2)%l_use) )THEN
2144	il_idx=jj
2145	EXIT
2146	ENDIF
2147	ENDDO
2148	ENDIF
2149
2150	IF( il_idx==0 )THEN
2151
2152	CALL logger_error("GRID GET EW OVERLAP: no suitable variable to compute "//&
2153	& " east west overlap in mppfile "//TRIM(td_mpp%c_name))
2154
2155	ELSE
2156	! read variable (full array)
2157	!tl_var=iom_mpp_read_var(td_mpp, il_idx)
2158
2159	! full domain dimension
2160	DO ji=1,ip_maxdim
2161	IF( td_mpp%t_proc(1)%t_var(il_idx)%t_dim(ji)%l_use )THEN
2162	!il_dim(ji)=td_mpp%t_proc(1)%t_var(il_idx)%t_dim(ji)%i_len
2163	il_dim(ji)=td_mpp%t_dim(ji)%i_len
2164	ELSE
2165	il_dim(ji)=1
2166	ENDIF
2167	ENDDO
2168
2169	! read variable (only usefull part)
2170	i1=1 ; j1=1
2171	i2=im_max_overlap ; j2=il_dim(jp_J)
2172	ic=i2-i1+1 ; jc=j2-j1+1
2173	tl_tmp=iom_mpp_read_var(td_mpp, &
2174	& td_mpp%t_proc(1)%t_var(il_idx)%c_name, &
2175	& id_start=(/i1,j1,1,1/), &
2176	& id_count=(/ic,jc,1,1/) )
2177
2178	! copy variable struct here, to get change done inside read_var too.
2179	tl_var=var_copy(tl_tmp,ld_value=.false.)
2180	! force dimension to be full domain dimension
2181	! (instead of proc dimension)
2182	tl_var%t_dim(:)%i_len=il_dim(:)
2183	ALLOCATE(tl_var%d_value(il_dim(jp_I), &
2184	& il_dim(jp_J), &
2185	& il_dim(jp_K), &
2186	& il_dim(jp_L)))
2187	! init array
2188	tl_var%d_value(:,:,:,:)=tl_var%d_fill
2189
2190	tl_var%d_value(i1:i2,:,1,1)=tl_tmp%d_value(:,:,1,1)
2191	! clean
2192	CALL var_clean(tl_tmp)
2193
2194	! read variable (only usefull part)
2195	i1=il_dim(jp_I)-im_max_overlap ; j1=1
2196	i2=il_dim(jp_I) ; j2=il_dim(jp_J)
2197	ic=i2-i1+1 ; jc=j2-j1+1
2198	tl_tmp=iom_mpp_read_var(td_mpp, &
2199	& td_mpp%t_proc(1)%t_var(il_idx)%c_name, &
2200	& id_start=(/i1,j1,1,1/), &
2201	& id_count=(/ic,jc,1,1/) )
2202
2203	tl_var%d_value(i1:i2,:,1,1)=tl_tmp%d_value(:,:,1,1)
2204	! clean
2205	CALL var_clean(tl_tmp)
2206
2207	il_ew=grid_get_ew_overlap(tl_var)
2208	IF( il_ew >= 0 )THEN
2209	if_overlap=il_ew
2210	ENDIF
2211
2212	! clean
2213	CALL var_clean(tl_var)
2214
2215	ENDIF
2216
2217	END FUNCTION grid__get_ew_overlap_mpp
2218	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
2219	FUNCTION grid_is_north_fold(td_lat) &
2220	& RESULT (lf_north)
2221	!-------------------------------------------------------------------
2222	!> @brief This subroutine check if there is north fold.
2223	!>
2224	!> @details
2225	!> check if maximum latitude greater than 88°N
2226	!>
2227	!> @author J.Paul
2228	!> @date November, 2013 - Initial Version
2229	!>
2230	!> @param[in] td_lat latitude variable structure
2231	!> @return true if there is north fold
2232	!-------------------------------------------------------------------
2233
2234	IMPLICIT NONE
2235
2236	! Argument
2237	TYPE(TVAR), INTENT(IN) :: td_lat
2238
2239	! function
2240	LOGICAL :: lf_north
2241
2242	! local variable
2243	! loop indices
2244	!----------------------------------------------------------------
2245
2246	! init
2247	lf_north=.FALSE.
2248
2249	IF( .NOT. ASSOCIATED(td_lat%d_value) )THEN
2250	CALL logger_error("GRID IS NORTH FOLD: "//&
2251	& " no value associated to latitude")
2252	ELSE
2253	IF( MAXVAL(td_lat%d_value(:,:,:,:), &
2254	& td_lat%d_value(:,:,:,:)/= td_lat%d_fill) >= 88.0 )THEN
2255
2256	lf_north=.TRUE.
2257
2258	ENDIF
2259	ENDIF
2260
2261	END FUNCTION grid_is_north_fold
2262	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
2263	SUBROUTINE grid_check_dom(td_coord, id_imin, id_imax, id_jmin, id_jmax)
2264	!-------------------------------------------------------------------
2265	!> @brief This subroutine check domain validity.
2266	!>
2267	!> @details
2268	!> If maximum latitude greater than 88°N, program will stop.
2269	!> @note Not able to manage north fold for now.
2270	!>
2271	!> @author J.Paul
2272	!> @date November, 2013 - Initial Version
2273	!> @date October, 2014
2274	!> - work on mpp file structure instead of file structure
2275	!>
2276	!> @param[in] cd_coord coordinate file
2277	!> @param[in] id_imin i-direction lower left point indice
2278	!> @param[in] id_imax i-direction upper right point indice
2279	!> @param[in] id_jmin j-direction lower left point indice
2280	!> @param[in] id_jmax j-direction upper right point indice
2281	!-------------------------------------------------------------------
2282
2283	IMPLICIT NONE
2284
2285	! Argument
2286	TYPE(TMPP) , INTENT(IN) :: td_coord
2287	INTEGER(i4), INTENT(IN) :: id_imin
2288	INTEGER(i4), INTENT(IN) :: id_imax
2289	INTEGER(i4), INTENT(IN) :: id_jmin
2290	INTEGER(i4), INTENT(IN) :: id_jmax
2291
2292	! local variable
2293	CHARACTER(LEN=lc) :: cl_name
2294
2295	INTEGER(i4) :: il_ind
2296
2297	TYPE(TVAR) :: tl_var
2298
2299	TYPE(TMPP) :: tl_coord
2300
2301	TYPE(TDOM) :: tl_dom
2302	! loop indices
2303	!----------------------------------------------------------------
2304
2305	IF( id_jmin > id_jmax .OR. id_jmax == 0 )THEN
2306
2307	CALL logger_fatal("GRID CHECK DOM: invalid domain. "//&
2308	& "can not create configuration with north pole.")
2309
2310	ELSE
2311
2312	IF( id_imin == id_imax .AND. td_coord%i_ew < 0 )THEN
2313	CALL logger_fatal("GRID CHECK DOM: invalid domain."//&
2314	& " can not create east-west cyclic fine grid"//&
2315	& " inside closed coarse grid")
2316	ENDIF
2317
2318	! copy structure
2319	tl_coord=mpp_copy(td_coord)
2320
2321	! compute domain
2322	tl_dom=dom_init( tl_coord, &
2323	& id_imin, id_imax,&
2324	& id_jmin, id_jmax )
2325
2326	! open mpp files to be used
2327	CALL iom_dom_open(tl_coord, tl_dom)
2328
2329	! read variable value on domain
2330	WRITE(cl_name,*) 'latitude'
2331	il_ind=var_get_id(tl_coord%t_proc(1)%t_var(:), cl_name)
2332	IF( il_ind == 0 )THEN
2333	CALL logger_warn("GRID CHECK DOM: no variable "//&
2334	& TRIM(cl_name)//" in file "//TRIM(tl_coord%c_name)//". &
2335	& try to use latitude_T.")
2336	WRITE(cl_name,*) 'latitude_T'
2337	ENDIF
2338	tl_var=iom_dom_read_var(tl_coord,TRIM(cl_name),tl_dom)
2339
2340	! close mpp files
2341	CALL iom_dom_close(tl_coord)
2342
2343	! clean structure
2344	CALL mpp_clean(tl_coord)
2345
2346	IF( MAXVAL(tl_var%d_value(:,:,:,:), &
2347	& tl_var%d_value(:,:,:,:)/= tl_var%d_fill) >= 88.0 )THEN
2348
2349	CALL logger_debug("GRID CHECK DOM: max latitude "//&
2350	& TRIM(fct_str(MAXVAL(tl_var%d_value(:,:,:,:)))) )
2351	CALL logger_fatal("GRID CHECK DOM: invalid domain. "//&
2352	& "can not create configuration too close from north pole.")
2353
2354	ENDIF
2355
2356	! clean
2357	CALL dom_clean(tl_dom)
2358	CALL var_clean(tl_var)
2359
2360	ENDIF
2361
2362	END SUBROUTINE grid_check_dom
2363	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
2364	FUNCTION grid__get_coarse_index_ff(td_coord0, td_coord1, &
2365	& id_rho, cd_point) &
2366	& RESULT (if_idx)
2367	!-------------------------------------------------------------------
2368	!> @brief This function get closest coarse grid indices of fine grid domain.
2369	!>
2370	!> @details
2371	!> it use coarse and fine grid coordinates files.
2372	!> optionally, you could specify the array of refinment factor (default 1.)
2373	!> optionally, you could specify on which Arakawa grid point you want to
2374	!> work (default 'T')
2375	!>
2376	!> @author J.Paul
2377	!> @date November, 2013 - Initial Version
2378	!> @date September, 2014
2379	!> - use grid point to read coordinates variable.
2380	!> @date October, 2014
2381	!> - work on mpp file structure instead of file structure
2382	!> @date February, 2015
2383	!> - use longitude or latitude as standard name, if can not find
2384	!> longitude_T, latitude_T...
2385	!>
2386	!> @param[in] td_coord0 coarse grid coordinate mpp structure
2387	!> @param[in] td_coord1 fine grid coordinate mpp structure
2388	!> @param[in] id_rho array of refinment factor (default 1.)
2389	!> @param[in] cd_point Arakawa grid point (default 'T').
2390	!> @return coarse grid indices(/(/imin0, imax0/), (/jmin0, jmax0/)/)
2391	!>
2392	!-------------------------------------------------------------------
2393
2394	IMPLICIT NONE
2395
2396	! Argument
2397	TYPE(TMPP) , INTENT(IN) :: td_coord0
2398	TYPE(TMPP) , INTENT(IN) :: td_coord1
2399	INTEGER(i4) , DIMENSION(:), INTENT(IN), OPTIONAL :: id_rho
2400	CHARACTER(LEN=*) , INTENT(IN), OPTIONAL :: cd_point
2401
2402	! function
2403	INTEGER(i4) , DIMENSION(2,2) :: if_idx
2404
2405	! local variable
2406	CHARACTER(LEN= 1) :: cl_point
2407	CHARACTER(LEN=lc) :: cl_name
2408
2409	INTEGER(i4) :: il_imin0
2410	INTEGER(i4) :: il_imax0
2411	INTEGER(i4) :: il_jmin0
2412	INTEGER(i4) :: il_jmax0
2413	INTEGER(i4) :: il_ind
2414
2415	INTEGER(i4), DIMENSION(2,2) :: il_xghost0
2416	INTEGER(i4), DIMENSION(2,2) :: il_xghost1
2417
2418	INTEGER(i4), DIMENSION(:) , ALLOCATABLE :: il_rho
2419
2420	TYPE(TVAR) :: tl_lon0
2421	TYPE(TVAR) :: tl_lat0
2422	TYPE(TVAR) :: tl_lon1
2423	TYPE(TVAR) :: tl_lat1
2424
2425	TYPE(TMPP) :: tl_coord0
2426	TYPE(TMPP) :: tl_coord1
2427
2428	! loop indices
2429	!----------------------------------------------------------------
2430
2431	! init
2432	if_idx(:,:)=0
2433
2434	ALLOCATE(il_rho(ip_maxdim))
2435	il_rho(:)=1
2436	IF( PRESENT(id_rho) ) il_rho(:)=id_rho(:)
2437
2438	cl_point='T'
2439	IF( PRESENT(cd_point) ) cl_point=TRIM(fct_upper(cd_point))
2440
2441	! copy structure
2442	tl_coord0=mpp_copy(td_coord0)
2443	tl_coord1=mpp_copy(td_coord1)
2444
2445	IF( .NOT. ASSOCIATED(tl_coord0%t_proc) .OR. &
2446	& .NOT. ASSOCIATED(tl_coord1%t_proc) )THEN
2447	CALL logger_error("GRID GET COARSE INDEX: can not get coarse "//&
2448	& "grid indices. decompsition of mpp file "//TRIM(tl_coord0%c_name)//&
2449	& " and/or "//TRIM(tl_coord1%c_name)//" not defined." )
2450	ELSE
2451	! Coarse grid
2452	! get ghost cell factor on coarse grid
2453	il_xghost0(:,:)=grid_get_ghost( tl_coord0 )
2454
2455	! open mpp files
2456	CALL iom_mpp_open(tl_coord0)
2457
2458	! read coarse longitue and latitude
2459	WRITE(cl_name,*) 'longitude_'//TRIM(cl_point)
2460	il_ind=var_get_id(tl_coord0%t_proc(1)%t_var(:), cl_name)
2461	IF( il_ind == 0 )THEN
2462	CALL logger_warn("GRID GET COARSE INDEX: no variable "//&
2463	& TRIM(cl_name)//" in file "//TRIM(tl_coord0%c_name)//". &
2464	& try to use longitude.")
2465	WRITE(cl_name,*) 'longitude'
2466	ENDIF
2467	tl_lon0=iom_mpp_read_var(tl_coord0, TRIM(cl_name))
2468
2469	WRITE(cl_name,*) 'latitude_'//TRIM(cl_point)
2470	il_ind=var_get_id(tl_coord0%t_proc(1)%t_var(:), cl_name)
2471	IF( il_ind == 0 )THEN
2472	CALL logger_warn("GRID GET COARSE INDEX: no variable "//&
2473	& TRIM(cl_name)//" in file "//TRIM(tl_coord0%c_name)//". &
2474	& try to use latitude.")
2475	WRITE(cl_name,*) 'latitude'
2476	ENDIF
2477	tl_lat0=iom_mpp_read_var(tl_coord0, TRIM(cl_name))
2478
2479	CALL grid_del_ghost(tl_lon0, il_xghost0(:,:))
2480	CALL grid_del_ghost(tl_lat0, il_xghost0(:,:))
2481
2482	! close mpp files
2483	CALL iom_mpp_close(tl_coord0)
2484
2485	! Fine grid
2486
2487	! get ghost cell factor on fine grid
2488	il_xghost1(:,:)=grid_get_ghost( tl_coord1 )
2489
2490	! open mpp files
2491	CALL iom_mpp_open(tl_coord1)
2492
2493	! read fine longitue and latitude
2494	WRITE(cl_name,*) 'longitude_'//TRIM(cl_point)
2495	il_ind=var_get_id(tl_coord1%t_proc(1)%t_var(:), cl_name)
2496	IF( il_ind == 0 )THEN
2497	CALL logger_warn("GRID GET COARSE INDEX: no variable "//&
2498	& TRIM(cl_name)//" in file "//TRIM(tl_coord1%c_name)//". &
2499	& try to use longitude.")
2500	WRITE(cl_name,*) 'longitude'
2501	ENDIF
2502	tl_lon1=iom_mpp_read_var(tl_coord1, TRIM(cl_name))
2503
2504	WRITE(cl_name,*) 'latitude_'//TRIM(cl_point)
2505	il_ind=var_get_id(tl_coord1%t_proc(1)%t_var(:), cl_name)
2506	IF( il_ind == 0 )THEN
2507	CALL logger_warn("GRID GET COARSE INDEX: no variable "//&
2508	& TRIM(cl_name)//" in file "//TRIM(tl_coord1%c_name)//". &
2509	& try to use latitude.")
2510	WRITE(cl_name,*) 'latitude'
2511	ENDIF
2512	tl_lat1=iom_mpp_read_var(tl_coord1, TRIM(cl_name))
2513
2514	CALL grid_del_ghost(tl_lon1, il_xghost1(:,:))
2515	CALL grid_del_ghost(tl_lat1, il_xghost1(:,:))
2516
2517	! close mpp files
2518	CALL iom_mpp_close(tl_coord1)
2519
2520	! compute
2521	if_idx(:,:)=grid_get_coarse_index(tl_lon0,tl_lat0,&
2522	& tl_lon1,tl_lat1,&
2523	& il_rho(:) )
2524
2525	! add ghost cell to indices
2526	il_imin0=if_idx(1,1)+il_xghost0(jp_I,1)*ip_ghost
2527	il_imax0=if_idx(1,2)+il_xghost0(jp_I,1)*ip_ghost
2528
2529	il_jmin0=if_idx(2,1)+il_xghost0(jp_J,1)*ip_ghost
2530	il_jmax0=if_idx(2,2)+il_xghost0(jp_J,1)*ip_ghost
2531
2532	if_idx(jp_I,1)=il_imin0
2533	if_idx(jp_I,2)=il_imax0
2534	if_idx(jp_J,1)=il_jmin0
2535	if_idx(jp_J,2)=il_jmax0
2536
2537	CALL var_clean(tl_lon0)
2538	CALL var_clean(tl_lat0)
2539	CALL var_clean(tl_lon1)
2540	CALL var_clean(tl_lat1)
2541
2542	ENDIF
2543
2544	! clean
2545	CALL mpp_clean(tl_coord0)
2546	CALL mpp_clean(tl_coord1)
2547	DEALLOCATE(il_rho)
2548
2549	END FUNCTION grid__get_coarse_index_ff
2550	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
2551	FUNCTION grid__get_coarse_index_cf(td_lon0, td_lat0, td_coord1, &
2552	& id_rho, cd_point) &
2553	& RESULT (if_id)
2554	!-------------------------------------------------------------------
2555	!> @brief This function get closest coarse grid indices of fine grid domain.
2556	!>
2557	!> @details
2558	!> it use coarse array of longitude and latitude and fine grid coordinates file.
2559	!> optionaly, you could specify the array of refinment factor (default 1.)
2560	!> optionally, you could specify on which Arakawa grid point you want to
2561	!> work (default 'T')
2562	!>
2563	!> @author J.Paul
2564	!> @date November, 2013 - Initial Version
2565	!> @date September, 2014
2566	!> - use grid point to read coordinates variable.
2567	!> @date October, 2014
2568	!> - work on mpp file structure instead of file structure
2569	!> @date February, 2015
2570	!> - use longitude or latitude as standard name, if can not find
2571	!> longitude_T, latitude_T...
2572	!>
2573	!> @param[in] td_longitude0 coarse grid longitude
2574	!> @param[in] td_latitude0 coarse grid latitude
2575	!> @param[in] td_coord1 fine grid coordinate mpp structure
2576	!> @param[in] id_rho array of refinment factor
2577	!> @param[in] cd_point Arakawa grid point (default 'T')
2578	!> @return coarse grid indices (/(/imin0, imax0/), (/jmin0, jmax0/)/)
2579	!-------------------------------------------------------------------
2580
2581	IMPLICIT NONE
2582
2583	! Argument
2584	TYPE(TVAR ) , INTENT(IN) :: td_lon0
2585	TYPE(TVAR ) , INTENT(IN) :: td_lat0
2586	TYPE(TMPP ) , INTENT(IN) :: td_coord1
2587	INTEGER(i4) , DIMENSION(:), INTENT(IN), OPTIONAL :: id_rho
2588	CHARACTER(LEN=*) , INTENT(IN), OPTIONAL :: cd_point
2589
2590	! function
2591	INTEGER(i4), DIMENSION(2,2) :: if_id
2592
2593	! local variable
2594	CHARACTER(LEN= 1) :: cl_point
2595	CHARACTER(LEN=lc) :: cl_name
2596
2597	INTEGER(i4) :: il_ind
2598
2599	INTEGER(i4), DIMENSION(:) , ALLOCATABLE :: il_rho
2600
2601	INTEGER(i4), DIMENSION(2,2) :: il_xghost
2602
2603	TYPE(TVAR) :: tl_lon1
2604	TYPE(TVAR) :: tl_lat1
2605
2606	TYPE(TMPP) :: tl_coord1
2607
2608	! loop indices
2609	!----------------------------------------------------------------
2610
2611	! init
2612	if_id(:,:)=0
2613
2614	ALLOCATE(il_rho(ip_maxdim) )
2615	il_rho(:)=1
2616	IF( PRESENT(id_rho) ) il_rho(:)=id_rho(:)
2617
2618	! copy structure
2619	tl_coord1=mpp_copy(td_coord1)
2620
2621	cl_point='T'
2622	IF( PRESENT(cd_point) ) cl_point=TRIM(fct_upper(cd_point))
2623
2624	IF( .NOT. ASSOCIATED(tl_coord1%t_proc) )THEN
2625	CALL logger_error("GRID GET COARSE INDEX: decompsition of mpp "//&
2626	& "file "//TRIM(tl_coord1%c_name)//" not defined." )
2627
2628	ELSE IF( .NOT. ASSOCIATED(td_lon0%d_value) .OR. &
2629	& .NOT. ASSOCIATED(td_lat0%d_value) )THEN
2630
2631	CALL logger_error("GRID GET COARSE INDEX: some coarse grid"//&
2632	& " coordinate value are not associated.")
2633
2634	ELSE
2635
2636	IF( TRIM(td_lon0%c_point)/='' )THEN
2637	cl_point=TRIM(td_lon0%c_point)
2638	ELSEIF( TRIM(td_lat0%c_point)/='' )THEN
2639	cl_point=TRIM(td_lat0%c_point)
2640	ENDIF
2641
2642	! Fine grid
2643	! get ghost cell factor on fine grid
2644	il_xghost(:,:)=grid_get_ghost( tl_coord1 )
2645
2646	! open mpp files
2647	CALL iom_mpp_open(tl_coord1)
2648
2649	! read fine longitue and latitude
2650	WRITE(cl_name,*) 'longitude_'//TRIM(cl_point)
2651	il_ind=var_get_id(tl_coord1%t_proc(1)%t_var(:), cl_name)
2652	IF( il_ind == 0 )THEN
2653	CALL logger_warn("GRID GET COARSE INDEX: no variable "//&
2654	& TRIM(cl_name)//"in file "//TRIM(tl_coord1%c_name)//". &
2655	& try to use longitude.")
2656	WRITE(cl_name,*) 'longitude'
2657	ENDIF
2658	tl_lon1=iom_mpp_read_var(tl_coord1, TRIM(cl_name))
2659
2660	WRITE(cl_name,*) 'latitude_'//TRIM(cl_point)
2661	il_ind=var_get_id(tl_coord1%t_proc(1)%t_var(:), cl_name)
2662	IF( il_ind == 0 )THEN
2663	CALL logger_warn("GRID GET COARSE INDEX: no variable "//&
2664	& TRIM(cl_name)//"in file "//TRIM(tl_coord1%c_name)//". &
2665	& try to use longitude.")
2666	WRITE(cl_name,*) 'latitude'
2667	ENDIF
2668	tl_lat1=iom_mpp_read_var(tl_coord1, TRIM(cl_name))
2669
2670	CALL grid_del_ghost(tl_lon1, il_xghost(:,:))
2671	CALL grid_del_ghost(tl_lat1, il_xghost(:,:))
2672
2673	! close mpp files
2674	CALL iom_mpp_close(tl_coord1)
2675
2676	! compute
2677	if_id(:,:)=grid_get_coarse_index(td_lon0,td_lat0,&
2678	& tl_lon1,tl_lat1,&
2679	& il_rho(:), cl_point )
2680
2681	CALL var_clean(tl_lon1)
2682	CALL var_clean(tl_lat1)
2683
2684	ENDIF
2685
2686	DEALLOCATE(il_rho)
2687	CALL mpp_clean(tl_coord1)
2688
2689	END FUNCTION grid__get_coarse_index_cf
2690	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
2691	FUNCTION grid__get_coarse_index_fc(td_coord0, td_lon1, td_lat1, &
2692	& id_rho, cd_point) &
2693	& RESULT (if_idx)
2694	!-------------------------------------------------------------------
2695	!> @brief This function get closest coarse grid indices of fine grid domain.
2696	!>
2697	!> @details
2698	!> it use coarse grid coordinates file and fine grid array of longitude and latitude.
2699	!> optionaly, you could specify the array of refinment factor (default 1.)
2700	!> optionally, you could specify on which Arakawa grid point you want to
2701	!> work (default 'T')
2702	!>
2703	!> @author J.Paul
2704	!> @date November, 2013 - Initial Version
2705	!> @date September, 2014
2706	!> - use grid point to read coordinates variable.
2707	!> @date October, 2014
2708	!> - work on mpp file structure instead of file structure
2709	!> @date February, 2015
2710	!> - use longitude or latitude as standard name, if can not find
2711	!> longitude_T, latitude_T...
2712	!>
2713	!> @param[in] td_coord0 coarse grid coordinate mpp structure
2714	!> @param[in] td_lon1 fine grid longitude
2715	!> @param[in] td_lat1 fine grid latitude
2716	!> @param[in] id_rho array of refinment factor (default 1.)
2717	!> @param[in] cd_point Arakawa grid point (default 'T')
2718	!> @return coarse grid indices (/(/imin0, imax0/), (/jmin0, jmax0/)/)
2719	!-------------------------------------------------------------------
2720
2721	IMPLICIT NONE
2722
2723	! Argument
2724	TYPE(TMPP ) , INTENT(IN) :: td_coord0
2725	TYPE(TVAR ) , INTENT(IN) :: td_lon1
2726	TYPE(TVAR ) , INTENT(IN) :: td_lat1
2727	INTEGER(i4) , DIMENSION(:), INTENT(IN), OPTIONAL :: id_rho
2728	CHARACTER(LEN=*) , INTENT(IN), OPTIONAL :: cd_point
2729
2730	! function
2731	INTEGER(i4) , DIMENSION(2,2) :: if_idx
2732
2733	! local variable
2734	CHARACTER(LEN= 1) :: cl_point
2735	CHARACTER(LEN=lc) :: cl_name
2736
2737	INTEGER(i4) :: il_imin0
2738	INTEGER(i4) :: il_imax0
2739	INTEGER(i4) :: il_jmin0
2740	INTEGER(i4) :: il_jmax0
2741	INTEGER(i4) :: il_ind
2742
2743	INTEGER(i4), DIMENSION(:), ALLOCATABLE :: il_rho
2744
2745	INTEGER(i4), DIMENSION(2,2) :: il_xghost
2746
2747	TYPE(TVAR) :: tl_lon0
2748	TYPE(TVAR) :: tl_lat0
2749
2750	TYPE(TMPP) :: tl_coord0
2751
2752	! loop indices
2753	!----------------------------------------------------------------
2754
2755	! init
2756	if_idx(:,:)=0
2757
2758	ALLOCATE(il_rho(ip_maxdim))
2759	il_rho(:)=1
2760	IF( PRESENT(id_rho) ) il_rho(:)=id_rho(:)
2761
2762	cl_point='T'
2763	IF( PRESENT(cd_point) ) cl_point=TRIM(fct_upper(cd_point))
2764
2765	! copy structure
2766	tl_coord0=mpp_copy(td_coord0)
2767
2768	IF( .NOT. ASSOCIATED(tl_coord0%t_proc) )THEN
2769	CALL logger_error("GRID GET COARSE INDEX: decompsition of mpp "//&
2770	& "file "//TRIM(tl_coord0%c_name)//" not defined." )
2771
2772	ELSE IF( .NOT. ASSOCIATED(td_lon1%d_value) .OR. &
2773	& .NOT. ASSOCIATED(td_lat1%d_value) )THEN
2774
2775	CALL logger_error("GRID GET COARSE INDEX: some fine grid"//&
2776	& " coordinate value are not associated.")
2777
2778	ELSE
2779
2780	IF( TRIM(td_lon1%c_point)/='' )THEN
2781	cl_point=TRIM(td_lon1%c_point)
2782	ELSEIF( TRIM(td_lat1%c_point)/='' )THEN
2783	cl_point=TRIM(td_lat1%c_point)
2784	ENDIF
2785
2786	! get ghost cell factor on coarse grid
2787	il_xghost(:,:)=grid_get_ghost( tl_coord0 )
2788
2789	! open mpp files
2790	CALL iom_mpp_open(tl_coord0)
2791
2792	! read coarse longitue and latitude
2793	WRITE(cl_name,*) 'longitude_'//TRIM(cl_point)
2794	il_ind=var_get_id(tl_coord0%t_proc(1)%t_var(:), cl_name)
2795	IF( il_ind == 0 )THEN
2796	CALL logger_warn("GRID GET COARSE INDEX: no variable "//&
2797	& TRIM(cl_name)//"in file "//TRIM(tl_coord0%c_name)//". &
2798	& try to use longitude.")
2799	WRITE(cl_name,*) 'longitude'
2800	ENDIF
2801	tl_lon0=iom_mpp_read_var(tl_coord0, TRIM(cl_name))
2802
2803	WRITE(cl_name,*) 'latitude_'//TRIM(cl_point)
2804	il_ind=var_get_id(tl_coord0%t_proc(1)%t_var(:), cl_name)
2805	IF( il_ind == 0 )THEN
2806	CALL logger_warn("GRID GET COARSE INDEX: no variable "//&
2807	& TRIM(cl_name)//"in file "//TRIM(tl_coord0%c_name)//". &
2808	& try to use latitude.")
2809	WRITE(cl_name,*) 'latitude'
2810	ENDIF
2811	tl_lat0=iom_mpp_read_var(tl_coord0, TRIM(cl_name))
2812
2813	CALL grid_del_ghost(tl_lon0, il_xghost(:,:))
2814	CALL grid_del_ghost(tl_lat0, il_xghost(:,:))
2815
2816	! close mpp files
2817	CALL iom_mpp_close(tl_coord0)
2818
2819	if_idx(:,:)=grid_get_coarse_index(tl_lon0,tl_lat0,&
2820	& td_lon1,td_lat1,&
2821	& il_rho(:), cl_point )
2822
2823	! remove ghost cell
2824	il_imin0=if_idx(1,1)+il_xghost(jp_I,1)*ip_ghost
2825	il_imax0=if_idx(1,2)+il_xghost(jp_I,1)*ip_ghost
2826
2827	il_jmin0=if_idx(2,1)+il_xghost(jp_J,1)*ip_ghost
2828	il_jmax0=if_idx(2,2)+il_xghost(jp_J,1)*ip_ghost
2829
2830	if_idx(1,1)=il_imin0
2831	if_idx(1,2)=il_imax0
2832	if_idx(2,1)=il_jmin0
2833	if_idx(2,2)=il_jmax0
2834
2835	CALL var_clean(tl_lon0)
2836	CALL var_clean(tl_lat0)
2837
2838	ENDIF
2839
2840	CALL mpp_clean(tl_coord0)
2841	DEALLOCATE(il_rho)
2842
2843	END FUNCTION grid__get_coarse_index_fc
2844	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
2845	FUNCTION grid__get_coarse_index_cc(td_lon0, td_lat0, td_lon1, td_lat1, &
2846	& id_rho, cd_point) &
2847	& RESULT (if_idx)
2848	!-------------------------------------------------------------------
2849	!> @brief This function get closest coarse grid indices of fine grid domain.
2850	!
2851	!> @details
2852	!> it use coarse and fine grid array of longitude and latitude.
2853	!> optionaly, you could specify the array of refinment factor (default 1.)
2854	!> optionally, you could specify on which Arakawa grid point you want to
2855	!> work (default 'T')
2856	!>
2857	!> @note do not use ghost cell
2858	!>
2859	!> @author J.Paul
2860	!> @date November, 2013 - Initial Version
2861	!> @date September, 2014
2862	!> - check grid point
2863	!> - take into account EW overlap
2864	!> @date February, 2016
2865	!> - use delta (lon or lat)
2866	!> - manage cases for T,U,V or F point, with even or odd refinment
2867	!>
2868	!> @param[in] td_lon0 coarse grid longitude
2869	!> @param[in] td_lat0 coarse grid latitude
2870	!> @param[in] td_lon1 fine grid longitude
2871	!> @param[in] td_lat1 fine grid latitude
2872	!> @param[in] id_rho array of refinment factor
2873	!> @param[in] cd_point Arakawa grid point ('T','U','V','F')
2874	!> @return coarse grid indices (/(/imin0, imax0/), (/jmin0, jmax0/)/)
2875	!>
2876	!> @todo
2877	!> -check case boundary domain on overlap band
2878	!-------------------------------------------------------------------
2879
2880	IMPLICIT NONE
2881
2882	! Argument
2883	TYPE(TVAR) , INTENT(IN) :: td_lon0
2884	TYPE(TVAR) , INTENT(IN) :: td_lat0
2885	TYPE(TVAR) , INTENT(IN) :: td_lon1
2886	TYPE(TVAR) , INTENT(IN) :: td_lat1
2887	INTEGER(i4) , DIMENSION(:), INTENT(IN), OPTIONAL :: id_rho
2888	CHARACTER(LEN=*) , INTENT(IN), OPTIONAL :: cd_point
2889
2890	! function
2891	INTEGER(i4) , DIMENSION(2,2) :: if_idx
2892
2893	! local variable
2894	CHARACTER(LEN= 1) :: cl_point0
2895	CHARACTER(LEN= 1) :: cl_point1
2896
2897	LOGICAL , DIMENSION(2) :: ll_even
2898
2899	REAL(dp) :: dl_lon1
2900	REAL(dp) :: dl_dlon
2901	REAL(dp) :: dl_lat1
2902	REAL(dp) :: dl_dlat
2903
2904	INTEGER(i4) :: il_ew0
2905	INTEGER(i4) :: il_imin0
2906	INTEGER(i4) :: il_imax0
2907	INTEGER(i4) :: il_jmin0
2908	INTEGER(i4) :: il_jmax0
2909
2910	INTEGER(i4) :: il_ew1
2911	INTEGER(i4) :: il_imin1
2912	INTEGER(i4) :: il_imax1
2913	INTEGER(i4) :: il_jmin1
2914	INTEGER(i4) :: il_jmax1
2915
2916	INTEGER(i4) :: il_imin
2917	INTEGER(i4) :: il_imax
2918	INTEGER(i4) :: il_jmin
2919	INTEGER(i4) :: il_jmax
2920
2921	INTEGER(i4), DIMENSION(:), ALLOCATABLE :: il_rho
2922
2923	INTEGER(i4), DIMENSION(2) :: il_ill
2924	INTEGER(i4), DIMENSION(2) :: il_ilr
2925	INTEGER(i4), DIMENSION(2) :: il_iul
2926	INTEGER(i4), DIMENSION(2) :: il_iur
2927
2928	INTEGER(i4), DIMENSION(2,2) :: il_xghost0
2929	INTEGER(i4), DIMENSION(2,2) :: il_yghost0
2930	INTEGER(i4), DIMENSION(2,2) :: il_xghost1
2931	INTEGER(i4), DIMENSION(2,2) :: il_yghost1
2932
2933	TYPE(TVAR) :: tl_lon0
2934	TYPE(TVAR) :: tl_lat0
2935	TYPE(TVAR) :: tl_lon1
2936	TYPE(TVAR) :: tl_lat1
2937
2938	! loop indices
2939	!----------------------------------------------------------------
2940	! init
2941	if_idx(:,:)=0
2942
2943	ALLOCATE( il_rho(ip_maxdim) )
2944	il_rho(:)=1
2945	IF( PRESENT(id_rho) ) il_rho(:)=id_rho(:)
2946
2947	ll_even(:)=(/ (MOD(id_rho(jp_I),2)==0), (MOD(id_rho(jp_J),2)==0) /)
2948
2949	cl_point0='T'
2950	cl_point1='T'
2951	IF( PRESENT(cd_point) )THEN
2952	cl_point0=TRIM(fct_upper(cd_point))
2953	cl_point1=TRIM(fct_upper(cd_point))
2954	ENDIF
2955
2956	IF( .NOT. ASSOCIATED(td_lon0%d_value) .OR. &
2957	& .NOT. ASSOCIATED(td_lat0%d_value) .OR. &
2958	& .NOT. ASSOCIATED(td_lon1%d_value) .OR. &
2959	& .NOT. ASSOCIATED(td_lat1%d_value) )THEN
2960	CALL logger_error("GRID GET COARSE INDEX: some fine or coarse grid"//&
2961	& " coordinate value not associated.")
2962	ELSE
2963
2964	IF( TRIM(td_lon0%c_point)/='' )THEN
2965	cl_point0=TRIM(td_lon0%c_point)
2966	ELSEIF( TRIM(td_lat0%c_point)/='' )THEN
2967	cl_point0=TRIM(td_lat0%c_point)
2968	ENDIF
2969	IF( TRIM(td_lon1%c_point)/='' )THEN
2970	cl_point1=TRIM(td_lon1%c_point)
2971	ELSEIF( TRIM(td_lat1%c_point)/='' )THEN
2972	cl_point1=TRIM(td_lat1%c_point)
2973	ENDIF
2974	IF( cl_point0 /= cl_point1 )THEN
2975	CALL logger_error("GRID GET COARSE INDEX: fine and coarse grid"//&
2976	& " coordinate not on same grid point.")
2977	ENDIF
2978
2979	IF( grid_is_global(td_lon1, td_lat1) )THEN
2980
2981	IF( grid_is_global(td_lon0, td_lat0) )THEN
2982	CALL logger_trace("GRID GET COARSE INDEX: fine grid is global ")
2983	if_idx(:,:) = 1
2984	if_idx(:,:) = 0
2985	ELSE
2986	CALL logger_error("GRID GET COARSE INDEX: fine grid is "//&
2987	& "global, coarse grid not.")
2988	ENDIF
2989
2990	ELSE
2991
2992	il_xghost0(:,:)=grid_get_ghost( td_lon0 )
2993	il_yghost0(:,:)=grid_get_ghost( td_lat0 )
2994	IF( ANY(il_xghost0(:,:) /= il_yghost0(:,:)) )THEN
2995	CALL logger_error("GRID GET COARSE INDEX: coarse grid "//&
2996	& "coordinate do not share same ghost cell")
2997	ENDIF
2998
2999	tl_lon0=var_copy(td_lon0)
3000	tl_lat0=var_copy(td_lat0)
3001	CALL grid_del_ghost(tl_lon0, il_xghost0(:,:))
3002	CALL grid_del_ghost(tl_lat0, il_xghost0(:,:))
3003
3004	! "global" coarse grid indice
3005	il_imin0=1
3006	il_jmin0=1
3007
3008	il_imax0=tl_lon0%t_dim(1)%i_len
3009	il_jmax0=tl_lon0%t_dim(2)%i_len
3010
3011	! get east west overlap for coarse grid
3012	il_ew0=tl_lon0%i_ew
3013	IF( il_ew0 >= 0 )THEN
3014	! last point before overlap
3015	il_imax0=il_imax0-il_ew0
3016	ENDIF
3017
3018	il_xghost1(:,:)=grid_get_ghost( td_lon1 )
3019	il_yghost1(:,:)=grid_get_ghost( td_lat1 )
3020	IF( ANY(il_xghost1(:,:) /= il_yghost1(:,:)) )THEN
3021	CALL logger_error("GRID GET COARSE INDEX: fine grid "//&
3022	& "coordinate do not share same ghost cell")
3023	ENDIF
3024
3025	tl_lon1=var_copy(td_lon1)
3026	tl_lat1=var_copy(td_lat1)
3027	CALL grid_del_ghost(tl_lon1, il_xghost1(:,:))
3028	CALL grid_del_ghost(tl_lat1, il_xghost1(:,:))
3029
3030	! "global" fine grid indice
3031	il_imin1=1
3032	il_jmin1=1
3033
3034	il_imax1=tl_lon1%t_dim(1)%i_len
3035	il_jmax1=tl_lon1%t_dim(2)%i_len
3036
3037	! get east west overlap for fine grid
3038	il_ew1=tl_lon1%i_ew
3039	IF( il_ew1 >= 0 )THEN
3040	! last point before overlap
3041	il_imax1=il_imax1-il_ew1
3042	ENDIF
3043
3044	! get indices for each corner
3045	!1- search lower left corner indices
3046	dl_lon1=tl_lon1%d_value( il_imin1, il_jmin1, 1, 1 )
3047	dl_lat1=tl_lat1%d_value( il_imin1, il_jmin1, 1, 1 )
3048
3049	IF( dl_lon1 == tl_lon1%d_fill .OR. &
3050	& dl_lat1 == tl_lat1%d_fill )THEN
3051	CALL logger_error("GRID GET COARSE INDEX: lower left corner "//&
3052	& "point is FillValue. remove ghost cell "//&
3053	& "before running grid_get_coarse_index.")
3054	ENDIF
3055
3056	!!!!! i-direction !!!!!
3057	IF( ll_even(jp_I) )THEN
3058	! even
3059	SELECT CASE(TRIM(cl_point1))
3060	CASE('F','U')
3061	dl_dlon= ( tl_lon1%d_value(il_imin1+1,il_jmin1,1,1) - &
3062	& tl_lon1%d_value(il_imin1 ,il_jmin1,1,1) ) / &
3063	& 2.
3064	CASE DEFAULT
3065	dl_dlon=0
3066	END SELECT
3067	ELSE
3068	! odd
3069	dl_dlon= ( tl_lon1%d_value(il_imin1+1,il_jmin1,1,1) - &
3070	& tl_lon1%d_value(il_imin1 ,il_jmin1,1,1) ) / &
3071	& 2.
3072	ENDIF
3073
3074	!!!!! j-direction !!!!!
3075	IF( ll_even(jp_J) )THEN
3076	! even
3077	SELECT CASE(TRIM(cl_point1))
3078	CASE('F','V')
3079	dl_dlat= ( tl_lat1%d_value(il_imin1,il_jmin1+1,1,1) - &
3080	& tl_lat1%d_value(il_imin1,il_jmin1 ,1,1) ) / &
3081	& 2.
3082	CASE DEFAULT
3083	dl_dlat=0
3084	END SELECT
3085	ELSE
3086	! odd
3087	dl_dlat= ( tl_lat1%d_value(il_imin1,il_jmin1+1,1,1) - &
3088	& tl_lat1%d_value(il_imin1,il_jmin1 ,1,1) ) / &
3089	& 2.
3090	ENDIF
3091
3092	dl_lon1 = dl_lon1 + dl_dlon
3093	dl_lat1 = dl_lat1 + dl_dlat
3094
3095	! look for closest point on coarse grid
3096	il_ill(:)= grid_get_closest(tl_lon0%d_value(il_imin0:il_imax0, &
3097	& il_jmin0:il_jmax0, &
3098	& 1,1), &
3099	& tl_lat0%d_value(il_imin0:il_imax0, &
3100	& il_jmin0:il_jmax0, &
3101	& 1,1), &
3102	& dl_lon1, dl_lat1, 'll' )
3103
3104
3105	!2- search upper left corner indices
3106	dl_lon1=tl_lon1%d_value( il_imin1, il_jmax1, 1, 1 )
3107	dl_lat1=tl_lat1%d_value( il_imin1, il_jmax1, 1, 1 )
3108
3109	IF( dl_lon1 == tl_lon1%d_fill .OR. &
3110	& dl_lat1 == tl_lat1%d_fill )THEN
3111	CALL logger_error("GRID GET COARSE INDEX: upper left corner "//&
3112	& "point is FillValue. remove ghost cell "//&
3113	& "running grid_get_coarse_index.")
3114	ENDIF
3115
3116	!!!!! i-direction !!!!!
3117	IF( ll_even(jp_I) )THEN
3118	! even
3119	SELECT CASE(TRIM(cl_point1))
3120	CASE('F','U')
3121	dl_dlon= ( tl_lon1%d_value(il_imin1+1,il_jmax1,1,1) - &
3122	& tl_lon1%d_value(il_imin1 ,il_jmax1,1,1) ) / &
3123	& 2.
3124	CASE DEFAULT
3125	dl_dlon=0
3126	END SELECT
3127	ELSE
3128	! odd
3129	dl_dlon= ( tl_lon1%d_value(il_imin1+1,il_jmax1,1,1) - &
3130	& tl_lon1%d_value(il_imin1 ,il_jmax1,1,1) ) / &
3131	& 2.
3132	ENDIF
3133
3134	!!!!! j-direction !!!!!
3135	IF( ll_even(jp_J) )THEN
3136	! even
3137	SELECT CASE(TRIM(cl_point1))
3138	CASE('F','V')
3139	dl_dlat= ( tl_lat1%d_value(il_imin1,il_jmax1 ,1,1) - &
3140	& tl_lat1%d_value(il_imin1,il_jmax1-1,1,1) ) / &
3141	& 2.
3142	CASE DEFAULT
3143	dl_dlat=0
3144	END SELECT
3145	ELSE
3146	! odd
3147	dl_dlat= ( tl_lat1%d_value(il_imin1,il_jmax1 ,1,1) - &
3148	& tl_lat1%d_value(il_imin1,il_jmax1-1,1,1) ) / &
3149	& 2.
3150	ENDIF
3151
3152	dl_lon1 = dl_lon1 + dl_dlon
3153	dl_lat1 = dl_lat1 - dl_dlat
3154
3155	! look for closest point on coarse grid
3156	il_iul(:)= grid_get_closest(tl_lon0%d_value(il_imin0:il_imax0, &
3157	& il_jmin0:il_jmax0, &
3158	& 1,1), &
3159	& tl_lat0%d_value(il_imin0:il_imax0, &
3160	& il_jmin0:il_jmax0, &
3161	& 1,1), &
3162	& dl_lon1, dl_lat1, 'ul' )
3163
3164	!3- search lower right corner indices
3165	dl_lon1=tl_lon1%d_value( il_imax1, il_jmin1, 1, 1 )
3166	dl_lat1=tl_lat1%d_value( il_imax1, il_jmin1, 1, 1 )
3167
3168	IF( dl_lon1 == tl_lon1%d_fill .OR. &
3169	& dl_lat1 == tl_lat1%d_fill )THEN
3170	CALL logger_error("GRID GET COARSE INDEX: lower right corner "//&
3171	& "point is FillValue. remove ghost cell "//&
3172	& "running grid_get_coarse_index.")
3173	ENDIF
3174
3175	!!!!! i-direction !!!!!
3176	IF( ll_even(jp_I) )THEN
3177	! even
3178	SELECT CASE(TRIM(cl_point1))
3179	CASE('F','U')
3180	dl_dlon= ( tl_lon1%d_value(il_imax1 ,il_jmin1,1,1) - &
3181	& tl_lon1%d_value(il_imax1-1,il_jmin1,1,1) ) / &
3182	& 2.
3183	CASE DEFAULT
3184	dl_dlon=0
3185	END SELECT
3186	ELSE
3187	! odd
3188	dl_dlon= ( tl_lon1%d_value(il_imax1 ,il_jmin1,1,1) - &
3189	& tl_lon1%d_value(il_imax1-1,il_jmin1,1,1) ) / &
3190	& 2.
3191	ENDIF
3192
3193	!!!!! j-direction !!!!!
3194	IF( ll_even(jp_J) )THEN
3195	! even
3196	SELECT CASE(TRIM(cl_point1))
3197	CASE('F','V')
3198	dl_dlat= ( tl_lat1%d_value(il_imax1,il_jmin1+1,1,1) - &
3199	& tl_lat1%d_value(il_imax1,il_jmin1 ,1,1) ) / &
3200	& 2.
3201	CASE DEFAULT
3202	dl_dlat=0
3203	END SELECT
3204	ELSE
3205	! odd
3206	dl_dlat= ( tl_lat1%d_value(il_imax1,il_jmin1+1,1,1) - &
3207	& tl_lat1%d_value(il_imax1,il_jmin1 ,1,1) ) / &
3208	& 2.
3209	ENDIF
3210
3211	dl_lon1 = dl_lon1 - dl_dlon
3212	dl_lat1 = dl_lat1 + dl_dlat
3213
3214	! look for closest point on coarse grid
3215	il_ilr(:)= grid_get_closest(tl_lon0%d_value(il_imin0:il_imax0, &
3216	& il_jmin0:il_jmax0, &
3217	& 1,1), &
3218	& tl_lat0%d_value(il_imin0:il_imax0, &
3219	& il_jmin0:il_jmax0, &
3220	& 1,1), &
3221	& dl_lon1, dl_lat1, 'lr' )
3222
3223	!4- search upper right corner indices
3224	dl_lon1=tl_lon1%d_value( il_imax1, il_jmax1, 1, 1 )
3225	dl_lat1=tl_lat1%d_value( il_imax1, il_jmax1, 1, 1 )
3226
3227	IF( dl_lon1 == tl_lon1%d_fill .OR. &
3228	& dl_lat1 == tl_lat1%d_fill )THEN
3229	CALL logger_error("GRID GET COARSE INDEX: upper right corner "//&
3230	& "point is FillValue. remove ghost cell "//&
3231	& "before running grid_get_coarse_index.")
3232	ENDIF
3233
3234	!!!!! i-direction !!!!!
3235	IF( ll_even(jp_I) )THEN
3236	! even
3237	SELECT CASE(TRIM(cl_point1))
3238	CASE('F','U')
3239	dl_dlon= ( tl_lon1%d_value(il_imax1 ,il_jmax1,1,1) - &
3240	& tl_lon1%d_value(il_imax1-1,il_jmax1,1,1) ) / &
3241	& 2.
3242	CASE DEFAULT
3243	dl_dlon=0
3244	END SELECT
3245	ELSE
3246	! odd
3247	dl_dlon= ( tl_lon1%d_value(il_imax1 ,il_jmax1,1,1) - &
3248	& tl_lon1%d_value(il_imax1-1,il_jmax1,1,1) ) / &
3249	& 2.
3250	ENDIF
3251
3252	!!!!! j-direction !!!!!
3253	IF( ll_even(jp_J) )THEN
3254	! even
3255	SELECT CASE(TRIM(cl_point1))
3256	CASE('F','V')
3257	dl_dlat= ( tl_lat1%d_value(il_imax1,il_jmax1 ,1,1) - &
3258	& tl_lat1%d_value(il_imax1,il_jmax1-1,1,1) ) / &
3259	& 2.
3260	CASE DEFAULT
3261	dl_dlat=0
3262	END SELECT
3263	ELSE
3264	! odd
3265	dl_dlat= ( tl_lat1%d_value(il_imax1,il_jmax1 ,1,1) - &
3266	& tl_lat1%d_value(il_imax1,il_jmax1-1,1,1) ) / &
3267	& 2.
3268	ENDIF
3269
3270	dl_lon1 = dl_lon1 - dl_dlon
3271	dl_lat1 = dl_lat1 - dl_dlat
3272
3273	! look for closest point on coarse grid
3274	il_iur(:)= grid_get_closest(tl_lon0%d_value(il_imin0:il_imax0, &
3275	& il_jmin0:il_jmax0, &
3276	& 1,1), &
3277	& tl_lat0%d_value(il_imin0:il_imax0, &
3278	& il_jmin0:il_jmax0, &
3279	& 1,1), &
3280	& dl_lon1, dl_lat1, 'ur' )
3281
3282	! coarse grid indices
3283	il_imin = il_imin0-1+MIN(il_ill(1), il_iul(1))
3284	il_imax = il_imin0-1+MAX(il_ilr(1), il_iur(1))
3285
3286	IF( il_imax <= il_ew0 )THEN
3287	!il_imin = 1
3288	il_imax = tl_lon0%t_dim(1)%i_len - il_ew0 + il_imax
3289	ENDIF
3290
3291	il_jmin = il_jmin0-1+MIN(il_ill(2), il_ilr(2))
3292	il_jmax = il_jmin0-1+MAX(il_iul(2), il_iur(2))
3293
3294	! special case if east west overlap
3295	IF( il_ew1 >= 0 )THEN
3296	CALL logger_debug("GRID GET COARSE INDEX: East-West overlap "//&
3297	& "found for fine grid " )
3298
3299	il_imin = 1
3300	il_imax = tl_lon0%t_dim(1)%i_len
3301
3302	ENDIF
3303	ENDIF
3304
3305	if_idx(1,1) = il_imin
3306	if_idx(1,2) = il_imax
3307
3308	if_idx(2,1) = il_jmin
3309	if_idx(2,2) = il_jmax
3310
3311	! clean
3312	CALL var_clean(tl_lon1)
3313	CALL var_clean(tl_lat1)
3314	CALL var_clean(tl_lon0)
3315	CALL var_clean(tl_lat0)
3316	ENDIF
3317
3318	DEALLOCATE( il_rho )
3319
3320	END FUNCTION grid__get_coarse_index_cc
3321	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
3322	FUNCTION grid_is_global(td_lon, td_lat) &
3323	& RESULT (lf_global)
3324	!-------------------------------------------------------------------
3325	!> @brief This function check if grid is global or not
3326	!>
3327	!> @details
3328	!>
3329	!> @author J.Paul
3330	!> @date November, 2013 - Initial Version
3331	!>
3332	!> @param[in] td_lon longitude structure
3333	!> @param[in] td_lat latitude structure
3334	!> @return true if grid is global
3335	!-------------------------------------------------------------------
3336
3337	IMPLICIT NONE
3338
3339	! Argument
3340	TYPE(TVAR), INTENT(IN) :: td_lon
3341	TYPE(TVAR), INTENT(IN) :: td_lat
3342
3343	! function
3344	LOGICAL :: lf_global
3345
3346	! local variable
3347	INTEGER(i4) :: il_ew
3348	INTEGER(i4) :: il_south
3349	INTEGER(i4) :: il_north
3350
3351	REAL(dp) :: dl_lat_min
3352	REAL(dp) :: dl_lat_max
3353
3354	! loop indices
3355	!----------------------------------------------------------------
3356
3357	! init
3358	lf_global=.FALSE.
3359
3360	IF( ANY( td_lon%t_dim(:)%i_len /= td_lat%t_dim(:)%i_len ) )THEN
3361	CALL logger_fatal("GRID IS GLOBAL: dimension of longitude and "//&
3362	& " latitude differ")
3363	ENDIF
3364
3365	IF( .NOT. ASSOCIATED(td_lon%d_value) .OR. &
3366	& .NOT. ASSOCIATED(td_lat%d_value) )THEN
3367	CALL logger_error("GRID IS GLOBAL: no value associated to "//&
3368	& " longitude or latitude strucutre")
3369	ELSE
3370
3371	il_south=1
3372	il_north=td_lon%t_dim(2)%i_len
3373
3374	dl_lat_min=MINVAL(td_lat%d_value(:,il_south,1,1))
3375	dl_lat_max=MAXVAL(td_lat%d_value(:,il_north,1,1))
3376
3377	IF( dl_lat_min < -77.0 .AND. dl_lat_max > 89.0 )THEN
3378
3379	il_ew=td_lon%i_ew
3380	IF( il_ew >= 0 )THEN
3381
3382	lf_global=.TRUE.
3383
3384	ENDIF
3385
3386	ENDIF
3387	ENDIF
3388
3389	END FUNCTION grid_is_global
3390	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
3391	FUNCTION grid__get_closest_str( td_coord0, dd_lon1, dd_lat1, cd_pos, dd_fill ) &
3392	& RESULT(if_idx)
3393	!-------------------------------------------------------------------
3394	!> @brief This function return grid indices of the closest point
3395	!> from point (lon1,lat1)
3396	!>
3397	!> @details
3398	!>
3399	!> @note overlap band should have been already removed from coarse grid array
3400	!> of longitude and latitude, before running this function
3401	!>
3402	!> if you add cd_pos argument, you could choice to return closest point at
3403	!> - lower left (ll) of the point
3404	!> - lower right (lr) of the point
3405	!> - upper left (ul) of the point
3406	!> - upper right (ur) of the point
3407	!> - lower (lo) of the point
3408	!> - upper (up) of the point
3409	!> - left (le) of the point
3410	!> - right (ri) of the point
3411	!>
3412	!> @author J.Paul
3413	!> @date April, 2016 - Initial Version
3414	!> @date October, 2016
3415	!> - use max of zero and east-west overlap instead of east-west overlap
3416	!>
3417	!> @param[in] td_coord0 coarse grid coordinate mpp structure
3418	!> @param[in] dd_lon1 fine grid longitude
3419	!> @param[in] dd_lat1 fine grid latitude
3420	!> @param[in] cd_pos relative position of grid point from point
3421	!> @param[in] dd_fill fill value
3422	!> @return coarse grid indices of closest point of fine grid point
3423	!-------------------------------------------------------------------
3424
3425	IMPLICIT NONE
3426
3427	! Argument
3428	TYPE(TMPP ) , INTENT(IN) :: td_coord0
3429	REAL(dp), INTENT(IN) :: dd_lon1
3430	REAL(dp), INTENT(IN) :: dd_lat1
3431	CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: cd_pos
3432	REAL(dp), INTENT(IN), OPTIONAL :: dd_fill
3433
3434	! function
3435	INTEGER(i4), DIMENSION(2) :: if_idx
3436
3437	! local variable
3438	CHARACTER(LEN=lc) :: cl_point
3439	CHARACTER(LEN=lc) :: cl_name
3440
3441	INTEGER(i4) :: il_ind
3442	INTEGER(i4) :: il_ew
3443
3444	REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_lon0
3445	REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_lat0
3446
3447	TYPE(TVAR) :: tl_lon0
3448	TYPE(TVAR) :: tl_lat0
3449	TYPE(TMPP) :: tl_coord0
3450	!----------------------------------------------------------------
3451
3452	if_idx(:)=-1
3453	cl_point='T'
3454
3455	! copy structure
3456	tl_coord0=mpp_copy(td_coord0)
3457
3458	IF( .NOT. ASSOCIATED(tl_coord0%t_proc) )THEN
3459
3460	CALL logger_error("GRID GET CLOSEST: decompsition of mpp "//&
3461	& "file "//TRIM(tl_coord0%c_name)//" not defined." )
3462
3463	ELSE
3464
3465	! open mpp files
3466	CALL iom_mpp_open(tl_coord0)
3467
3468	! read coarse longitue and latitude
3469	WRITE(cl_name,*) 'longitude_'//TRIM(cl_point)
3470	il_ind=var_get_id(tl_coord0%t_proc(1)%t_var(:), cl_name)
3471	IF( il_ind == 0 )THEN
3472	CALL logger_warn("GRID GET CLOSEST: no variable "//&
3473	& TRIM(cl_name)//"in file "//TRIM(tl_coord0%c_name)//". &
3474	& try to use longitude.")
3475	WRITE(cl_name,*) 'longitude'
3476	ENDIF
3477	tl_lon0=iom_mpp_read_var(tl_coord0, TRIM(cl_name))
3478
3479	WRITE(cl_name,*) 'latitude_'//TRIM(cl_point)
3480	il_ind=var_get_id(tl_coord0%t_proc(1)%t_var(:), cl_name)
3481	IF( il_ind == 0 )THEN
3482	CALL logger_warn("GRID GET CLOSEST: no variable "//&
3483	& TRIM(cl_name)//"in file "//TRIM(tl_coord0%c_name)//". &
3484	& try to use latitude.")
3485	WRITE(cl_name,*) 'latitude'
3486	ENDIF
3487	tl_lat0=iom_mpp_read_var(tl_coord0, TRIM(cl_name))
3488
3489	! close mpp files
3490	CALL iom_mpp_close(tl_coord0)
3491
3492	il_ew=MAX(0,tl_coord0%i_ew)
3493	ALLOCATE(dl_lon0(tl_coord0%t_dim(jp_I)%i_len-il_ew, &
3494	& tl_coord0%t_dim(jp_J)%i_len) )
3495	ALLOCATE(dl_lat0(tl_coord0%t_dim(jp_I)%i_len-il_ew, &
3496	& tl_coord0%t_dim(jp_J)%i_len) )
3497
3498	dl_lon0(:,:)=tl_lon0%d_value(il_ew+1:,:,1,1)
3499	dl_lat0(:,:)=tl_lat0%d_value(il_ew+1:,:,1,1)
3500
3501	if_idx(:)=grid_get_closest( dl_lon0, dl_lat0, dd_lon1, dd_lat1, cd_pos, dd_fill )
3502
3503	DEALLOCATE(dl_lon0, dl_lat0)
3504	CALL var_clean(tl_lon0)
3505	CALL var_clean(tl_lat0)
3506	CALL mpp_clean(tl_coord0)
3507
3508	ENDIF
3509
3510	END FUNCTION grid__get_closest_str
3511	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
3512	FUNCTION grid__get_closest_arr(dd_lon0, dd_lat0, dd_lon1, dd_lat1, cd_pos, dd_fill) &
3513	& RESULT (if_idx)
3514	!-------------------------------------------------------------------
3515	!> @brief This function return grid indices of the closest point
3516	!> from point (lon1,lat1)
3517	!>
3518	!> @details
3519	!>
3520	!> @note overlap band should have been already removed from coarse grid array
3521	!> of longitude and latitude, before running this function
3522	!>
3523	!> if you add cd_pos argument, you could choice to return closest point at
3524	!> - lower left (ll) of the point
3525	!> - lower right (lr) of the point
3526	!> - upper left (ul) of the point
3527	!> - upper right (ur) of the point
3528	!> - lower (lo) of the point
3529	!> - upper (up) of the point
3530	!> - left (le) of the point
3531	!> - right (ri) of the point
3532	!>
3533	!> @author J.Paul
3534	!> @date November, 2013 - Initial Version
3535	!> @date February, 2015
3536	!> - change dichotomy method to manage ORCA grid
3537	!> @date February, 2016
3538	!> - add optional use of relative position
3539	!>
3540	!> @param[in] dd_lon0 coarse grid array of longitude
3541	!> @param[in] dd_lat0 coarse grid array of latitude
3542	!> @param[in] dd_lon1 fine grid longitude
3543	!> @param[in] dd_lat1 fine grid latitude
3544	!> @param[in] cd_pos relative position of grid point from point
3545	!> @param[in] dd_fill fill value
3546	!> @return coarse grid indices of closest point of fine grid point
3547	!-------------------------------------------------------------------
3548
3549	IMPLICIT NONE
3550
3551	! Argument
3552	REAL(dp), DIMENSION(:,:), INTENT(IN) :: dd_lon0
3553	REAL(dp), DIMENSION(:,:), INTENT(IN) :: dd_lat0
3554	REAL(dp), INTENT(IN) :: dd_lon1
3555	REAL(dp), INTENT(IN) :: dd_lat1
3556	CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: cd_pos
3557	REAL(dp), INTENT(IN), OPTIONAL :: dd_fill
3558
3559	! function
3560	INTEGER(i4), DIMENSION(2) :: if_idx
3561
3562	! local variable
3563	INTEGER(i4) :: il_iinf
3564	INTEGER(i4) :: il_imid
3565	INTEGER(i4) :: il_isup
3566	INTEGER(i4) :: il_jinf
3567	INTEGER(i4) :: il_jmid
3568	INTEGER(i4) :: il_jsup
3569	INTEGER(i4), DIMENSION(2) :: il_shape
3570	INTEGER(i4), DIMENSION(1) :: il_ind
3571
3572	LOGICAL :: ll_north
3573	LOGICAL :: ll_continue
3574
3575	REAL(dp) :: dl_lon1
3576	REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_dist
3577	REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_lon0
3578
3579	! loop indices
3580	!----------------------------------------------------------------
3581
3582	IF( ANY( SHAPE(dd_lon0(:,:)) /= SHAPE(dd_lat0(:,:)) ) )THEN
3583	CALL logger_fatal("GRID GET CLOSEST: dimension of longitude and "//&
3584	& " latitude differ")
3585	ENDIF
3586
3587	il_shape(:)=SHAPE(dd_lon0(:,:))
3588
3589	ALLOCATE( dl_lon0(il_shape(1),il_shape(2)) )
3590
3591	dl_lon0(:,:) = dd_lon0(:,:)
3592	WHERE(dd_lon0(:,:) < 0 ) dl_lon0(:,:) = dd_lon0(:,:) + 360.
3593
3594	dl_lon1 = dd_lon1
3595	IF( dd_lon1 < 0 ) dl_lon1 = dd_lon1 + 360.
3596
3597	! first, use dichotomy to reduce domain
3598	il_iinf = 1 ; il_jinf = 1
3599	il_isup = il_shape(1) ; il_jsup = il_shape(2)
3600
3601	il_shape(1)= il_isup - il_iinf + 1
3602	il_shape(2)= il_jsup - il_jinf + 1
3603
3604	ll_north=.FALSE.
3605	ll_continue=.FALSE.
3606
3607	! avoid to use fillvalue for reduce domain on first time
3608	IF( PRESENT(dd_fill) )THEN
3609	DO WHILE( ALL(dl_lon0(il_isup,:) == dd_fill) )
3610	il_isup=il_isup-1
3611	ENDDO
3612	DO WHILE( ALL(dl_lon0(il_iinf,:) == dd_fill) )
3613	il_iinf=il_iinf+1
3614	ENDDO
3615	DO WHILE( ALL(dd_lat0(:,il_jsup) == dd_fill) )
3616	il_jsup=il_jsup-1
3617	ENDDO
3618	DO WHILE( ALL(dd_lat0(:,il_jinf) == dd_fill) )
3619	il_jinf=il_jinf+1
3620	ENDDO
3621
3622	il_shape(1)= il_isup - il_iinf + 1
3623	il_shape(2)= il_jsup - il_jinf + 1
3624
3625	ENDIF
3626
3627	! special case for north ORCA grid
3628	IF( dd_lat1 > 19. .AND. dl_lon1 < 74. )THEN
3629	ll_north=.TRUE.
3630	ENDIF
3631
3632	IF( .NOT. ll_north )THEN
3633	! look for meridian 0°/360°
3634	il_jmid = il_jinf + INT(il_shape(2)/2)
3635	il_ind(:) = MAXLOC( dl_lon0(il_iinf:il_isup,il_jmid), &
3636	& dl_lon0(il_iinf:il_isup,il_jmid) <= 360._dp )
3637
3638	il_imid=il_ind(1)
3639
3640	IF( dl_lon1 == dl_lon0(il_imid,il_jmid) .AND. &
3641	& dd_lat1 == dd_lat0(il_imid,il_jmid) )THEN
3642
3643	il_iinf = il_imid ; il_isup = il_imid
3644	il_jinf = il_jmid ; il_jsup = il_jmid
3645
3646	ELSE
3647	IF( ALL(dl_lon0(il_isup,il_jinf:il_jsup) > dl_lon1 ) .AND. &
3648	& il_imid /= il_isup )THEN
3649	! 0 < lon1 < lon0(isup)
3650	! point east
3651	il_iinf = il_imid+1
3652	ll_continue=.TRUE.
3653
3654	ELSE IF( ALL(dl_lon0(il_iinf,il_jinf:il_jsup) < dl_lon1 ) .AND. &
3655	& il_imid /= il_iinf )THEN
3656	! lon0(iinf) < lon1 < 360
3657	! point west
3658	il_isup = il_imid
3659	ll_continue=.TRUE.
3660
3661	ENDIF
3662
3663	il_shape(1)= il_isup - il_iinf + 1
3664	il_shape(2)= il_jsup - il_jinf + 1
3665
3666	il_imid = il_iinf + INT(il_shape(1)/2)
3667	il_jmid = il_jinf + INT(il_shape(2)/2)
3668
3669	! exit when close enough of point
3670	IF( ANY(il_shape(:) < 10 ) ) ll_continue=.FALSE.
3671	ENDIF
3672	ENDIF
3673
3674	!
3675	DO WHILE( ll_continue .AND. .NOT. ll_north )
3676
3677	ll_continue=.FALSE.
3678	IF( dl_lon1 == dl_lon0(il_imid,il_jmid) .AND. &
3679	& dd_lat1 == dd_lat0(il_imid,il_jmid) )THEN
3680
3681	il_iinf = il_imid ; il_isup = il_imid
3682	il_jinf = il_jmid ; il_jsup = il_jmid
3683
3684	ELSE
3685	IF( ALL(dl_lon0(il_imid,il_jinf:il_jsup) < dl_lon1) )THEN
3686
3687	! point east
3688	il_iinf = il_imid
3689	ll_continue=.TRUE.
3690
3691	ELSE IF( ALL(dl_lon0(il_imid,il_jinf:il_jsup) > dl_lon1) )THEN
3692
3693	! point west
3694	il_isup = il_imid
3695	ll_continue=.TRUE.
3696
3697	ENDIF
3698
3699	IF( ALL(dd_lat0(il_iinf:il_isup,il_jmid) < dd_lat1) )THEN
3700
3701	! point north
3702	il_jinf = il_jmid
3703	ll_continue=.TRUE.
3704
3705	ELSE IF( ALL(dd_lat0(il_iinf:il_isup,il_jmid) > dd_lat1) )THEN
3706
3707	! point south
3708	il_jsup = il_jmid
3709	ll_continue=.TRUE.
3710
3711	ENDIF
3712
3713	il_shape(1)= il_isup - il_iinf + 1
3714	il_shape(2)= il_jsup - il_jinf + 1
3715
3716	il_imid = il_iinf + INT(il_shape(1)/2)
3717	il_jmid = il_jinf + INT(il_shape(2)/2)
3718
3719	! exit when close enough of point
3720	IF( ANY(il_shape(:) < 10 ) ) ll_continue=.FALSE.
3721	ENDIF
3722
3723	ENDDO
3724
3725	! then find closest point by computing distances
3726	il_shape(1)= il_isup - il_iinf + 1
3727	il_shape(2)= il_jsup - il_jinf + 1
3728
3729	ALLOCATE( dl_dist(il_shape(1), il_shape(2)) )
3730
3731	dl_dist(:,:)=grid_distance(dl_lon0(il_iinf:il_isup,il_jinf:il_jsup), &
3732	& dd_lat0(il_iinf:il_isup,il_jinf:il_jsup), &
3733	& dl_lon1, dd_lat1 )
3734
3735	IF( PRESENT(cd_pos) )THEN
3736	!
3737	SELECT CASE(TRIM(cd_pos))
3738	CASE('le')
3739	WHERE( dd_lat0(il_iinf:il_isup,il_jinf:il_jsup) > dd_lat1 )
3740	dl_dist(:,:)=NF90_FILL_DOUBLE
3741	END WHERE
3742	CASE('ri')
3743	WHERE( dd_lat0(il_iinf:il_isup,il_jinf:il_jsup) < dd_lat1 )
3744	dl_dist(:,:)=NF90_FILL_DOUBLE
3745	END WHERE
3746	CASE('up')
3747	WHERE( dl_lon0(il_iinf:il_isup,il_jinf:il_jsup) < dl_lon1 )
3748	dl_dist(:,:)=NF90_FILL_DOUBLE
3749	END WHERE
3750	CASE('lo')
3751	WHERE( dl_lon0(il_iinf:il_isup,il_jinf:il_jsup) > dl_lon1 )
3752	dl_dist(:,:)=NF90_FILL_DOUBLE
3753	END WHERE
3754	CASE('ll')
3755	WHERE( dl_lon0(il_iinf:il_isup,il_jinf:il_jsup) > dl_lon1 .OR. &
3756	& dd_lat0(il_iinf:il_isup,il_jinf:il_jsup) > dd_lat1 )
3757	dl_dist(:,:)=NF90_FILL_DOUBLE
3758	END WHERE
3759	CASE('lr')
3760	WHERE( dl_lon0(il_iinf:il_isup,il_jinf:il_jsup) < dl_lon1 .OR. &
3761	& dd_lat0(il_iinf:il_isup,il_jinf:il_jsup) > dd_lat1 )
3762	dl_dist(:,:)=NF90_FILL_DOUBLE
3763	END WHERE
3764	CASE('ul')
3765	WHERE( dl_lon0(il_iinf:il_isup,il_jinf:il_jsup) > dl_lon1 .OR. &
3766	& dd_lat0(il_iinf:il_isup,il_jinf:il_jsup) < dd_lat1 )
3767	dl_dist(:,:)=NF90_FILL_DOUBLE
3768	END WHERE
3769	CASE('ur')
3770	WHERE( dl_lon0(il_iinf:il_isup,il_jinf:il_jsup) < dl_lon1 .OR. &
3771	& dd_lat0(il_iinf:il_isup,il_jinf:il_jsup) < dd_lat1 )
3772	dl_dist(:,:)=NF90_FILL_DOUBLE
3773	END WHERE
3774	END SELECT
3775	ENDIF
3776	if_idx(:)=MINLOC(dl_dist(:,:),dl_dist(:,:)/=NF90_FILL_DOUBLE)
3777
3778	if_idx(1)=if_idx(1)+il_iinf-1
3779	if_idx(2)=if_idx(2)+il_jinf-1
3780
3781	DEALLOCATE( dl_dist )
3782	DEALLOCATE( dl_lon0 )
3783
3784	END FUNCTION grid__get_closest_arr
3785	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
3786	FUNCTION grid_distance(dd_lon, dd_lat, dd_lonA, dd_latA) &
3787	& RESULT (df_dist)
3788	!-------------------------------------------------------------------
3789	!> @brief This function compute the distance between a point A and grid points.
3790	!>
3791	!> @details
3792	!>
3793	!> @author J.Paul
3794	!> @date November, 2013 - Initial Version
3795	!>
3796	!> @param[in] dd_lon grid longitude array
3797	!> @param[in] dd_lat grid latitude array
3798	!> @param[in] dd_lonA longitude of point A
3799	!> @param[in] dd_latA latitude of point A
3800	!> @param[in] dd_fill
3801	!> @return array of distance between point A and grid points.
3802	!-------------------------------------------------------------------
3803
3804	IMPLICIT NONE
3805
3806	! Argument
3807	REAL(dp), DIMENSION(:,:), INTENT(IN) :: dd_lon
3808	REAL(dp), DIMENSION(:,:), INTENT(IN) :: dd_lat
3809	REAL(dp), INTENT(IN) :: dd_lonA
3810	REAL(dp), INTENT(IN) :: dd_latA
3811
3812	! function
3813	REAL(dp), DIMENSION(SIZE(dd_lon(:,:),DIM=1),&
3814	& SIZE(dd_lon(:,:),DIM=2)) :: df_dist
3815
3816	! local variable
3817	INTEGER(i4), DIMENSION(2) :: il_shape
3818
3819	REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_lon
3820	REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_lat
3821	REAL(dp) :: dl_lonA
3822	REAL(dp) :: dl_latA
3823
3824	REAL(dp) :: dl_tmp
3825
3826	! loop indices
3827	INTEGER(i4) :: ji
3828	INTEGER(i4) :: jj
3829	!----------------------------------------------------------------
3830
3831	IF( ANY( SHAPE(dd_lon(:,:)) /= SHAPE(dd_lat(:,:)) ) )THEN
3832	CALL logger_fatal("GRID DISTANCE: dimension of longitude and "//&
3833	& " latitude differ")
3834	ENDIF
3835	il_shape(:)=SHAPE(dd_lon(:,:))
3836
3837	ALLOCATE(dl_lon(il_shape(1),il_shape(2)))
3838	ALLOCATE(dl_lat(il_shape(1),il_shape(2)))
3839
3840	dl_lon(:,:) = dd_lon(:,:)
3841	dl_lonA = dd_lonA
3842
3843	WHERE(dd_lon(:,:) < 0 ) dl_lon(:,:) = dd_lon(:,:) + 360.
3844	IF( dd_lonA < 0 ) dl_lonA = dd_lonA + 360.
3845
3846	dl_lonA = dd_lonA * dp_deg2rad
3847	dl_latA = dd_latA * dp_deg2rad
3848
3849	dl_lon(:,:) = dl_lon(:,:) * dp_deg2rad
3850	dl_lat(:,:) = dd_lat(:,:) * dp_deg2rad
3851
3852	df_dist(:,:)=NF90_FILL_DOUBLE
3853
3854	DO jj=1,il_shape(2)
3855	DO ji=1,il_shape(1)
3856	IF( dl_lon(ji,jj) == dl_lonA .AND. &
3857	& dl_lat(ji,jj) == dl_latA )THEN
3858	df_dist(ji,jj)=0.0
3859	ELSE
3860	dl_tmp= SIN(dl_latA)*SIN(dl_lat(ji,jj)) + &
3861	& COS(dl_latA)COS(dl_lat(ji,jj)) &
3862	& COS(dl_lon(ji,jj)-dl_lonA)
3863	! check to avoid mistake with ACOS
3864	IF( dl_tmp < -1.0 ) dl_tmp = -1.0
3865	IF( dl_tmp > 1.0 ) dl_tmp = 1.0
3866	df_dist(ji,jj)=ACOS(dl_tmp)*dp_rearth
3867	ENDIF
3868	ENDDO
3869	ENDDO
3870
3871	DEALLOCATE(dl_lon)
3872	DEALLOCATE(dl_lat)
3873
3874	END FUNCTION grid_distance
3875	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
3876	FUNCTION grid__get_fine_offset_ff(td_coord0, &
3877	& id_imin0, id_jmin0, id_imax0, id_jmax0, &
3878	& td_coord1, id_rho, cd_point) &
3879	& RESULT (if_offset)
3880	!-------------------------------------------------------------------
3881	!> @brief This function get offset between fine grid and coarse grid.
3882	!>
3883	!> @details
3884	!> optionally, you could specify on which Arakawa grid point you want to
3885	!> work (default 'T')
3886	!> offset value could be 0,1,..,rho-1
3887	!>
3888	!> @author J.Paul
3889	!> @date September, 2014 - Initial Version
3890	!> @date October, 2014
3891	!> - work on mpp file structure instead of file structure
3892	!>
3893	!> @param[in] td_coord0 coarse grid coordinate
3894	!> @param[in] id_imin0 coarse grid lower left corner i-indice of fine grid domain
3895	!> @param[in] id_jmin0 coarse grid lower left corner j-indice of fine grid domain
3896	!> @param[in] id_imax0 coarse grid upper right corner i-indice of fine grid domain
3897	!> @param[in] id_jmax0 coarse grid upper right corner j-indice of fine grid domain
3898	!> @param[in] td_coord1 fine grid coordinate
3899	!> @param[in] id_rho array of refinement factor
3900	!> @param[in] cd_point Arakawa grid point
3901	!> @return offset array (/ (/i_offset_left,i_offset_right/),(/j_offset_lower,j_offset_upper/) /)
3902	!-------------------------------------------------------------------
3903
3904	IMPLICIT NONE
3905
3906	! Argument
3907	TYPE(TMPP) , INTENT(IN) :: td_coord0
3908	TYPE(TMPP) , INTENT(IN) :: td_coord1
3909
3910	INTEGER(i4) , INTENT(IN) :: id_imin0
3911	INTEGER(i4) , INTENT(IN) :: id_jmin0
3912	INTEGER(i4) , INTENT(IN) :: id_imax0
3913	INTEGER(i4) , INTENT(IN) :: id_jmax0
3914
3915	INTEGER(i4) , DIMENSION(:), INTENT(IN), OPTIONAL :: id_rho
3916	CHARACTER(LEN=*) , INTENT(IN), OPTIONAL :: cd_point
3917
3918	! function
3919	INTEGER(i4) , DIMENSION(2,2) :: if_offset
3920
3921	! local variable
3922	INTEGER(i4) :: il_imin0
3923	INTEGER(i4) :: il_jmin0
3924	INTEGER(i4) :: il_imax0
3925	INTEGER(i4) :: il_jmax0
3926	INTEGER(i4) :: il_ind
3927
3928	INTEGER(i4), DIMENSION(:), ALLOCATABLE :: il_rho
3929
3930	INTEGER(i4), DIMENSION(2,2) :: il_xghost0
3931	INTEGER(i4), DIMENSION(2,2) :: il_xghost1
3932
3933	CHARACTER(LEN= 1) :: cl_point
3934	CHARACTER(LEN=lc) :: cl_name
3935
3936	REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_lon0
3937	REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_lat0
3938	REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_lon1
3939	REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_lat1
3940
3941	TYPE(TVAR) :: tl_lon0
3942	TYPE(TVAR) :: tl_lat0
3943	TYPE(TVAR) :: tl_lon1
3944	TYPE(TVAR) :: tl_lat1
3945
3946	TYPE(TMPP) :: tl_coord0
3947	TYPE(TMPP) :: tl_coord1
3948
3949	! loop indices
3950	!----------------------------------------------------------------
3951	! init
3952	if_offset(:,:)=-1
3953
3954	ALLOCATE(il_rho(ip_maxdim))
3955	il_rho(:)=1
3956	IF( PRESENT(id_rho) ) il_rho(:)=id_rho(:)
3957
3958	cl_point='T'
3959	IF( PRESENT(cd_point) ) cl_point=TRIM(fct_upper(cd_point))
3960
3961	! copy structure
3962	tl_coord0=mpp_copy(td_coord0)
3963	tl_coord1=mpp_copy(td_coord1)
3964
3965	IF( .NOT. ASSOCIATED(tl_coord0%t_proc) .OR. &
3966	& .NOT. ASSOCIATED(tl_coord1%t_proc) )THEN
3967	CALL logger_error("GRID GET FINE OFFSET: can not get coarse "//&
3968	& "grid indices. decompsition of mpp file "//TRIM(tl_coord0%c_name)//&
3969	& " and/or "//TRIM(tl_coord1%c_name)//" not defined." )
3970	ELSE
3971	!1- Coarse grid
3972	! get ghost cell factor on coarse grid
3973	il_xghost0(:,:)=grid_get_ghost( tl_coord0 )
3974
3975	! open mpp files
3976	CALL iom_mpp_open(tl_coord0)
3977
3978	! read coarse longitue and latitude
3979	WRITE(cl_name,*) 'longitude_'//TRIM(cl_point)
3980	il_ind=var_get_id(tl_coord0%t_proc(1)%t_var(:), cl_name)
3981	IF( il_ind == 0 )THEN
3982	CALL logger_warn("GRID GET FINE OFFSET: no variable "//&
3983	& TRIM(cl_name)//" in file "//TRIM(tl_coord0%c_name)//". &
3984	& try to use longitude.")
3985	WRITE(cl_name,*) 'longitude'
3986	ENDIF
3987	tl_lon0=iom_mpp_read_var(tl_coord0, TRIM(cl_name))
3988
3989	WRITE(cl_name,*) 'latitude_'//TRIM(cl_point)
3990	il_ind=var_get_id(tl_coord0%t_proc(1)%t_var(:), cl_name)
3991	IF( il_ind == 0 )THEN
3992	CALL logger_warn("GRID GET FINE OFFSET: no variable "//&
3993	& TRIM(cl_name)//" in file "//TRIM(tl_coord0%c_name)//". &
3994	& try to use latitude.")
3995	WRITE(cl_name,*) 'latitude'
3996	ENDIF
3997	tl_lat0=iom_mpp_read_var(tl_coord0, TRIM(cl_name))
3998
3999	! close mpp files
4000	CALL iom_mpp_close(tl_coord0)
4001
4002	CALL grid_del_ghost(tl_lon0, il_xghost0(:,:))
4003	CALL grid_del_ghost(tl_lat0, il_xghost0(:,:))
4004
4005	ALLOCATE(dl_lon0(tl_lon0%t_dim(jp_I)%i_len, &
4006	& tl_lon0%t_dim(jp_J)%i_len ))
4007
4008	dl_lon0(:,:)=tl_lon0%d_value(:,:,1,1)
4009
4010	ALLOCATE(dl_lat0(tl_lat0%t_dim(jp_I)%i_len, &
4011	& tl_lat0%t_dim(jp_J)%i_len ))
4012
4013	dl_lat0(:,:)=tl_lat0%d_value(:,:,1,1)
4014
4015	! clean
4016	CALL var_clean(tl_lon0)
4017	CALL var_clean(tl_lat0)
4018
4019	! adjust coarse grid indices
4020	il_imin0=id_imin0-il_xghost0(jp_I,1)
4021	il_imax0=id_imax0-il_xghost0(jp_I,1)
4022
4023	il_jmin0=id_jmin0-il_xghost0(jp_J,1)
4024	il_jmax0=id_jmax0-il_xghost0(jp_J,1)
4025
4026	!2- Fine grid
4027	! get ghost cell factor on fine grid
4028	il_xghost1(:,:)=grid_get_ghost( tl_coord1 )
4029
4030	! open mpp files
4031	CALL iom_mpp_open(tl_coord1)
4032
4033	! read fine longitue and latitude
4034	WRITE(cl_name,*) 'longitude_'//TRIM(cl_point)
4035	il_ind=var_get_id(tl_coord1%t_proc(1)%t_var(:), cl_name)
4036	IF( il_ind == 0 )THEN
4037	CALL logger_warn("GRID GET FINE OFFSET: no variable "//&
4038	& TRIM(cl_name)//" in file "//TRIM(tl_coord1%c_name)//". &
4039	& try to use longitude.")
4040	WRITE(cl_name,*) 'longitude'
4041	ENDIF
4042	tl_lon1=iom_mpp_read_var(tl_coord1, TRIM(cl_name))
4043
4044	WRITE(cl_name,*) 'latitude_'//TRIM(cl_point)
4045	il_ind=var_get_id(tl_coord1%t_proc(1)%t_var(:), cl_name)
4046	IF( il_ind == 0 )THEN
4047	CALL logger_warn("GRID GET FINE OFFSET: no variable "//&
4048	& TRIM(cl_name)//" in file "//TRIM(tl_coord1%c_name)//". &
4049	& try to use latitude.")
4050	WRITE(cl_name,*) 'latitude'
4051	ENDIF
4052	tl_lat1=iom_mpp_read_var(tl_coord1, TRIM(cl_name))
4053
4054	! close mpp files
4055	CALL iom_mpp_close(tl_coord1)
4056
4057	CALL grid_del_ghost(tl_lon1, il_xghost1(:,:))
4058	CALL grid_del_ghost(tl_lat1, il_xghost1(:,:))
4059
4060	ALLOCATE(dl_lon1(tl_lon1%t_dim(jp_I)%i_len, &
4061	& tl_lon1%t_dim(jp_J)%i_len ))
4062
4063	dl_lon1(:,:)=tl_lon1%d_value(:,:,1,1)
4064
4065	ALLOCATE(dl_lat1(tl_lat1%t_dim(jp_I)%i_len, &
4066	& tl_lat1%t_dim(jp_J)%i_len ))
4067
4068	dl_lat1(:,:)=tl_lat1%d_value(:,:,1,1)
4069
4070	! clean
4071	CALL var_clean(tl_lon1)
4072	CALL var_clean(tl_lat1)
4073
4074	!3- compute
4075	if_offset(:,:)=grid_get_fine_offset( dl_lon0(:,:), dl_lat0(:,:),&
4076	& il_imin0, il_jmin0, &
4077	& il_imax0, il_jmax0, &
4078	& dl_lon1(:,:), dl_lat1(:,:),&
4079	& id_rho(:), cl_point )
4080
4081	DEALLOCATE(dl_lon0, dl_lat0)
4082	DEALLOCATE(dl_lon1, dl_lat1)
4083	ENDIF
4084
4085	! clean
4086	CALL mpp_clean(tl_coord0)
4087	CALL mpp_clean(tl_coord1)
4088	DEALLOCATE(il_rho)
4089
4090	END FUNCTION grid__get_fine_offset_ff
4091	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
4092	FUNCTION grid__get_fine_offset_cf(dd_lon0, dd_lat0, &
4093	& id_imin0, id_jmin0, id_imax0, id_jmax0, &
4094	& td_coord1, id_rho, cd_point) &
4095	& RESULT (if_offset)
4096	!-------------------------------------------------------------------
4097	!> @brief This function get offset between fine grid and coarse grid.
4098	!>
4099	!> @details
4100	!> optionally, you could specify on which Arakawa grid point you want to
4101	!> work (default 'T')
4102	!> offset value could be 0,1,..,rho-1
4103	!>
4104	!> @author J.Paul
4105	!> @date September, 2014 - Initial Version
4106	!> @date October, 2014
4107	!> - work on mpp file structure instead of file structure
4108	!>
4109	!> @param[in] dd_lon0 coarse grid longitude array
4110	!> @param[in] dd_lat0 coarse grid latitude array
4111	!> @param[in] id_imin0 coarse grid lower left corner i-indice of fine grid domain
4112	!> @param[in] id_jmin0 coarse grid lower left corner j-indice of fine grid domain
4113	!> @param[in] id_imax0 coarse grid upper right corner i-indice of fine grid domain
4114	!> @param[in] id_jmax0 coarse grid upper right corner j-indice of fine grid domain
4115	!> @param[in] td_coord1 fine grid coordinate
4116	!> @param[in] id_rho array of refinement factor
4117	!> @param[in] cd_point Arakawa grid point
4118	!> @return offset array (/ (/i_offset_left,i_offset_right/),(/j_offset_lower,j_offset_upper/) /)
4119	!-------------------------------------------------------------------
4120
4121	IMPLICIT NONE
4122
4123	! Argument
4124	REAL(dp) , DIMENSION(:,:), INTENT(IN) :: dd_lon0
4125	REAL(dp) , DIMENSION(:,:), INTENT(IN) :: dd_lat0
4126	TYPE(TMPP) , INTENT(IN) :: td_coord1
4127
4128	INTEGER(i4) , INTENT(IN) :: id_imin0
4129	INTEGER(i4) , INTENT(IN) :: id_jmin0
4130	INTEGER(i4) , INTENT(IN) :: id_imax0
4131	INTEGER(i4) , INTENT(IN) :: id_jmax0
4132
4133	INTEGER(i4) , DIMENSION(:) , INTENT(IN), OPTIONAL :: id_rho
4134	CHARACTER(LEN=*) , INTENT(IN), OPTIONAL :: cd_point
4135
4136	! function
4137	INTEGER(i4) , DIMENSION(2,2) :: if_offset
4138
4139	! local variable
4140	INTEGER(i4) :: il_ind
4141	INTEGER(i4), DIMENSION(2,2) :: il_xghost1
4142	INTEGER(i4), DIMENSION(:), ALLOCATABLE :: il_rho
4143
4144	CHARACTER(LEN= 1) :: cl_point
4145	CHARACTER(LEN=lc) :: cl_name
4146
4147	REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_lon1
4148	REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_lat1
4149
4150	TYPE(TVAR) :: tl_lon1
4151	TYPE(TVAR) :: tl_lat1
4152
4153	TYPE(TMPP) :: tl_coord1
4154	! loop indices
4155	!----------------------------------------------------------------
4156	! init
4157	if_offset(:,:)=-1
4158
4159	ALLOCATE(il_rho(ip_maxdim))
4160	il_rho(:)=1
4161	IF( PRESENT(id_rho) ) il_rho(:)=id_rho(:)
4162
4163	cl_point='T'
4164	IF( PRESENT(cd_point) ) cl_point=TRIM(fct_upper(cd_point))
4165
4166	! copy structure
4167	tl_coord1=mpp_copy(td_coord1)
4168
4169	IF( .NOT. ASSOCIATED(tl_coord1%t_proc) )THEN
4170	CALL logger_error("GRID GET FINE OFFSET: decompsition of mpp "//&
4171	& "file "//TRIM(tl_coord1%c_name)//" not defined." )
4172	ELSE
4173
4174	! Fine grid
4175	! get ghost cell factor on fine grid
4176	il_xghost1(:,:)=grid_get_ghost( tl_coord1 )
4177
4178	! open mpp files
4179	CALL iom_mpp_open(tl_coord1)
4180
4181	! read fine longitue and latitude
4182	WRITE(cl_name,*) 'longitude_'//TRIM(cl_point)
4183	il_ind=var_get_id(tl_coord1%t_proc(1)%t_var(:), cl_name)
4184	IF( il_ind == 0 )THEN
4185	CALL logger_warn("GRID GET FINE OFFSET: no variable "//&
4186	& TRIM(cl_name)//" in file "//TRIM(tl_coord1%c_name)//". &
4187	& try to use longitude.")
4188	WRITE(cl_name,*) 'longitude'
4189	ENDIF
4190	tl_lon1=iom_mpp_read_var(tl_coord1, TRIM(cl_name))
4191
4192	WRITE(cl_name,*) 'latitude_'//TRIM(cl_point)
4193	il_ind=var_get_id(tl_coord1%t_proc(1)%t_var(:), cl_name)
4194	IF( il_ind == 0 )THEN
4195	CALL logger_warn("GRID GET FINE OFFSET: no variable "//&
4196	& TRIM(cl_name)//" in file "//TRIM(tl_coord1%c_name)//". &
4197	& try to use latitude.")
4198	WRITE(cl_name,*) 'latitude'
4199	ENDIF
4200	tl_lat1=iom_mpp_read_var(tl_coord1, TRIM(cl_name))
4201
4202	! close mpp files
4203	CALL iom_mpp_close(tl_coord1)
4204
4205	CALL grid_del_ghost(tl_lon1, il_xghost1(:,:))
4206	CALL grid_del_ghost(tl_lat1, il_xghost1(:,:))
4207
4208	ALLOCATE(dl_lon1(tl_lon1%t_dim(jp_I)%i_len, &
4209	& tl_lon1%t_dim(jp_J)%i_len ))
4210
4211	dl_lon1(:,:)=tl_lon1%d_value(:,:,1,1)
4212
4213	ALLOCATE(dl_lat1(tl_lat1%t_dim(jp_I)%i_len, &
4214	& tl_lat1%t_dim(jp_J)%i_len ))
4215
4216	dl_lat1(:,:)=tl_lat1%d_value(:,:,1,1)
4217
4218	! clean
4219	CALL var_clean(tl_lon1)
4220	CALL var_clean(tl_lat1)
4221
4222	! compute
4223	if_offset(:,:)=grid_get_fine_offset( dd_lon0(:,:), dd_lat0(:,:),&
4224	& id_imin0, id_jmin0, &
4225	& id_imax0, id_jmax0, &
4226	& dl_lon1(:,:), dl_lat1(:,:),&
4227	& id_rho(:), cl_point )
4228
4229	DEALLOCATE(dl_lon1, dl_lat1)
4230	ENDIF
4231
4232	! clean
4233	CALL mpp_clean(tl_coord1)
4234	DEALLOCATE(il_rho)
4235
4236	END FUNCTION grid__get_fine_offset_cf
4237	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
4238	FUNCTION grid__get_fine_offset_fc(td_coord0, &
4239	& id_imin0, id_jmin0, id_imax0, id_jmax0, &
4240	& dd_lon1, dd_lat1, &
4241	& id_rho, cd_point) &
4242	& RESULT (if_offset)
4243	!-------------------------------------------------------------------
4244	!> @brief This function get offset between fine grid and coarse grid.
4245	!>
4246	!> @details
4247	!> optionally, you could specify on which Arakawa grid point you want to
4248	!> work (default 'T')
4249	!> offset value could be 0,1,..,rho-1
4250	!>
4251	!> @author J.Paul
4252	!> @date September, 2014 - Initial Version
4253	!> @date October, 2014
4254	!> - work on mpp file structure instead of file structure
4255	!>
4256	!> @param[in] td_coord0 coarse grid coordinate
4257	!> @param[in] id_imin0 coarse grid lower left corner i-indice of fine grid domain
4258	!> @param[in] id_jmin0 coarse grid lower left corner j-indice of fine grid domain
4259	!> @param[in] id_imax0 coarse grid upper right corner i-indice of fine grid domain
4260	!> @param[in] id_jmax0 coarse grid upper right corner j-indice of fine grid domain
4261	!> @param[in] dd_lon1 fine grid longitude array
4262	!> @param[in] dd_lat1 fine grid latitude array
4263	!> @param[in] id_rho array of refinement factor
4264	!> @param[in] cd_point Arakawa grid point
4265	!> @return offset array (/ (/i_offset_left,i_offset_right/),(/j_offset_lower,j_offset_upper/) /)
4266	!-------------------------------------------------------------------
4267
4268	IMPLICIT NONE
4269
4270	! Argument
4271	TYPE(TMPP) , INTENT(IN) :: td_coord0
4272	REAL(dp) , DIMENSION(:,:), INTENT(IN) :: dd_lon1
4273	REAL(dp) , DIMENSION(:,:), INTENT(IN) :: dd_lat1
4274
4275	INTEGER(i4) , INTENT(IN) :: id_imin0
4276	INTEGER(i4) , INTENT(IN) :: id_jmin0
4277	INTEGER(i4) , INTENT(IN) :: id_imax0
4278	INTEGER(i4) , INTENT(IN) :: id_jmax0
4279
4280	INTEGER(i4) , DIMENSION(:) , INTENT(IN), OPTIONAL :: id_rho
4281	CHARACTER(LEN=*) , INTENT(IN), OPTIONAL :: cd_point
4282
4283	! function
4284	INTEGER(i4) , DIMENSION(2,2) :: if_offset
4285
4286	! local variable
4287	INTEGER(i4) :: il_imin0
4288	INTEGER(i4) :: il_jmin0
4289	INTEGER(i4) :: il_imax0
4290	INTEGER(i4) :: il_jmax0
4291	INTEGER(i4) :: il_ind
4292
4293	INTEGER(i4), DIMENSION(:), ALLOCATABLE :: il_rho
4294
4295	INTEGER(i4), DIMENSION(2,2) :: il_xghost0
4296
4297	CHARACTER(LEN= 1) :: cl_point
4298	CHARACTER(LEN=lc) :: cl_name
4299
4300	REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_lon0
4301	REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_lat0
4302
4303	TYPE(TVAR) :: tl_lon0
4304	TYPE(TVAR) :: tl_lat0
4305
4306	TYPE(TMPP) :: tl_coord0
4307	! loop indices
4308	!----------------------------------------------------------------
4309	! init
4310	if_offset(:,:)=-1
4311	ALLOCATE(il_rho(ip_maxdim))
4312	il_rho(:)=1
4313	IF( PRESENT(id_rho) ) il_rho(:)=id_rho(:)
4314
4315	cl_point='T'
4316	IF( PRESENT(cd_point) ) cl_point=TRIM(fct_upper(cd_point))
4317
4318	! copy structure
4319	tl_coord0=mpp_copy(td_coord0)
4320
4321	IF( .NOT. ASSOCIATED(tl_coord0%t_proc) )THEN
4322	CALL logger_error("GRID GET FINE OFFSET: decompsition of mpp "//&
4323	& "file "//TRIM(tl_coord0%c_name)//" not defined." )
4324	ELSE
4325	!1- Coarse grid
4326	! get ghost cell factor on coarse grid
4327	il_xghost0(:,:)=grid_get_ghost( tl_coord0 )
4328
4329	! open mpp files
4330	CALL iom_mpp_open(tl_coord0)
4331
4332	! read coarse longitude and latitude
4333	WRITE(cl_name,*) 'longitude_'//TRIM(cl_point)
4334	il_ind=var_get_id(tl_coord0%t_proc(1)%t_var(:), cl_name)
4335	IF( il_ind == 0 )THEN
4336	CALL logger_warn("GRID GET FINE OFFSET: no variable "//&
4337	& TRIM(cl_name)//" in file "//TRIM(tl_coord0%c_name)//". &
4338	& try to use longitude.")
4339	WRITE(cl_name,*) 'longitude'
4340	ENDIF
4341	tl_lon0=iom_mpp_read_var(tl_coord0, TRIM(cl_name))
4342
4343	WRITE(cl_name,*) 'latitude_'//TRIM(cl_point)
4344	il_ind=var_get_id(tl_coord0%t_proc(1)%t_var(:), cl_name)
4345	IF( il_ind == 0 )THEN
4346	CALL logger_warn("GRID GET FINE OFFSET: no variable "//&
4347	& TRIM(cl_name)//" in file "//TRIM(tl_coord0%c_name)//". &
4348	& try to use latitude.")
4349	WRITE(cl_name,*) 'latitude'
4350	ENDIF
4351	tl_lat0=iom_mpp_read_var(tl_coord0, TRIM(cl_name))
4352
4353	! close mpp files
4354	CALL iom_mpp_close(tl_coord0)
4355
4356	CALL grid_del_ghost(tl_lon0, il_xghost0(:,:))
4357	CALL grid_del_ghost(tl_lat0, il_xghost0(:,:))
4358
4359
4360	ALLOCATE(dl_lon0(tl_lon0%t_dim(jp_I)%i_len, &
4361	& tl_lon0%t_dim(jp_J)%i_len ))
4362
4363	dl_lon0(:,:)=tl_lon0%d_value(:,:,1,1)
4364
4365	ALLOCATE(dl_lat0(tl_lat0%t_dim(jp_I)%i_len, &
4366	& tl_lat0%t_dim(jp_J)%i_len ))
4367
4368	dl_lat0(:,:)=tl_lat0%d_value(:,:,1,1)
4369
4370	! clean
4371	CALL var_clean(tl_lon0)
4372	CALL var_clean(tl_lat0)
4373
4374	! adjust coarse grid indices
4375	il_imin0=id_imin0-il_xghost0(jp_I,1)
4376	il_imax0=id_imax0-il_xghost0(jp_I,1)
4377
4378	il_jmin0=id_jmin0-il_xghost0(jp_J,1)
4379	il_jmax0=id_jmax0-il_xghost0(jp_J,1)
4380
4381	!3- compute
4382	if_offset(:,:)=grid_get_fine_offset( dl_lon0(:,:), dl_lat0(:,:),&
4383	& il_imin0, il_jmin0, &
4384	& il_imax0, il_jmax0, &
4385	& dd_lon1(:,:), dd_lat1(:,:),&
4386	& id_rho(:), cl_point )
4387
4388	DEALLOCATE(dl_lon0, dl_lat0)
4389	ENDIF
4390
4391	! clean
4392	CALL mpp_clean(tl_coord0)
4393	DEALLOCATE(il_rho)
4394
4395	END FUNCTION grid__get_fine_offset_fc
4396	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
4397	FUNCTION grid__get_fine_offset_cc(dd_lon0, dd_lat0, &
4398	& id_imin0, id_jmin0, id_imax0, id_jmax0, &
4399	& dd_lon1, dd_lat1, id_rho, cd_point) &
4400	& RESULT (if_offset)
4401	!-------------------------------------------------------------------
4402	!> @brief This function get offset between fine grid and coarse grid.
4403	!>
4404	!> @details
4405	!> offset value could be 0,1,..,rho-1
4406	!>
4407	!> @author J.Paul
4408	!> @date November, 2013 - Initial Version
4409	!> @date September, 2014
4410	!> - rename from grid_get_fine_offset
4411	!> @date May, 2015
4412	!> - improve way to find offset
4413	!> @date July, 2015
4414	!> - manage case close to greenwich meridian
4415	!> @date February, 2016
4416	!> - use grid_get_closest to assess offset
4417	!> - use delta (lon or lat)
4418	!> - manage cases for T,U,V or F point, with even or odd refinment
4419	!> - check lower left(upper right) fine grid point inside lower left(upper
4420	!> right) coarse grid cell.
4421	!>
4422	!> @todo check case close from North fold.
4423	!>
4424	!> @param[in] dd_lon0 coarse grid longitude array
4425	!> @param[in] dd_lat0 coarse grid latitude array
4426	!> @param[in] id_imin0 coarse grid lower left corner i-indice of fine grid domain
4427	!> @param[in] id_jmin0 coarse grid lower left corner j-indice of fine grid domain
4428	!> @param[in] id_imax0 coarse grid upper right corner i-indice of fine grid domain
4429	!> @param[in] id_jmax0 coarse grid upper right corner j-indice of fine grid domain
4430	!> @param[in] dd_lon1 fine grid longitude array
4431	!> @param[in] dd_lat1 fine grid latitude array
4432	!> @param[in] id_rho array of refinement factor
4433	!> @param[in] cd_point Arakawa grid point
4434	!> @return offset array (/ (/i_offset_left,i_offset_right/),(/j_offset_lower,j_offset_upper/) /)
4435	!-------------------------------------------------------------------
4436
4437	IMPLICIT NONE
4438
4439	! Argument
4440	REAL(dp) , DIMENSION(:,:), INTENT(IN) :: dd_lon0
4441	REAL(dp) , DIMENSION(:,:), INTENT(IN) :: dd_lat0
4442	REAL(dp) , DIMENSION(:,:), INTENT(IN) :: dd_lon1
4443	REAL(dp) , DIMENSION(:,:), INTENT(IN) :: dd_lat1
4444
4445	INTEGER(i4) , INTENT(IN) :: id_imin0
4446	INTEGER(i4) , INTENT(IN) :: id_jmin0
4447	INTEGER(i4) , INTENT(IN) :: id_imax0
4448	INTEGER(i4) , INTENT(IN) :: id_jmax0
4449
4450	INTEGER(i4) , DIMENSION(:) , INTENT(IN) :: id_rho
4451	CHARACTER(LEN=*) , INTENT(IN), OPTIONAL :: cd_point
4452
4453	! function
4454	INTEGER(i4) , DIMENSION(2,2) :: if_offset
4455
4456	! local variable
4457	CHARACTER(LEN= 1) :: cl_point
4458
4459	INTEGER(i4) :: i1
4460	INTEGER(i4) :: i2
4461	INTEGER(i4) :: j1
4462	INTEGER(i4) :: j2
4463
4464	INTEGER(i4), DIMENSION(2) :: il_shape0
4465	INTEGER(i4), DIMENSION(2) :: il_shape1
4466
4467	INTEGER(i4), DIMENSION(2) :: il_ind
4468
4469	REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_lon0
4470	REAL(dp) , DIMENSION(:,:), ALLOCATABLE :: dl_lon1
4471
4472	REAL(dp) :: dl_lonmax0
4473	REAL(dp) :: dl_latmax0
4474	REAL(dp) :: dl_lonmin0
4475	REAL(dp) :: dl_latmin0
4476
4477	REAL(dp) :: dl_lon0F
4478	REAL(dp) :: dl_lat0F
4479	REAL(dp) :: dl_dlon
4480	REAL(dp) :: dl_dlat
4481
4482	LOGICAL , DIMENSION(2) :: ll_even
4483	LOGICAL :: ll_greenwich
4484
4485	! loop indices
4486	INTEGER(i4) :: ii
4487	INTEGER(i4) :: ij
4488	!----------------------------------------------------------------
4489	IF( ANY( SHAPE(dd_lon0(:,:)) /= SHAPE(dd_lat0(:,:)) ) )THEN
4490	CALL logger_fatal("GRID GET FINE OFFSET: dimension of coarse "//&
4491	& "longitude and latitude differ")
4492	ENDIF
4493
4494	IF( ANY( SHAPE(dd_lon1(:,:)) /= SHAPE(dd_lat1(:,:)) ) )THEN
4495	CALL logger_fatal("GRID GET FINE OFFSET: dimension of fine "//&
4496	& "longitude and latitude differ")
4497	ENDIF
4498
4499	ll_even(:)=(/ (MOD(id_rho(jp_I),2)==0), (MOD(id_rho(jp_J),2)==0) /)
4500
4501	cl_point='T'
4502	IF( PRESENT(cd_point) ) cl_point=TRIM(fct_upper(cd_point))
4503
4504	il_shape0(:)=SHAPE(dd_lon0(:,:))
4505	ALLOCATE( dl_lon0(il_shape0(1),il_shape0(2)) )
4506
4507	il_shape1(:)=SHAPE(dd_lon1(:,:))
4508	ALLOCATE( dl_lon1(il_shape1(1),il_shape1(2)) )
4509
4510	dl_lon0(:,:)=dd_lon0(:,:)
4511	WHERE( dd_lon0(:,:) < 0 ) dl_lon0(:,:)=dd_lon0(:,:)+360.
4512
4513	dl_lon1(:,:)=dd_lon1(:,:)
4514	WHERE( dd_lon1(:,:) < 0 ) dl_lon1(:,:)=dd_lon1(:,:)+360.
4515
4516	! init
4517	if_offset(:,:)=-1
4518	ll_greenwich=.FALSE.
4519
4520	IF( il_shape1(jp_J) == 1 )THEN
4521
4522	if_offset(jp_J,:)=((id_rho(jp_J)-1)/2)
4523
4524	!!! work on i-direction
4525	!!! look for i-direction left offset
4526	i1=1 ; i2=MIN((id_rho(jp_I)+2),il_shape1(jp_I))
4527	j1=1 ; j2=1
4528
4529	! check if cross greenwich meridien
4530	IF( minval(dl_lon0(id_imin0:id_imin0+1,id_jmin0))<5. .OR. &
4531	& maxval(dl_lon0(id_imin0:id_imin0+1,id_jmin0))>355. )THEN
4532	! close to greenwich meridien
4533	ll_greenwich=.TRUE.
4534	! 0:360 => -180:180
4535	WHERE( dl_lon0(id_imin0:id_imin0+1,id_jmin0) > 180. )
4536	dl_lon0(id_imin0:id_imin0+1,id_jmin0) = &
4537	& dl_lon0(id_imin0:id_imin0+1,id_jmin0)-360.
4538	END WHERE
4539
4540	WHERE( dl_lon1(i1:i2,j1:j2) > 180. )
4541	dl_lon1(i1:i2,j1:j2)=dl_lon1(i1:i2,j1:j2)-360.
4542	END WHERE
4543	ENDIF
4544
4545	! max lognitude of the left cell
4546	dl_lonmax0=dl_lon0(id_imin0+1,id_jmin0)
4547	IF( dl_lon1(1,1) < dl_lonmax0 )THEN
4548
4549	!!!!! i-direction !!!!!
4550	IF( ll_even(jp_I) )THEN
4551	! even
4552	SELECT CASE(TRIM(cl_point))
4553	CASE('F','U')
4554	dl_dlon= ( dl_lon0(id_imin0+1,id_jmin0) - &
4555	& dl_lon0(id_imin0 ,id_jmin0) ) / &
4556	& ( 2.*id_rho(jp_I) )
4557	CASE DEFAULT
4558	dl_dlon=0
4559	END SELECT
4560	ELSE
4561	! odd
4562	dl_dlon= ( dl_lon0(id_imin0+1,id_jmin0) - &
4563	& dl_lon0(id_imin0 ,id_jmin0) ) / &
4564	& ( 2.*id_rho(jp_I) )
4565	ENDIF
4566
4567	dl_lon0F= dl_lon0(id_imin0+1,id_jmin0) + dl_dlon
4568	dl_lat0F= dd_lat0(id_imin0+1,id_jmin0)
4569
4570	il_ind(:)=grid_get_closest( dl_lon1(i1:i2,j1:j2), dd_lat1(i1:i2,j1:j2), &
4571	& dl_lon0F, dl_lat0F, 'le' )
4572
4573	ii=il_ind(1)
4574
4575	!!!!! i-direction !!!!!
4576	IF( ll_even(jp_I) )THEN
4577	! even
4578	SELECT CASE(TRIM(cl_point))
4579	CASE('T','V')
4580	if_offset(jp_I,1)=id_rho(jp_I)-ii
4581	CASE DEFAULT !'F','U'
4582	if_offset(jp_I,1)=(id_rho(jp_I)+1)-ii
4583	END SELECT
4584	ELSE
4585	! odd
4586	if_offset(jp_I,1)=(id_rho(jp_I)+1)-ii
4587	ENDIF
4588
4589	ELSE
4590	CALL logger_error("GRID GET FINE OFFSET: coarse grid indices do "//&
4591	& " not match fine grid left corner.")
4592	ENDIF
4593
4594	IF( ll_greenwich )THEN
4595	! close to greenwich meridien
4596	ll_greenwich=.FALSE.
4597	! -180:180 => 0:360
4598	WHERE( dl_lon0(id_imin0:id_imin0+1,id_jmin0) < 0. )
4599	dl_lon0(id_imin0:id_imin0+1,id_jmin0) = &
4600	& dl_lon0(id_imin0:id_imin0+1,id_jmin0)+360.
4601	END WHERE
4602
4603	WHERE( dl_lon1(i1:i2,j1:j2) < 0. )
4604	dl_lon1(i1:i2,j1:j2)=dl_lon1(i1:i2,j1:j2)+360.
4605	END WHERE
4606	ENDIF
4607
4608	!!!!!! look for i-direction right offset !!!!!!
4609	i1=MAX(1,il_shape1(jp_I)-(id_rho(jp_I)+2)+1) ; i2=il_shape1(jp_I)
4610	j1=1 ; j2=1
4611
4612	! check if cross greenwich meridien
4613	IF( minval(dl_lon0(id_imax0-1:id_imax0,id_jmin0))<5. .OR. &
4614	& maxval(dl_lon0(id_imax0-1:id_imax0,id_jmin0))>355. )THEN
4615	! close to greenwich meridien
4616	ll_greenwich=.TRUE.
4617	! 0:360 => -180:180
4618	WHERE( dl_lon0(id_imax0-1:id_imax0,id_jmin0) > 180. )
4619	dl_lon0(id_imax0-1:id_imax0,id_jmin0) = &
4620	& dl_lon0(id_imax0-1:id_imax0,id_jmin0)-360.
4621	END WHERE
4622
4623	WHERE( dl_lon1(i1:i2,j1:j2) > 180. )
4624	dl_lon1(i1:i2,j1:j2)=dl_lon1(i1:i2,j1:j2)-360.
4625	END WHERE
4626	ENDIF
4627
4628	! min lognitude of the right cell
4629	dl_lonmin0=dl_lon0(id_imax0-1,id_jmin0)
4630	IF( dl_lon1(il_shape1(jp_I),il_shape1(jp_J)) > dl_lonmin0 )THEN
4631
4632	!!!!! i-direction !!!!!
4633	IF( ll_even(jp_I) )THEN
4634	! even
4635	SELECT CASE(TRIM(cl_point))
4636	CASE('F','U')
4637	dl_dlon= ( dl_lon0(id_imax0 ,id_jmin0) - &
4638	& dl_lon0(id_imax0-1,id_jmin0) ) / &
4639	& ( 2.*id_rho(jp_I) )
4640	CASE DEFAULT
4641	dl_dlon=0
4642	END SELECT
4643	ELSE
4644	! odd
4645	dl_dlon= ( dl_lon0(id_imax0 ,id_jmin0) - &
4646	& dl_lon0(id_imax0-1,id_jmin0) ) / &
4647	& ( 2.*id_rho(jp_I) )
4648	ENDIF
4649
4650	dl_lon0F= dl_lon0(id_imax0-1,id_jmin0) - dl_dlon
4651	dl_lat0F= dd_lat0(id_imax0-1,id_jmin0)
4652
4653	il_ind(:)=grid_get_closest( dl_lon1(i1:i2,j1:j2), dd_lat1(i1:i2,j1:j2), &
4654	& dl_lon0F, dl_lat0F, 'ri' )
4655
4656	ii=(MIN(il_shape1(jp_I),(id_rho(jp_I)+2))-il_ind(1)+1)
4657
4658	!!!!! i-direction !!!!!
4659	IF( ll_even(jp_I) )THEN
4660	! even
4661	SELECT CASE(TRIM(cl_point))
4662	CASE('T','V')
4663	if_offset(jp_I,2)=id_rho(jp_I)-ii
4664	CASE DEFAULT !'F','U'
4665	if_offset(jp_I,2)=(id_rho(jp_I)+1)-ii
4666	END SELECT
4667	ELSE
4668	! odd
4669	if_offset(jp_I,2)=(id_rho(jp_I)+1)-ii
4670	ENDIF
4671
4672	ELSE
4673	CALL logger_error("GRID GET FINE OFFSET: coarse grid indices do "//&
4674	& " not match fine grid right corner.")
4675	ENDIF
4676
4677	IF( ll_greenwich )THEN
4678	! close to greenwich meridien
4679	ll_greenwich=.FALSE.
4680	! -180:180 => 0:360
4681	WHERE( dl_lon0(id_imax0-1:id_imax0,id_jmin0) < 0. )
4682	dl_lon0(id_imax0-1:id_imax0,id_jmin0) = &
4683	& dl_lon0(id_imax0-1:id_imax0,id_jmin0)+360.
4684	END WHERE
4685
4686	WHERE( dl_lon1(i1:i2,j1:j2) < 0. )
4687	dl_lon1(i1:i2,j1:j2)=dl_lon1(i1:i2,j1:j2)+360.
4688	END WHERE
4689	ENDIF
4690
4691	ELSEIF( il_shape1(jp_I) == 1 )THEN
4692
4693	if_offset(jp_I,:)=((id_rho(jp_I)-1)/2)
4694
4695	!!! work on j-direction
4696	!!! look for j-direction lower offset
4697	i1=1 ; i2=1
4698	j1=1 ; j2=MIN((id_rho(jp_J)+2),il_shape1(jp_J))
4699
4700
4701	! max latitude of the lower cell
4702	dl_latmax0=dd_lat0(id_imin0,id_jmin0+1)
4703	IF( dd_lat1(1,1) < dl_latmax0 )THEN
4704
4705	IF( ll_even(jp_J) )THEN
4706	! even
4707	SELECT CASE(TRIM(cl_point))
4708	CASE('F','V')
4709	dl_dlat= ( dd_lat0(id_imin0,id_jmin0+1) - &
4710	& dd_lat0(id_imin0,id_jmin0 ) ) / &
4711	& ( 2.*id_rho(jp_J) )
4712	CASE DEFAULT
4713	dl_dlat=0
4714	END SELECT
4715	ELSE
4716	! odd
4717	dl_dlat= ( dd_lat0(id_imin0,id_jmin0+1) - &
4718	& dd_lat0(id_imin0,id_jmin0 ) ) / &
4719	& ( 2.*id_rho(jp_J) )
4720	ENDIF
4721
4722	dl_lon0F= dl_lon0(id_imin0,id_jmin0+1)
4723	dl_lat0F= dd_lat0(id_imin0,id_jmin0+1) + dl_dlat
4724
4725	il_ind(:)=grid_get_closest( dl_lon1(i1:i2,j1:j2), dd_lat1(i1:i2,j1:j2), &
4726	& dl_lon0F, dl_lat0F, 'lo' )
4727
4728	ij=il_ind(2)
4729
4730	!!!!! j-direction !!!!!
4731	IF( ll_even(jp_J) )THEN
4732	! even
4733	SELECT CASE(TRIM(cl_point))
4734	CASE('T','V')
4735	if_offset(jp_J,1)=id_rho(jp_J)-ij
4736	CASE DEFAULT !'F','U'
4737	if_offset(jp_J,1)=(id_rho(jp_J)+1)-ij
4738	END SELECT
4739	ELSE
4740	! odd
4741	if_offset(jp_J,1)=(id_rho(jp_J)+1)-ij
4742	ENDIF
4743
4744	ELSE
4745	CALL logger_error("GRID GET FINE OFFSET: coarse grid indices do "//&
4746	& " not match fine grid lower corner.")
4747	ENDIF
4748
4749	!!! look for j-direction upper offset
4750	i1=1 ; i2=1
4751	j1=MAX(1,il_shape1(jp_J)-(id_rho(jp_J)+2)+1) ; j2=il_shape1(jp_J)
4752
4753	! min latitude of the upper cell
4754	dl_latmin0=dd_lat0(id_imin0,id_jmax0-1)
4755	IF( dd_lat1(il_shape1(jp_I),il_shape1(jp_J)) > dl_latmin0 )THEN
4756
4757	IF( ll_even(jp_J) )THEN
4758	! even
4759	SELECT CASE(TRIM(cl_point))
4760	CASE('F','V')
4761	dl_dlat= ( dd_lat0(id_imin0,id_jmax0 ) - &
4762	& dd_lat0(id_imin0,id_jmax0-1) ) / &
4763	& ( 2.*id_rho(jp_J) )
4764	CASE DEFAULT
4765	dl_dlat=0
4766	END SELECT
4767	ELSE
4768	! odd
4769	dl_dlat= ( dd_lat0(id_imin0,id_jmax0 ) - &
4770	& dd_lat0(id_imin0,id_jmax0-1) ) / &
4771	& ( 2*id_rho(jp_J) )
4772	ENDIF
4773
4774	dl_lon0F= dl_lon0(id_imin0,id_jmax0-1)
4775	dl_lat0F= dd_lat0(id_imin0,id_jmax0-1) - dl_dlat
4776
4777	il_ind(:)=grid_get_closest( dl_lon1(i1:i2,j1:j2), dd_lat1(i1:i2,j1:j2), &
4778	& dl_lon0F, dl_lat0F, 'up' )
4779
4780	ij=(MIN(il_shape1(jp_J),(id_rho(jp_J)+2))-il_ind(2)+1)
4781
4782	!!!!! j-direction !!!!!
4783	IF( ll_even(jp_J) )THEN
4784	! even
4785	SELECT CASE(TRIM(cl_point))
4786	CASE('T','U')
4787	if_offset(jp_J,2)=id_rho(jp_J)-ij
4788	CASE DEFAULT !'F','V'
4789	if_offset(jp_J,2)=(id_rho(jp_J)+1)-ij
4790	END SELECT
4791	ELSE
4792	! odd
4793	if_offset(jp_J,2)=(id_rho(jp_J)+1)-ij
4794	ENDIF
4795
4796	ELSE
4797	CALL logger_error("GRID GET FINE OFFSET: coarse grid indices do "//&
4798	& " not match fine grid upper corner.")
4799	ENDIF
4800
4801	ELSE ! il_shape1(1) > 1 .AND. il_shape1(2) > 1
4802
4803	!!!!!! look for lower left offset !!!!!!
4804	i1=1 ; i2=MIN((id_rho(jp_I)+2),il_shape1(jp_I))
4805	j1=1 ; j2=MIN((id_rho(jp_J)+2),il_shape1(jp_J))
4806
4807	! check if cross greenwich meridien
4808	IF( minval(dl_lon0(id_imin0:id_imin0+1,id_jmin0:id_jmin0+1))<5. .OR. &
4809	& maxval(dl_lon0(id_imin0:id_imin0+1,id_jmin0:id_jmin0+1))>355. )THEN
4810	! close to greenwich meridien
4811	ll_greenwich=.TRUE.
4812	! 0:360 => -180:180
4813	WHERE( dl_lon0(id_imin0:id_imin0+1,id_jmin0:id_jmin0+1) > 180. )
4814	dl_lon0(id_imin0:id_imin0+1,id_jmin0:id_jmin0+1) = &
4815	& dl_lon0(id_imin0:id_imin0+1,id_jmin0:id_jmin0+1)-360.
4816	END WHERE
4817
4818	WHERE( dl_lon1(i1:i2,j1:j2) > 180. )
4819	dl_lon1(i1:i2,j1:j2)=dl_lon1(i1:i2,j1:j2)-360.
4820	END WHERE
4821	ENDIF
4822
4823	! max longitude of the lower left cell
4824	dl_lonmax0=MAX(dl_lon0(id_imin0+1,id_jmin0),dl_lon0(id_imin0+1,id_jmin0+1))
4825	! max latitude of the lower left cell
4826	dl_latmax0=MAX(dd_lat0(id_imin0,id_jmin0+1),dd_lat0(id_imin0+1,id_jmin0+1))
4827	IF( dl_lon1(1,1) < dl_lonmax0 .AND. &
4828	& dd_lat1(1,1) < dl_latmax0 )THEN
4829
4830	!!!!! i-direction !!!!!
4831	IF( ll_even(jp_I) )THEN
4832	! even
4833	SELECT CASE(TRIM(cl_point))
4834	CASE('F','U')
4835	dl_dlon= ( dl_lon0(id_imin0+1,id_jmin0+1) - &
4836	& dl_lon0(id_imin0 ,id_jmin0+1) ) / &
4837	& ( 2.*id_rho(jp_I) )
4838	CASE DEFAULT
4839	dl_dlon=0
4840	END SELECT
4841	ELSE
4842	! odd
4843	dl_dlon= ( dl_lon0(id_imin0+1,id_jmin0+1) - &
4844	& dl_lon0(id_imin0 ,id_jmin0+1) ) / &
4845	& ( 2.*id_rho(jp_I) )
4846	ENDIF
4847
4848	!!!!! j-direction !!!!!
4849	IF( ll_even(jp_J) )THEN
4850	! even
4851	SELECT CASE(TRIM(cl_point))
4852	CASE('F','V')
4853	dl_dlat= ( dd_lat0(id_imin0+1,id_jmin0+1) - &
4854	& dd_lat0(id_imin0+1,id_jmin0 ) ) / &
4855	& ( 2.*id_rho(jp_J) )
4856	CASE DEFAULT
4857	dl_dlat=0
4858	END SELECT
4859	ELSE
4860	! odd
4861	dl_dlat= ( dd_lat0(id_imin0+1,id_jmin0+1) - &
4862	& dd_lat0(id_imin0+1,id_jmin0 ) ) / &
4863	& ( 2.*id_rho(jp_J) )
4864	ENDIF
4865
4866	dl_lon0F= dl_lon0(id_imin0+1,id_jmin0+1) + dl_dlon
4867	dl_lat0F= dd_lat0(id_imin0+1,id_jmin0+1) + dl_dlat
4868
4869	il_ind(:)=grid_get_closest( dl_lon1(i1:i2,j1:j2), dd_lat1(i1:i2,j1:j2), &
4870	& dl_lon0F, dl_lat0F, 'll' )
4871
4872	ii=il_ind(1)
4873	ij=il_ind(2)
4874
4875	!!!!! i-direction !!!!!
4876	IF( ll_even(jp_I) )THEN
4877	! even
4878	SELECT CASE(TRIM(cl_point))
4879	CASE('T','V')
4880	if_offset(jp_I,1)=id_rho(jp_I)-ii
4881	CASE DEFAULT !'F','U'
4882	if_offset(jp_I,1)=(id_rho(jp_I)+1)-ii
4883	END SELECT
4884	ELSE
4885	! odd
4886	if_offset(jp_I,1)=(id_rho(jp_I)+1)-ii
4887	ENDIF
4888
4889	!!!!! j-direction !!!!!
4890	IF( ll_even(jp_J) )THEN
4891	! even
4892	SELECT CASE(TRIM(cl_point))
4893	CASE('T','U')
4894	if_offset(jp_J,1)=id_rho(jp_J)-ij
4895	CASE DEFAULT !'F','V'
4896	if_offset(jp_J,1)=(id_rho(jp_J)+1)-ij
4897	END SELECT
4898	ELSE
4899	! odd
4900	if_offset(jp_J,1)=(id_rho(jp_J)+1)-ij
4901	ENDIF
4902
4903	ELSE
4904	CALL logger_error("GRID GET FINE OFFSET: coarse grid indices do"//&
4905	& " not match fine grid lower left corner.")
4906	ENDIF
4907
4908	IF( ll_greenwich )THEN
4909	! close to greenwich meridien
4910	ll_greenwich=.FALSE.
4911	! -180:180 => 0:360
4912	WHERE( dl_lon0(id_imin0:id_imin0+1,id_jmin0:id_jmin0+1) < 0. )
4913	dl_lon0(id_imin0:id_imin0+1,id_jmin0:id_jmin0+1) = &
4914	& dl_lon0(id_imin0:id_imin0+1,id_jmin0:id_jmin0+1)+360.
4915	END WHERE
4916
4917	WHERE( dl_lon1(i1:i2,j1:j2) < 0. )
4918	dl_lon1(i1:i2,j1:j2)=dl_lon1(i1:i2,j1:j2)+360.
4919	END WHERE
4920	ENDIF
4921
4922	!!!!!! look for upper right offset !!!!!!
4923	i1=MAX(1,il_shape1(jp_I)-(id_rho(jp_I)+2)+1) ; i2=il_shape1(jp_I)
4924	j1=MAX(1,il_shape1(jp_J)-(id_rho(jp_J)+2)+1) ; j2=il_shape1(jp_J)
4925
4926	! check if cross greenwich meridien
4927	IF( minval(dl_lon0(id_imax0-1:id_imax0,id_jmax0-1:id_jmax0))<5. .OR. &
4928	& maxval(dl_lon0(id_imax0-1:id_imax0,id_jmax0-1:id_jmax0))>355. )THEN
4929	! close to greenwich meridien
4930	ll_greenwich=.TRUE.
4931	! 0:360 => -180:180
4932	WHERE( dl_lon0(id_imax0-1:id_imax0,id_jmax0-1:id_jmax0) > 180. )
4933	dl_lon0(id_imax0-1:id_imax0,id_jmax0-1:id_jmax0) = &
4934	& dl_lon0(id_imax0-1:id_imax0,id_jmax0-1:id_jmax0)-360.
4935	END WHERE
4936
4937	WHERE( dl_lon1(i1:i2,j1:j2) > 180. )
4938	dl_lon1(i1:i2,j1:j2)=dl_lon1(i1:i2,j1:j2)-360.
4939	END WHERE
4940	ENDIF
4941
4942	! min latitude of the upper right cell
4943	dl_lonmin0=MIN(dl_lon0(id_imax0-1,id_jmax0-1),dl_lon0(id_imax0-1,id_jmax0))
4944	! min latitude of the upper right cell
4945	dl_latmin0=MIN(dd_lat0(id_imax0-1,id_jmax0-1),dd_lat0(id_imax0,id_jmax0-1))
4946	IF( dl_lon1(il_shape1(jp_I),il_shape1(jp_J)) > dl_lonmin0 .AND. &
4947	& dd_lat1(il_shape1(jp_I),il_shape1(jp_J)) > dl_latmin0 )THEN
4948
4949	!!!!! i-direction !!!!!
4950	IF( ll_even(jp_I) )THEN
4951	! even
4952	SELECT CASE(TRIM(cl_point))
4953	CASE('F','U')
4954	dl_dlon= ( dl_lon0(id_imax0 ,id_jmax0-1) - &
4955	& dl_lon0(id_imax0-1,id_jmax0-1) ) / &
4956	& ( 2.*id_rho(jp_I) )
4957	CASE DEFAULT
4958	dl_dlon=0
4959	END SELECT
4960	ELSE
4961	! odd
4962	dl_dlon= ( dl_lon0(id_imax0 ,id_jmax0-1) - &
4963	& dl_lon0(id_imax0-1,id_jmax0-1) ) / &
4964	& ( 2*id_rho(jp_I) )
4965	ENDIF
4966
4967	!!!!! j-direction !!!!!
4968	IF( ll_even(jp_J) )THEN
4969	! even
4970	SELECT CASE(TRIM(cl_point))
4971	CASE('F','V')
4972	dl_dlat= ( dd_lat0(id_imax0-1,id_jmax0 ) - &
4973	& dd_lat0(id_imax0-1,id_jmax0-1) ) / &
4974	& ( 2.*id_rho(jp_J) )
4975	CASE DEFAULT
4976	dl_dlat=0
4977	END SELECT
4978	ELSE
4979	! odd
4980	dl_dlat= ( dd_lat0(id_imax0-1,id_jmax0 ) - &
4981	& dd_lat0(id_imax0-1,id_jmax0-1) ) / &
4982	& ( 2*id_rho(jp_J) )
4983	ENDIF
4984
4985	dl_lon0F= dl_lon0(id_imax0-1,id_jmax0-1) - dl_dlon
4986	dl_lat0F= dd_lat0(id_imax0-1,id_jmax0-1) - dl_dlat
4987
4988	il_ind(:)=grid_get_closest( dl_lon1(i1:i2,j1:j2), dd_lat1(i1:i2,j1:j2), &
4989	& dl_lon0F, dl_lat0F, 'ur' )
4990
4991	ii=(MIN(il_shape1(jp_I),(id_rho(jp_I)+2))-il_ind(1)+1)
4992	ij=(MIN(il_shape1(jp_J),(id_rho(jp_J)+2))-il_ind(2)+1)
4993
4994	!!!!! i-direction !!!!!
4995	IF( ll_even(jp_I) )THEN
4996	! even
4997	SELECT CASE(TRIM(cl_point))
4998	CASE('T','V')
4999	if_offset(jp_I,2)=id_rho(jp_I)-ii
5000	CASE DEFAULT !'F','U'
5001	if_offset(jp_I,2)=(id_rho(jp_I)+1)-ii
5002	END SELECT
5003	ELSE
5004	! odd
5005	if_offset(jp_I,2)=(id_rho(jp_I)+1)-ii
5006	ENDIF
5007
5008	!!!!! j-direction !!!!!
5009	IF( ll_even(jp_J) )THEN
5010	! even
5011	SELECT CASE(TRIM(cl_point))
5012	CASE('T','U')
5013	if_offset(jp_J,2)=id_rho(jp_J)-ij
5014	CASE DEFAULT !'F','V'
5015	if_offset(jp_J,2)=(id_rho(jp_J)+1)-ij
5016	END SELECT
5017	ELSE
5018	! odd
5019	if_offset(jp_J,2)=(id_rho(jp_J)+1)-ij
5020	ENDIF
5021
5022	ELSE
5023	CALL logger_error("GRID GET FINE OFFSET: coarse grid indices do"//&
5024	& " not match fine grid upper right corner.")
5025	ENDIF
5026
5027	IF( ll_greenwich )THEN
5028	! close to greenwich meridien
5029	ll_greenwich=.FALSE.
5030	! -180:180 => 0:360
5031	WHERE( dl_lon0(id_imax0-1:id_imax0,id_jmax0-1:id_jmax0) < 0. )
5032	dl_lon0(id_imax0-1:id_imax0,id_jmax0-1:id_jmax0) = &
5033	& dl_lon0(id_imax0-1:id_imax0,id_jmax0-1:id_jmax0)+360.
5034	END WHERE
5035
5036	WHERE( dl_lon1(i1:i2,j1:j2) < 0. )
5037	dl_lon1(i1:i2,j1:j2)=dl_lon1(i1:i2,j1:j2)+360.
5038	END WHERE
5039	ENDIF
5040
5041	ENDIF
5042
5043	DEALLOCATE( dl_lon0 )
5044	DEALLOCATE( dl_lon1 )
5045
5046	IF( ANY(if_offset(:,:)==-1) )THEN
5047	CALL logger_fatal("GRID GET FINE OFFSET: can not found "//&
5048	& " offset between coarse and fine grid.")
5049	ENDIF
5050
5051	END FUNCTION grid__get_fine_offset_cc
5052	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
5053	SUBROUTINE grid_check_coincidence(td_coord0, td_coord1, &
5054	& id_imin0, id_imax0, &
5055	& id_jmin0, id_jmax0, &
5056	& id_rho)
5057	!-------------------------------------------------------------------
5058	!> @brief This subroutine check fine and coarse grid coincidence.
5059	!>
5060	!> @details
5061	!>
5062	!> @author J.Paul
5063	!> @date November, 2013- Initial Version
5064	!> @date October, 2014
5065	!> - work on mpp file structure instead of file structure
5066	!> @date February, 2016
5067	!> - use F-point to check coincidence for even refinment
5068	!> - use F-point estimation, if can not read it.
5069	!>
5070	!> @param[in] td_coord0 coarse grid coordinate file structure
5071	!> @param[in] td_coord1 fine grid coordinate file structure
5072	!> @param[in] id_imin0 coarse grid lower left corner i-indice of fine grid domain
5073	!> @param[in] id_imax0 coarse grid upper right corner i-indice of fine grid domain
5074	!> @param[in] id_jmin0 coarse grid lower left corner j-indice of fine grid domain
5075	!> @param[in] id_jmax0 coarse grid upper right corner j-indice of fine grid domain
5076	!> @param[in] id_rho array of refinement factor
5077	!-------------------------------------------------------------------
5078
5079	IMPLICIT NONE
5080
5081	! Argument
5082	TYPE(TMPP) , INTENT(IN) :: td_coord0
5083	TYPE(TMPP) , INTENT(IN) :: td_coord1
5084	INTEGER(i4) , INTENT(IN) :: id_imin0
5085	INTEGER(i4) , INTENT(IN) :: id_imax0
5086	INTEGER(i4) , INTENT(IN) :: id_jmin0
5087	INTEGER(i4) , INTENT(IN) :: id_jmax0
5088	INTEGER(i4), DIMENSION(:), INTENT(IN) :: id_rho
5089
5090	! local variable
5091	INTEGER(i4) :: il_imid1
5092	INTEGER(i4) :: il_jmid1
5093
5094	INTEGER(i4) :: il_ew0
5095	INTEGER(i4) :: il_ew1
5096
5097	INTEGER(i4) :: il_ind
5098
5099	INTEGER(i4) :: il_imin1
5100	INTEGER(i4) :: il_imax1
5101	INTEGER(i4) :: il_jmin1
5102	INTEGER(i4) :: il_jmax1
5103
5104	INTEGER(i4), DIMENSION(2) :: il_ind0
5105	INTEGER(i4), DIMENSION(2) :: il_ind1
5106
5107	INTEGER(i4), DIMENSION(2) :: il_ill1
5108	INTEGER(i4), DIMENSION(2) :: il_ilr1
5109	INTEGER(i4), DIMENSION(2) :: il_iul1
5110	INTEGER(i4), DIMENSION(2) :: il_iur1
5111
5112	REAL(dp) :: dl_lon0F
5113	REAL(dp) :: dl_lat0F
5114	REAL(dp) :: dl_lon0
5115	REAL(dp) :: dl_lat0
5116	REAL(dp) :: dl_lon1F
5117	REAL(dp) :: dl_lat1F
5118	REAL(dp) :: dl_lon1
5119	REAL(dp) :: dl_lat1
5120
5121	REAL(dp) :: dl_delta
5122
5123	LOGICAL :: ll_coincidence
5124	LOGICAL :: ll_even
5125	LOGICAL :: ll_grid0F
5126	LOGICAL :: ll_grid1F
5127
5128	TYPE(TVAR) :: tl_lon0
5129	TYPE(TVAR) :: tl_lat0
5130	TYPE(TVAR) :: tl_lon0F
5131	TYPE(TVAR) :: tl_lat0F
5132	TYPE(TVAR) :: tl_lon1
5133	TYPE(TVAR) :: tl_lat1
5134	TYPE(TVAR) :: tl_lon1F
5135	TYPE(TVAR) :: tl_lat1F
5136
5137	TYPE(TMPP) :: tl_coord0
5138	TYPE(TMPP) :: tl_coord1
5139
5140	TYPE(TDOM) :: tl_dom0
5141
5142	! loop indices
5143	INTEGER(i4) :: ji
5144	INTEGER(i4) :: jj
5145	!----------------------------------------------------------------
5146
5147	ll_coincidence=.TRUE.
5148
5149	ll_even=.FALSE.
5150	IF( MOD(id_rho(jp_I)*id_rho(jp_J),2) == 0 )THEN
5151	ll_even=.TRUE.
5152	ENDIF
5153
5154	! copy structure
5155	tl_coord0=mpp_copy(td_coord0)
5156
5157	! compute domain
5158	tl_dom0=dom_init( tl_coord0, &
5159	& id_imin0, id_imax0,&
5160	& id_jmin0, id_jmax0 )
5161
5162	! open mpp files
5163	CALL iom_dom_open(tl_coord0, tl_dom0)
5164
5165	! read variable value on domain
5166	il_ind=var_get_index(tl_coord0%t_proc(1)%t_var(:), 'longitude_T')
5167	IF( il_ind /= 0 )THEN
5168	tl_lon0=iom_dom_read_var(tl_coord0,'longitude_T',tl_dom0)
5169	ELSE
5170	tl_lon0=iom_dom_read_var(tl_coord0,'longitude',tl_dom0)
5171	ENDIF
5172
5173	il_ind=var_get_index(tl_coord0%t_proc(1)%t_var(:), 'latitude_T')
5174	IF( il_ind /= 0 )THEN
5175	tl_lat0=iom_dom_read_var(tl_coord0,'latitude_T' ,tl_dom0)
5176	ELSE
5177	tl_lat0=iom_dom_read_var(tl_coord0,'latitude' ,tl_dom0)
5178	ENDIF
5179
5180	IF( ll_even )THEN
5181
5182	! look for variable value on domain for F point
5183	il_ind=var_get_index(tl_coord0%t_proc(1)%t_var(:), 'longitude_F')
5184	IF( il_ind /= 0 )THEN
5185	tl_lon0F=iom_dom_read_var(tl_coord0,'longitude_F',tl_dom0)
5186	ENDIF
5187
5188	il_ind=var_get_index(tl_coord0%t_proc(1)%t_var(:), 'latitude_F')
5189	IF( il_ind /= 0 )THEN
5190	tl_lat0F=iom_dom_read_var(tl_coord0,'latitude_F' ,tl_dom0)
5191	ENDIF
5192
5193	ll_grid0F=.FALSE.
5194	IF( ASSOCIATED(tl_lon0F%d_value) .AND. &
5195	& ASSOCIATED(tl_lat0F%d_value) )THEN
5196	ll_grid0F=.TRUE.
5197	ENDIF
5198
5199	ENDIF
5200
5201	! close mpp files
5202	CALL iom_dom_close(tl_coord0)
5203
5204	! clean structure
5205	CALL mpp_clean(tl_coord0)
5206	CALL dom_clean(tl_dom0)
5207
5208	! copy structure
5209	tl_coord1=mpp_copy(td_coord1)
5210
5211	! open mpp files
5212	CALL iom_mpp_open(tl_coord1)
5213
5214	! read fine longitue and latitude
5215	il_ind=var_get_index(tl_coord1%t_proc(1)%t_var(:), TRIM(tl_lon0%c_longname))
5216	IF( il_ind /= 0 )THEN
5217	tl_lon1=iom_mpp_read_var(tl_coord1,TRIM(tl_lon0%c_longname))
5218	ELSE
5219	tl_lon1=iom_mpp_read_var(tl_coord1,'longitude')
5220	ENDIF
5221	il_ind=var_get_index(tl_coord1%t_proc(1)%t_var(:), TRIM(tl_lat0%c_longname))
5222	IF( il_ind /= 0 )THEN
5223	tl_lat1=iom_mpp_read_var(tl_coord1,TRIM(tl_lat0%c_longname))
5224	ELSE
5225	tl_lat1=iom_mpp_read_var(tl_coord1,'latitude')
5226	ENDIF
5227
5228	IF( ll_even )THEN
5229
5230	! look for variable value on domain for F point
5231	il_ind=var_get_index(tl_coord1%t_proc(1)%t_var(:), 'longitude_F')
5232	IF( il_ind /= 0 )THEN
5233	tl_lon1F=iom_mpp_read_var(tl_coord1,'longitude_F')
5234	ENDIF
5235
5236	il_ind=var_get_index(tl_coord1%t_proc(1)%t_var(:), 'latitude_F')
5237	IF( il_ind /= 0 )THEN
5238	tl_lat1F=iom_mpp_read_var(tl_coord1,'latitude_F')
5239	ENDIF
5240
5241	ll_grid1F=.FALSE.
5242	IF( ASSOCIATED(tl_lon1F%d_value) .AND. &
5243	& ASSOCIATED(tl_lat1F%d_value) )THEN
5244	ll_grid1F=.TRUE.
5245	ENDIF
5246
5247	ENDIF
5248
5249	! close mpp files
5250	CALL iom_mpp_close(tl_coord1)
5251	! clean structure
5252	CALL mpp_clean(tl_coord1)
5253
5254	CALL logger_debug("GRID CHECK COINCIDENCE:"//&
5255	& " fine grid "//TRIM(td_coord1%c_name) )
5256	CALL logger_debug("GRID CHECK COINCIDENCE:"//&
5257	& " coarse grid "//TRIM(td_coord0%c_name) )
5258
5259	! check domain
5260	! check global grid
5261	IF( .NOT. grid_is_global(tl_lon0, tl_lat0) )THEN
5262	IF( grid_is_global(tl_lon1, tl_lat1) )THEN
5263
5264	ll_coincidence=.FALSE.
5265	CALL logger_fatal("GRID CHECK COINCIDENCE:"//&
5266	& " fine grid is global,"//&
5267	& " coarse grid is not ")
5268
5269	ELSE
5270	il_ew1=tl_lon1%i_ew
5271	IF( il_ew1 >= 0 )THEN
5272	! ew overlap
5273
5274	il_ew0=tl_lon0%i_ew
5275	IF( il_ew0 < 0 )THEN
5276	CALL logger_fatal("GRID CHECK COINCIDENCE: "//&
5277	& "fine grid has east west overlap,"//&
5278	& " coarse grid not ")
5279	ENDIF
5280
5281	il_jmin1=1+ip_ghost
5282	il_jmax1=tl_lon1%t_dim(2)%i_len-ip_ghost
5283
5284	ll_coincidence=grid__check_lat(&
5285	& tl_lat0%d_value(1,:,1,1),&
5286	& tl_lat1%d_value(1,il_jmin1:il_jmax1,1,1))
5287
5288	ELSE
5289	! other case
5290	il_imin1=1+ip_ghost
5291	il_jmin1=1+ip_ghost
5292
5293	il_imax1=tl_lon1%t_dim(1)%i_len-ip_ghost
5294	il_jmax1=tl_lon1%t_dim(2)%i_len-ip_ghost
5295
5296	ll_coincidence=grid__check_corner(&
5297	& tl_lon0%d_value(:,:,1,1),&
5298	& tl_lat0%d_value(:,:,1,1),&
5299	& tl_lon1%d_value(il_imin1:il_imax1, &
5300	& il_jmin1:il_jmax1, &
5301	& 1,1),&
5302	& tl_lat1%d_value(il_imin1:il_imax1, &
5303	& il_jmin1:il_jmax1, &
5304	& 1,1) )
5305
5306	ENDIF
5307
5308	ENDIF
5309
5310	IF( .NOT. ll_coincidence )THEN
5311	CALL logger_fatal("GRID CHECK COINCIDENCE: no coincidence "//&
5312	& "between fine grid and coarse grid: invalid domain." )
5313	ENDIF
5314
5315	ENDIF
5316
5317	! check refinement factor
5318	! select point in middle of fine grid
5319	il_imid1=INT(tl_lon1%t_dim(1)%i_len*0.5)
5320	il_jmid1=INT(tl_lon1%t_dim(2)%i_len*0.5)
5321
5322	dl_lon1=tl_lon1%d_value(il_imid1, il_jmid1,1,1)
5323	dl_lat1=tl_lat1%d_value(il_imid1, il_jmid1,1,1)
5324
5325	! select closest point on coarse grid
5326	il_ind0(:)=grid_get_closest(tl_lon0%d_value(:,:,1,1),&
5327	& tl_lat0%d_value(:,:,1,1),&
5328	& dl_lon1, dl_lat1 )
5329
5330	IF( ANY(il_ind0(:)==0) )THEN
5331	CALL logger_fatal("GRID CHECK COINCIDENCE: can not find valid "//&
5332	& "coarse grid indices: invalid domain." )
5333	ENDIF
5334
5335	IF( .NOT. ll_even )THEN
5336	! case odd refinment in both direction
5337	! work on T-point
5338
5339	dl_lon0=tl_lon0%d_value(il_ind0(1),il_ind0(2),1,1)
5340	dl_lat0=tl_lat0%d_value(il_ind0(1),il_ind0(2),1,1)
5341
5342	il_ind1(:)=grid_get_closest(tl_lon1%d_value(:,:,1,1),&
5343	& tl_lat1%d_value(:,:,1,1),&
5344	& dl_lon0, dl_lat0 )
5345
5346	! check i-direction refinement factor
5347	DO ji=0,MIN(3,il_imid1)
5348
5349	IF( il_ind1(1)+ji*id_rho(jp_I)+1 > tl_lon1%t_dim(1)%i_len )THEN
5350	CALL logger_warn("GRID CHECK COINCIDENCE: domain to small "//&
5351	& " to check i-direction refinement factor ")
5352	EXIT
5353	ELSE
5354	dl_lon0=tl_lon0%d_value(il_ind0(1)+ji ,il_ind0(2),1,1)
5355	dl_lon1=tl_lon1%d_value(il_ind1(1)+ji*id_rho(jp_I),il_ind1(2),1,1)
5356
5357	! assume there could be little difference due to interpolation
5358	IF( ABS(dl_lon1 - dl_lon0) > dp_delta )THEN
5359	ll_coincidence=.FALSE.
5360	CALL logger_debug("GRID CHECK COINCIDENCE: invalid "//&
5361	& "i-direction refinement factor ("//&
5362	& TRIM(fct_str(id_rho(jp_I)))//&
5363	& ") between fine grid and coarse grid ")
5364	ENDIF
5365	ENDIF
5366
5367	ENDDO
5368
5369	! check j-direction refinement factor
5370	DO jj=0,MIN(3,il_jmid1)
5371
5372	IF( il_ind1(2)+jj*id_rho(jp_J)+1 > tl_lat1%t_dim(2)%i_len )THEN
5373	CALL logger_warn("GRID CHECK COINCIDENCE: domain to small "//&
5374	& " to check j-direction refinement factor ")
5375	EXIT
5376	ELSE
5377	dl_lat0=tl_lat0%d_value(il_ind0(1),il_ind0(2)+jj ,1,1)
5378	dl_lat1=tl_lat1%d_value(il_ind1(1),il_ind1(2)+jj*id_rho(jp_J),1,1)
5379
5380	! assume there could be little difference due to interpolation
5381	IF( ABS(dl_lat1-dl_lat0) > dp_delta )THEN
5382	ll_coincidence=.FALSE.
5383	CALL logger_debug("GRID CHECK COINCIDENCE: invalid "//&
5384	& "j-direction refinement factor ("//&
5385	& TRIM(fct_str(id_rho(jp_J)))//&
5386	& ") between fine grid and coarse grid ")
5387	ENDIF
5388	ENDIF
5389
5390	ENDDO
5391
5392	ELSE
5393	! case even refinment at least in one direction
5394	! work on F-point
5395
5396	dl_delta=dp_delta
5397	! look for lower left fine point in coarse cell.
5398	IF( ll_grid0F )THEN
5399
5400	! lower left corner of coarse cell
5401	dl_lon0F=tl_lon0F%d_value(il_ind0(1)-1,il_ind0(2)-1,1,1)
5402	dl_lat0F=tl_lat0F%d_value(il_ind0(1)-1,il_ind0(2)-1,1,1)
5403
5404	ELSE
5405
5406	! approximate lower left corner of coarse cell (with T point)
5407	dl_lon0F=( tl_lon0%d_value(il_ind0(1) ,il_ind0(2) ,1,1) + &
5408	& tl_lon0%d_value(il_ind0(1) ,il_ind0(2)-1,1,1) + &
5409	& tl_lon0%d_value(il_ind0(1)-1,il_ind0(2) ,1,1) + &
5410	& tl_lon0%d_value(il_ind0(1)-1,il_ind0(2)-1,1,1) ) * 0.25
5411
5412	dl_lat0F=( tl_lat0%d_value(il_ind0(1) ,il_ind0(2) ,1,1) + &
5413	& tl_lat0%d_value(il_ind0(1) ,il_ind0(2)-1,1,1) + &
5414	& tl_lat0%d_value(il_ind0(1)-1,il_ind0(2) ,1,1) + &
5415	& tl_lat0%d_value(il_ind0(1)-1,il_ind0(2)-1,1,1) ) * 0.25
5416
5417	! as we use approximation of F-point we relax condition
5418	dl_delta=100*dp_delta
5419
5420	ENDIF
5421
5422	IF( ll_grid1F )THEN
5423
5424	il_ind1(:)=grid_get_closest(tl_lon1F%d_value(:,:,1,1),&
5425	& tl_lat1F%d_value(:,:,1,1),&
5426	& dl_lon0F, dl_lat0F )
5427
5428	ELSE
5429
5430	il_ill1(:)=grid_get_closest(tl_lon1%d_value(:,:,1,1),&
5431	& tl_lat1%d_value(:,:,1,1),&
5432	& dl_lon0F, dl_lat0F, 'll' )
5433
5434	il_ilr1(:)=grid_get_closest(tl_lon1%d_value(:,:,1,1),&
5435	& tl_lat1%d_value(:,:,1,1),&
5436	& dl_lon0F, dl_lat0F, 'lr' )
5437
5438	il_iul1(:)=grid_get_closest(tl_lon1%d_value(:,:,1,1),&
5439	& tl_lat1%d_value(:,:,1,1),&
5440	& dl_lon0F, dl_lat0F, 'ul' )
5441
5442	il_iur1(:)=grid_get_closest(tl_lon1%d_value(:,:,1,1),&
5443	& tl_lat1%d_value(:,:,1,1),&
5444	& dl_lon0F, dl_lat0F, 'ur' )
5445
5446	! as we use approximation of F-point we relax condition
5447	dl_delta=100*dp_delta
5448
5449	ENDIF
5450
5451	! check i-direction refinement factor
5452	DO ji=0,MIN(3,il_imid1)
5453
5454	IF( il_ind1(1)+ji*id_rho(jp_I)+1 > tl_lon1%t_dim(1)%i_len )THEN
5455	CALL logger_warn("GRID CHECK COINCIDENCE: domain to small "//&
5456	& " to check i-direction refinement factor ")
5457	EXIT
5458	ELSE
5459	IF( ll_grid0F )THEN
5460	dl_lon0F=tl_lon0F%d_value(il_ind0(1)+ji-1, il_ind0(2)-1,1,1)
5461	ELSE
5462	dl_lon0F= 0.25 * &
5463	& ( tl_lon0%d_value(il_ind0(1)+ji , il_ind0(2) ,1,1) + &
5464	& tl_lon0%d_value(il_ind0(1)+ji-1, il_ind0(2) ,1,1) + &
5465	& tl_lon0%d_value(il_ind0(1)+ji , il_ind0(2)-1,1,1) + &
5466	& tl_lon0%d_value(il_ind0(1)+ji-1, il_ind0(2)-1,1,1) )
5467	ENDIF
5468
5469	IF( ll_grid1F )THEN
5470	dl_lon1F= tl_lon1F%d_value( il_ind1(1)+ji*id_rho(jp_I), &
5471	& il_ind1(2),1,1)
5472	ELSE
5473	dl_lon1F= 0.25 * &
5474	& ( tl_lon1%d_value( il_ill1(1)+ji*id_rho(jp_I), &
5475	& il_ill1(2),1,1) + &
5476	& tl_lon1%d_value( il_ilr1(1)+ji*id_rho(jp_I), &
5477	& il_ilr1(2),1,1) + &
5478	& tl_lon1%d_value( il_iul1(1)+ji*id_rho(jp_I), &
5479	& il_iul1(2),1,1) + &
5480	& tl_lon1%d_value( il_iur1(1)+ji*id_rho(jp_I), &
5481	& il_iur1(2),1,1) )
5482
5483	ENDIF
5484
5485	! assume there could be little difference due to interpolation
5486	IF( ABS(dl_lon1F - dl_lon0F) > dl_delta )THEN
5487	ll_coincidence=.FALSE.
5488	CALL logger_debug("GRID CHECK COINCIDENCE: invalid "//&
5489	& "i-direction refinement factor ("//&
5490	& TRIM(fct_str(id_rho(jp_I)))//&
5491	& ") between fine grid and coarse grid ")
5492	ENDIF
5493	ENDIF
5494
5495	ENDDO
5496
5497	! check j-direction refinement factor
5498	DO jj=0,MIN(3,il_jmid1)
5499
5500	IF( il_ind1(2)+jj*id_rho(jp_J)+1 > tl_lat1%t_dim(2)%i_len )THEN
5501	CALL logger_warn("GRID CHECK COINCIDENCE: domain to small "//&
5502	& " to check j-direction refinement factor ")
5503	EXIT
5504	ELSE
5505	IF( ll_grid0F )THEN
5506	dl_lat0F=tl_lat0F%d_value(il_ind0(1)-1, il_ind0(2)+jj-1,1,1)
5507	ELSE
5508	dl_lat0F= 0.25 * &
5509	& ( tl_lat0%d_value(il_ind0(1) , il_ind0(2)+jj ,1,1) + &
5510	& tl_lat0%d_value(il_ind0(1)-1, il_ind0(2)+jj ,1,1) + &
5511	& tl_lat0%d_value(il_ind0(1) , il_ind0(2)+jj-1,1,1) + &
5512	& tl_lat0%d_value(il_ind0(1)-1, il_ind0(2)+jj-1,1,1) )
5513	ENDIF
5514
5515	IF( ll_grid1F )THEN
5516	dl_lat1F= tl_lat1F%d_value( il_ind1(1), &
5517	& il_ind1(2)+jj*id_rho(jp_J),1,1)
5518	ELSE
5519	dl_lat1F= 0.25 * &
5520	& ( tl_lat1%d_value( il_ill1(1), &
5521	& il_ill1(2)+jj*id_rho(jp_J),1,1) + &
5522	& tl_lat1%d_value( il_ilr1(1), &
5523	& il_ilr1(2)+jj*id_rho(jp_J),1,1) + &
5524	& tl_lat1%d_value( il_iul1(1), &
5525	& il_iul1(2)+jj*id_rho(jp_J),1,1) + &
5526	& tl_lat1%d_value( il_iur1(1), &
5527	& il_iur1(2)+jj*id_rho(jp_J),1,1) )
5528
5529	ENDIF
5530
5531	! assume there could be little difference due to interpolation
5532	IF( ABS(dl_lat1F - dl_lat0F) > dl_delta )THEN
5533	ll_coincidence=.FALSE.
5534	CALL logger_debug("GRID CHECK COINCIDENCE: invalid "//&
5535	& "i-direction refinement factor ("//&
5536	& TRIM(fct_str(id_rho(jp_I)))//&
5537	& ") between fine grid and coarse grid ")
5538	ENDIF
5539	ENDIF
5540
5541	ENDDO
5542	ENDIF
5543
5544	! clean
5545	CALL var_clean(tl_lon1)
5546	CALL var_clean(tl_lat1)
5547	CALL var_clean(tl_lon0)
5548	CALL var_clean(tl_lat0)
5549
5550	IF( .NOT. ll_coincidence )THEN
5551	CALL logger_fatal("GRID CHECK COINCIDENCE: no coincidence "//&
5552	& "between fine and coarse grid: "//&
5553	& "invalid refinement factor" )
5554	ENDIF
5555
5556	END SUBROUTINE grid_check_coincidence
5557	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
5558	FUNCTION grid__check_corner(dd_lon0, dd_lat0, &
5559	& dd_lon1, dd_lat1) &
5560	& RESULT (lf_inside)
5561	!-------------------------------------------------------------------
5562	!> @brief This function check that fine grid is
5563	!> inside coarse grid
5564	!>
5565	!> @details
5566	!>
5567	!> @author J.Paul
5568	!> @date November, 2013 - Initial Version
5569	!>
5570	!> @param[in] dd_lon0 array of coarse grid longitude
5571	!> @param[in] dd_lat0 array of coarse grid latitude
5572	!> @param[in] dd_lon1 array of fine grid longitude
5573	!> @param[in] dd_lat1 array of fine grid latitude
5574	!> @return true if fine grid is inside coarse grid
5575	!-------------------------------------------------------------------
5576
5577	IMPLICIT NONE
5578
5579	! Argument
5580	REAL(dp), DIMENSION(:,:), INTENT(IN) :: dd_lon0
5581	REAL(dp), DIMENSION(:,:), INTENT(IN) :: dd_lat0
5582	REAL(dp), DIMENSION(:,:), INTENT(IN) :: dd_lon1
5583	REAL(dp), DIMENSION(:,:), INTENT(IN) :: dd_lat1
5584
5585	! function
5586	LOGICAL :: lf_inside
5587
5588	! local variable
5589	INTEGER(i4), DIMENSION(2) :: il_shape0
5590	INTEGER(i4), DIMENSION(2) :: il_shape1
5591
5592	INTEGER(i4) :: il_imin0
5593	INTEGER(i4) :: il_jmin0
5594	INTEGER(i4) :: il_imax0
5595	INTEGER(i4) :: il_jmax0
5596
5597	INTEGER(i4) :: il_imin1
5598	INTEGER(i4) :: il_jmin1
5599	INTEGER(i4) :: il_imax1
5600	INTEGER(i4) :: il_jmax1
5601
5602	REAL(dp) :: dl_lon0
5603	REAL(dp) :: dl_lat0
5604
5605	REAL(dp) :: dl_lon1
5606	REAL(dp) :: dl_lat1
5607	! loop indices
5608	!----------------------------------------------------------------
5609
5610	! init
5611	lf_inside=.TRUE.
5612
5613	il_shape0=SHAPE(dd_lon0(:,:))
5614	il_shape1=SHAPE(dd_lon1(:,:))
5615
5616	!1- check if fine grid inside coarse grid domain
5617	il_imin0=1 ; il_imax0=il_shape0(1)
5618	il_jmin0=1 ; il_jmax0=il_shape0(2)
5619
5620	il_imin1=1 ; il_imax1=il_shape1(1)
5621	il_jmin1=1 ; il_jmax1=il_shape1(2)
5622
5623	! check lower left corner
5624	dl_lon0 = dd_lon0(il_imin0, il_jmin0)
5625	dl_lat0 = dd_lat0(il_imin0, il_jmin0)
5626
5627	dl_lon1 = dd_lon1(il_imin1, il_jmin1)
5628	dl_lat1 = dd_lat1(il_imin1, il_jmin1)
5629
5630	IF( (ABS(dl_lon1-dl_lon0)>dp_delta) .AND. (dl_lon1 < dl_lon0 ) .OR. &
5631	&(ABS(dl_lat1-dl_lat0)>dp_delta) .AND. (dl_lat1 < dl_lat0 ) )THEN
5632
5633	CALL logger_error("GRID CHECK COINCIDENCE: fine grid lower left "//&
5634	& "corner not north east of coarse grid (imin,jmin) ")
5635	CALL logger_debug(" fine grid lower left ( "//&
5636	& TRIM(fct_str(dl_lon1))//","//&
5637	& TRIM(fct_str(dl_lat1))//")" )
5638	CALL logger_debug(" coarse grid lower left ( "//&
5639	& TRIM(fct_str(dl_lon0))//","//&
5640	& TRIM(fct_str(dl_lat0))//")" )
5641	lf_inside=.FALSE.
5642
5643	ENDIF
5644
5645	! check upper left corner
5646	dl_lon0 = dd_lon0(il_imin0, il_jmax0)
5647	dl_lat0 = dd_lat0(il_imin0, il_jmax0)
5648
5649	dl_lon1 = dd_lon1(il_imin1, il_jmax1)
5650	dl_lat1 = dd_lat1(il_imin1, il_jmax1)
5651
5652	IF( (ABS(dl_lon1-dl_lon0)>dp_delta) .AND. (dl_lon1 < dl_lon0) .OR. &
5653	& (ABS(dl_lat1-dl_lat0)>dp_delta) .AND. (dl_lat1 > dl_lat0) )THEN
5654
5655	CALL logger_error("GRID CHECK COINCIDENCE: fine grid upper left "//&
5656	& "corner not south east of coarse grid (imin,jmax) ")
5657	CALL logger_debug(" fine grid upper left ("//&
5658	& TRIM(fct_str(dl_lon1))//","//&
5659	& TRIM(fct_str(dl_lat1))//")")
5660	CALL logger_debug(" coasre grid upper left ("//&
5661	& TRIM(fct_str(dl_lon0))//","//&
5662	& TRIM(fct_str(dl_lat0))//")")
5663	lf_inside=.FALSE.
5664
5665	ENDIF
5666
5667	! check lower right corner
5668	dl_lon0 = dd_lon0(il_imax0, il_jmin0)
5669	dl_lat0 = dd_lat0(il_imax0, il_jmin0)
5670
5671	dl_lon1 = dd_lon1(il_imax1, il_jmin1)
5672	dl_lat1 = dd_lat1(il_imax1, il_jmin1)
5673
5674
5675	IF( (ABS(dl_lon1-dl_lon0)>dp_delta) .AND. (dl_lon1 > dl_lon0) .OR. &
5676	&(ABS(dl_lat1-dl_lat0)>dp_delta) .AND. (dl_lat1 < dl_lat0) )THEN
5677
5678	CALL logger_error("GRID CHECK COINCIDENCE: fine grid lower right "//&
5679	& "corner not north west of coarse grid (imax,jmin) ")
5680	CALL logger_debug(" fine grid lower right ( "//&
5681	& TRIM(fct_str(dl_lon1))//","//&
5682	& TRIM(fct_str(dl_lat1))//")" )
5683	CALL logger_debug(" coarse grid lower right ( "//&
5684	& TRIM(fct_str(dl_lon0))//","//&
5685	& TRIM(fct_str(dl_lat0))//")" )
5686	lf_inside=.FALSE.
5687
5688	ENDIF
5689
5690	! check upper right corner
5691	dl_lon0 = dd_lon0(il_imax0, il_jmax0)
5692	dl_lat0 = dd_lat0(il_imax0, il_jmax0)
5693
5694	dl_lon1 = dd_lon1(il_imax1, il_jmax1)
5695	dl_lat1 = dd_lat1(il_imax1, il_jmax1)
5696
5697	IF( (ABS(dl_lon1-dl_lon0)>dp_delta) .AND. (dl_lon1 > dl_lon0) .OR. &
5698	& (ABS(dl_lat1-dl_lat0)>dp_delta) .AND. (dl_lat1 > dl_lat0) )THEN
5699
5700	CALL logger_error("GRID CHECK COINCIDENCE: fine grid upper right "//&
5701	& "corner not south west of coarse grid (imax,jmax) ")
5702	CALL logger_debug(" fine grid upper right ( "//&
5703	& TRIM(fct_str(dl_lon1))//","//&
5704	& TRIM(fct_str(dl_lat1))//")" )
5705	CALL logger_debug(" fine imax1 jmax1 ( "//&
5706	& TRIM(fct_str(il_imax1))//","//&
5707	& TRIM(fct_str(il_jmax1))//")" )
5708	CALL logger_debug(" coarse grid upper right ( "//&
5709	& TRIM(fct_str(dl_lon0))//","//&
5710	& TRIM(fct_str(dl_lat0))//")" )
5711	CALL logger_debug(" fine imax0 jmax0 ( "//&
5712	& TRIM(fct_str(il_imax0))//","//&
5713	& TRIM(fct_str(il_jmax0))//")" )
5714	lf_inside=.FALSE.
5715
5716	ENDIF
5717
5718	END FUNCTION grid__check_corner
5719	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
5720	FUNCTION grid__check_lat(dd_lat0, dd_lat1) &
5721	& RESULT (lf_inside)
5722	!-------------------------------------------------------------------
5723	!> @brief This function check that fine grid latitude are
5724	!> inside coarse grid latitude
5725	!>
5726	!> @details
5727	!>
5728	!> @author J.Paul
5729	!> @date November, 2013 - Initial Version
5730	!>
5731	!> @param[in] dd_lat0 array of coarse grid latitude
5732	!> @param[in] dd_lat1 array of fine grid latitude
5733	!> @return true if fine grid is inside coarse grid
5734	!-------------------------------------------------------------------
5735
5736	IMPLICIT NONE
5737
5738	! Argument
5739	REAL(dp), DIMENSION(:), INTENT(IN) :: dd_lat0
5740	REAL(dp), DIMENSION(:), INTENT(IN) :: dd_lat1
5741
5742	! function
5743	LOGICAL :: lf_inside
5744
5745	! local variable
5746	INTEGER(i4), DIMENSION(1) :: il_shape0
5747	INTEGER(i4), DIMENSION(1) :: il_shape1
5748
5749	INTEGER(i4) :: il_jmin0
5750	INTEGER(i4) :: il_jmax0
5751
5752	INTEGER(i4) :: il_jmin1
5753	INTEGER(i4) :: il_jmax1
5754	! loop indices
5755	!----------------------------------------------------------------
5756
5757	! init
5758	lf_inside=.TRUE.
5759
5760	il_shape0(:)=SHAPE(dd_lat0(:))
5761	il_shape1(:)=SHAPE(dd_lat1(:))
5762
5763	!1- check if fine grid inside coarse grid domain
5764	il_jmin0=1 ; il_jmax0=il_shape0(1)
5765	il_jmin1=1 ; il_jmax1=il_shape1(1)
5766
5767	! check lower left fine grid
5768	IF( ABS(dd_lat1(il_jmin1)-dd_lat0(il_jmin0)) > dp_delta .AND. &
5769	&dd_lat1(il_jmin1) < dd_lat0(il_jmin0) )THEN
5770
5771	CALL logger_error("GRID CHECK COINCIDENCE: fine grid lower point"//&
5772	& " not north of coarse grid (jmin) ")
5773	CALL logger_debug(" fine grid lower point ( "//&
5774	& TRIM(fct_str(dd_lat1(il_jmin1)))//")" )
5775	CALL logger_debug(" coarse grid lower point ( "//&
5776	& TRIM(fct_str(dd_lat0(il_jmin0)))//")" )
5777	lf_inside=.FALSE.
5778
5779	ENDIF
5780
5781	! check upper left fine grid
5782	IF( ABS(dd_lat1(il_jmax1)-dd_lat0(il_jmax0)) > dp_delta .AND. &
5783	&dd_lat1(il_jmax1) > dd_lat0(il_jmax0) )THEN
5784
5785	CALL logger_error("GRID CHECK COINCIDENCE: fine grid upper point"//&
5786	& " not south of coarse grid (jmax) ")
5787	CALL logger_debug(" fine grid upper point ("//&
5788	& TRIM(fct_str(dd_lat1(il_jmax1)))//")")
5789	CALL logger_debug(" coasre grid upper point ("//&
5790	& TRIM(fct_str(dd_lat0(il_jmax0)))//")")
5791	lf_inside=.FALSE.
5792
5793	ENDIF
5794
5795	END FUNCTION grid__check_lat
5796	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
5797	SUBROUTINE grid_add_ghost(td_var, id_ghost)
5798	!-------------------------------------------------------------------
5799	!> @brief
5800	!> This subroutine add ghost cell at boundaries.
5801	!>
5802	!> @author J.Paul
5803	!> @date November, 2013 - Initial version
5804	!>
5805	!> @param[inout] td_var array of variable structure
5806	!> @param[in] id_ghost array of ghost cell factor
5807	!-------------------------------------------------------------------
5808
5809	IMPLICIT NONE
5810
5811	! Argument
5812	TYPE(TVAR) , INTENT(INOUT) :: td_var
5813	INTEGER(i4), DIMENSION(2,2), INTENT(IN ) :: id_ghost
5814
5815	! local variable
5816	INTEGER(i4) :: il_imin
5817	INTEGER(i4) :: il_jmin
5818	INTEGER(i4) :: il_imax
5819	INTEGER(i4) :: il_jmax
5820
5821	REAL(dp), DIMENSION(:,:,:,:) , ALLOCATABLE :: dl_value
5822
5823	TYPE(TVAR) :: tl_var
5824
5825	! loop indices
5826	!----------------------------------------------------------------
5827
5828	IF( ALL(td_var%t_dim(1:2)%l_use) )THEN
5829
5830	CALL logger_warn( "ADD GHOST: dimension change in variable "//&
5831	& TRIM(td_var%c_name) )
5832
5833	! copy variable
5834	tl_var=var_copy(td_var)
5835
5836	CALL var_del_value(td_var)
5837
5838	! compute indice to fill center
5839	il_imin=1+id_ghost(jp_I,1)*ip_ghost
5840	il_jmin=1+id_ghost(jp_J,1)*ip_ghost
5841
5842	il_imax=tl_var%t_dim(1)%i_len+id_ghost(jp_I,1)*ip_ghost
5843	il_jmax=tl_var%t_dim(2)%i_len+id_ghost(jp_J,1)*ip_ghost
5844
5845	! compute new dimension
5846	td_var%t_dim(1)%i_len = tl_var%t_dim(1)%i_len + &
5847	& SUM(id_ghost(jp_I,:))*ip_ghost
5848	td_var%t_dim(2)%i_len = tl_var%t_dim(2)%i_len + &
5849	& SUM(id_ghost(jp_J,:))*ip_ghost
5850
5851	ALLOCATE(dl_value(td_var%t_dim(1)%i_len, &
5852	& td_var%t_dim(2)%i_len, &
5853	& td_var%t_dim(3)%i_len, &
5854	& td_var%t_dim(4)%i_len) )
5855
5856	dl_value(:,:,:,:)=tl_var%d_fill
5857
5858	dl_value(il_imin:il_imax, &
5859	& il_jmin:il_jmax, &
5860	& :,:) = tl_var%d_value(:,:,:,:)
5861
5862	! add variable value
5863	CALL var_add_value(td_var,dl_value(:,:,:,:))
5864
5865	! save variable type
5866	td_var%i_type=tl_var%i_type
5867
5868	DEALLOCATE( dl_value )
5869
5870	CALL var_clean(tl_var)
5871
5872	ENDIF
5873
5874	END SUBROUTINE grid_add_ghost
5875	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
5876	SUBROUTINE grid_del_ghost(td_var, id_ghost)
5877	!-------------------------------------------------------------------
5878	!> @brief
5879	!> This subroutine delete ghost cell at boundaries.
5880	!>
5881	!> @author J.Paul
5882	!> @date November, 2013 - Initial version
5883	!>
5884	!> @param[inout] td_var array of variable structure
5885	!> @param[in] id_ghost array of ghost cell factor
5886	!-------------------------------------------------------------------
5887
5888	IMPLICIT NONE
5889
5890	! Argument
5891	TYPE(TVAR) , INTENT(INOUT) :: td_var
5892	INTEGER(i4), DIMENSION(2,2), INTENT(IN ) :: id_ghost
5893
5894	! local variable
5895	INTEGER(i4) :: il_imin
5896	INTEGER(i4) :: il_jmin
5897	INTEGER(i4) :: il_imax
5898	INTEGER(i4) :: il_jmax
5899
5900	REAL(dp), DIMENSION(:,:,:,:) , ALLOCATABLE :: dl_value
5901
5902	TYPE(TVAR) :: tl_var
5903
5904	! loop indices
5905	!----------------------------------------------------------------
5906
5907	IF( ALL(td_var%t_dim(1:2)%l_use) )THEN
5908
5909	IF( ANY(id_ghost(:,:)/=0) )THEN
5910	CALL logger_warn( "GRID DEL GHOST: dimension change in variable "//&
5911	& TRIM(td_var%c_name) )
5912	ENDIF
5913
5914	! copy variable
5915	tl_var=var_copy(td_var)
5916
5917	CALL var_del_value(td_var)
5918
5919	! compute indice to get center
5920	il_imin=1+id_ghost(jp_I,1)*ip_ghost
5921	il_jmin=1+id_ghost(jp_J,1)*ip_ghost
5922
5923	il_imax=tl_var%t_dim(1)%i_len-id_ghost(jp_I,2)*ip_ghost
5924	il_jmax=tl_var%t_dim(2)%i_len-id_ghost(jp_J,2)*ip_ghost
5925
5926	! compute new dimension
5927	td_var%t_dim(1)%i_len = il_imax - il_imin +1
5928	td_var%t_dim(2)%i_len = il_jmax - il_jmin +1
5929
5930	ALLOCATE(dl_value(td_var%t_dim(1)%i_len, &
5931	& td_var%t_dim(2)%i_len, &
5932	& td_var%t_dim(3)%i_len, &
5933	& td_var%t_dim(4)%i_len) )
5934
5935	dl_value(:,:,:,:)=tl_var%d_fill
5936
5937	dl_value(:,:,:,:) = tl_var%d_value(il_imin:il_imax, &
5938	& il_jmin:il_jmax, &
5939	& :,:)
5940
5941	! add variable value
5942	CALL var_add_value(td_var,dl_value(:,:,:,:))
5943
5944	! save variable type
5945	td_var%i_type=tl_var%i_type
5946
5947	DEALLOCATE( dl_value )
5948
5949	CALL var_clean(tl_var)
5950
5951	ENDIF
5952
5953	END SUBROUTINE grid_del_ghost
5954	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
5955	FUNCTION grid__get_ghost_var(td_var) &
5956	& RESULT (if_ghost)
5957	!-------------------------------------------------------------------
5958	!> @brief This function check if ghost cell are used or not, and return ghost
5959	!> cell factor (0,1) in horizontal plan.
5960	!>
5961	!> @details
5962	!> check if domain is global, and if there is an East-West overlap.
5963	!>
5964	!> @author J.Paul
5965	!> @date September, 2014 - Initial Version
5966	!>
5967	!> @param[in] td_var variable sturcture
5968	!> @return array of ghost cell factor
5969	!-------------------------------------------------------------------
5970
5971	IMPLICIT NONE
5972
5973	! Argument
5974	TYPE(TVAR), INTENT(IN) :: td_var
5975
5976	! function
5977	INTEGER(i4), DIMENSION(2,2) :: if_ghost
5978
5979	! local variable
5980	INTEGER(i4), DIMENSION(ip_maxdim) :: il_dim
5981
5982	! loop indices
5983	!----------------------------------------------------------------
5984	! init
5985	if_ghost(:,:)=0
5986
5987	IF( .NOT. ALL(td_var%t_dim(1:2)%l_use) )THEN
5988	CALL logger_error("GRID GET GHOST: "//TRIM(td_var%c_name)//" is not a suitable"//&
5989	& " variable to look for ghost cell (not 2D).")
5990	ELSE
5991	IF( .NOT. ASSOCIATED(td_var%d_value) )THEN
5992	CALL logger_error("GRID GET GHOST: no value associated to "//TRIM(td_var%c_name)//&
5993	& ". can't look for ghost cell.")
5994	ELSE
5995	il_dim(:)=td_var%t_dim(:)%i_len
5996
5997	IF(ALL(td_var%d_value( 1 , : ,1,1)/=td_var%d_fill).AND.&
5998	& ALL(td_var%d_value(il_dim(1), : ,1,1)/=td_var%d_fill).AND.&
5999	& ALL(td_var%d_value( : , 1 ,1,1)/=td_var%d_fill).AND.&
6000	& ALL(td_var%d_value( : ,il_dim(2),1,1)/=td_var%d_fill))THEN
6001	! no boundary closed
6002	CALL logger_warn("GRID GET GHOST: can't determined ghost cell. "//&
6003	& "there is no boundary closed for variable "//&
6004	& TRIM(td_var%c_name))
6005
6006	ELSE
6007	! check periodicity
6008	IF(ANY(td_var%d_value( 1 ,:,1,1)/=td_var%d_fill).OR.&
6009	& ANY(td_var%d_value(il_dim(1),:,1,1)/=td_var%d_fill))THEN
6010	! East-West cyclic (1,4,6)
6011	CALL logger_info("GRID GET GHOST: East West cyclic")
6012	if_ghost(jp_I,:)=0
6013
6014	IF( ANY(td_var%d_value(:, 1, 1, 1) /= td_var%d_fill) )THEN
6015	! South boundary not closed
6016
6017	CALL logger_debug("GRID GET GHOST: East_West cyclic")
6018	CALL logger_debug("GRID GET GHOST: South boundary not closed")
6019	CALL logger_error("GRID GET GHOST: should have been an "//&
6020	& "impossible case")
6021
6022	ELSE
6023	! South boundary closed (1,4,6)
6024	CALL logger_info("GRID GET GHOST: South boundary closed")
6025	if_ghost(jp_J,1)=1
6026
6027	IF(ANY(td_var%d_value(:,il_dim(2),1,1)/=td_var%d_fill) )THEN
6028	! North boundary not closed (4,6)
6029	CALL logger_info("GRID GET GHOST: North boundary not closed")
6030	if_ghost(jp_J,2)=0
6031	ELSE
6032	! North boundary closed
6033	CALL logger_info("GRID GET GHOST: North boundary closed")
6034	if_ghost(jp_J,2)=1
6035	ENDIF
6036
6037	ENDIF
6038
6039	ELSE
6040	! East-West boundaries closed (0,2,3,5)
6041	CALL logger_info("GRID GET GHOST: East West boundaries closed")
6042	if_ghost(jp_I,:)=1
6043
6044	IF( ANY(td_var%d_value(:, 1, 1, 1) /= td_var%d_fill) )THEN
6045	! South boundary not closed (2)
6046	CALL logger_info("GRID GET GHOST: South boundary not closed")
6047	if_ghost(jp_J,1)=0
6048
6049	IF(ANY(td_var%d_value(:,il_dim(2),1,1)/=td_var%d_fill))THEN
6050	! North boundary not closed
6051	CALL logger_debug("GRID GET GHOST: East West boundaries "//&
6052	& "closed")
6053	CALL logger_debug("GRID GET GHOST: South boundary not closed")
6054	CALL logger_debug("GRID GET GHOST: North boundary not closed")
6055	CALL logger_error("GRID GET GHOST: should have been "//&
6056	& "an impossible case")
6057	ELSE
6058	! North boundary closed
6059	if_ghost(jp_J,2)=1
6060	ENDIF
6061
6062	ELSE
6063	! South boundary closed (0,3,5)
6064	CALL logger_info("GRID GET GHOST: South boundary closed")
6065	if_ghost(jp_J,1)=1
6066
6067	IF(ANY(td_var%d_value(:,il_dim(2),1,1)/=td_var%d_fill))THEN
6068	! North boundary not closed (3,5)
6069	CALL logger_info("GRID GET GHOST: North boundary not closed")
6070	if_ghost(jp_J,2)=0
6071	ELSE
6072	! North boundary closed
6073	CALL logger_info("GRID GET GHOST: North boundary closed")
6074	if_ghost(jp_J,2)=1
6075	ENDIF
6076
6077	ENDIF
6078
6079	ENDIF
6080
6081	ENDIF
6082
6083	ENDIF
6084	ENDIF
6085
6086	END FUNCTION grid__get_ghost_var
6087	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
6088	FUNCTION grid__get_ghost_mpp(td_mpp) &
6089	& RESULT (if_ghost)
6090	!-------------------------------------------------------------------
6091	!> @brief This function check if ghost cell are used or not, and return ghost
6092	!> cell factor (0,1) in i- and j-direction.
6093	!>
6094	!> @details
6095	!> get longitude an latitude array, then
6096	!> check if domain is global, and if there is an East-West overlap
6097	!>
6098	!> @author J.Paul
6099	!> @date September, 2014 - Initial Version
6100	!> @date October, 2014
6101	!> - work on mpp file structure instead of file structure
6102	!>
6103	!> @param[in] td_file file sturcture
6104	!> @return array of ghost cell factor
6105	!-------------------------------------------------------------------
6106
6107	IMPLICIT NONE
6108
6109	! Argument
6110	TYPE(TMPP), INTENT(IN) :: td_mpp
6111
6112	! function
6113	INTEGER(i4), DIMENSION(2,2) :: if_ghost
6114
6115	! local variable
6116	!TYPE(TVAR) :: tl_lon
6117	!TYPE(TVAR) :: tl_lat
6118
6119	TYPE(TMPP) :: tl_mpp
6120
6121	!INTEGER(i4) :: il_lonid
6122	!INTEGER(i4) :: il_latid
6123	! loop indices
6124	!----------------------------------------------------------------
6125	! init
6126	if_ghost(:,:)=0
6127
6128	IF( .NOT. ASSOCIATED(td_mpp%t_proc) )THEN
6129	CALL logger_error("GRID GET GHOST: decomposition of mpp file "//&
6130	& TRIM(td_mpp%c_name)//" not defined." )
6131
6132	ELSE
6133
6134	! copy structure
6135	tl_mpp=mpp_copy(td_mpp)
6136
6137	CALL logger_info("GRID GET FINE GHOST perio"//TRIM(fct_str(tl_mpp%i_perio)))
6138	IF( tl_mpp%i_perio < 0 )THEN
6139	! compute NEMO periodicity index
6140	CALL grid_get_info(tl_mpp)
6141	ENDIF
6142
6143	SELECT CASE(tl_mpp%i_perio)
6144	CASE(0)
6145	if_ghost(:,:)=1
6146	CASE(1)
6147	if_ghost(jp_J,:)=1
6148	CASE(2)
6149	if_ghost(jp_I,:)=1
6150	if_ghost(jp_J,2)=1
6151	CASE(3,5)
6152	if_ghost(jp_I,:)=1
6153	if_ghost(jp_J,1)=1
6154	CASE(4,6)
6155	if_ghost(jp_J,1)=1
6156	CASE DEFAULT
6157	END SELECT
6158
6159	! clean
6160	CALL mpp_clean(tl_mpp)
6161
6162	ENDIF
6163
6164	END FUNCTION grid__get_ghost_mpp
6165	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
6166	FUNCTION grid_split_domain(td_var, id_level) &
6167	& RESULT (if_mask)
6168	!-------------------------------------------------------------------
6169	!> @brief This subroutine compute closed sea domain.
6170	!>
6171	!> @details
6172	!> to each domain is associated a negative value id (from -1 to ...)<br/>
6173	!> optionaly you could specify which level use (default 1)
6174	!>
6175	!> @author J.Paul
6176	!> @date November, 2013 - Initial Version
6177	!>
6178	!> @param[in] td_var variable strucutre
6179	!> @param[in] id_level level
6180	!> @return domain mask
6181	!-------------------------------------------------------------------
6182
6183	IMPLICIT NONE
6184
6185	! Argument
6186	TYPE(TVAR) , INTENT(IN) :: td_var
6187	INTEGER(i4), INTENT(IN), OPTIONAL :: id_level
6188
6189	! function
6190	INTEGER(i4), DIMENSION(td_var%t_dim(1)%i_len, &
6191	& td_var%t_dim(2)%i_len ) :: if_mask
6192
6193	! local variable
6194	INTEGER(i4) :: il_domid
6195	INTEGER(i4) :: il_level
6196	INTEGER(i4), DIMENSION(2) :: il_shape
6197	INTEGER(i4), DIMENSION(2) :: il_ind
6198	INTEGER(i4), DIMENSION(:,:), ALLOCATABLE :: il_tmp
6199
6200	LOGICAL :: ll_full
6201
6202	! loop indices
6203	INTEGER(i4) :: ji
6204	INTEGER(i4) :: jim
6205	INTEGER(i4) :: jip
6206	INTEGER(i4) :: jj
6207	INTEGER(i4) :: jjm
6208	INTEGER(i4) :: jjp
6209	!----------------------------------------------------------------
6210	il_level=1
6211	IF( PRESENT(id_level) ) il_level=id_level
6212
6213	! init
6214	il_domid=-1
6215
6216	il_shape(:)=td_var%t_dim(1:2)%i_len
6217	if_mask(:,:)=0
6218	WHERE( td_var%d_value(:,:,il_level,1)/=td_var%d_fill ) if_mask(:,:)=1
6219
6220	il_ind(:)=MAXLOC( if_mask(:,:) )
6221	DO WHILE( if_mask(il_ind(1), il_ind(2)) == 1 )
6222
6223	if_mask(il_ind(1),il_ind(2))=il_domid
6224	ll_full=.FALSE.
6225
6226	ALLOCATE( il_tmp(il_shape(1),il_shape(2)) )
6227
6228	DO WHILE( .NOT. ll_full )
6229	il_tmp(:,:)=0
6230
6231	ll_full=.TRUE.
6232	DO jj=1,il_shape(2)
6233	DO ji=1,il_shape(1)
6234	IF( if_mask(ji,jj)==il_domid )THEN
6235	jim=MAX(1,ji-1) ; jip=MIN(il_shape(1),ji+1)
6236	jjm=MAX(1,jj-1) ; jjp=MIN(il_shape(2),jj+1)
6237
6238	WHERE( if_mask(jim:jip,jjm:jjp)==1 )
6239	if_mask(jim:jip,jjm:jjp)=il_domid
6240	il_tmp(jim:jip,jjm:jjp)=1
6241	END WHERE
6242
6243	ENDIF
6244	ENDDO
6245	ENDDO
6246	IF( SUM(il_tmp(:,:))/=0 ) ll_full=.FALSE.
6247
6248	ENDDO
6249	DEALLOCATE( il_tmp )
6250
6251	il_ind(:)=MAXLOC( if_mask(:,:) )
6252	il_domid=il_domid-1
6253
6254	ENDDO
6255
6256	CALL logger_info("GRID SPLIT DOMAIN: "//TRIM( fct_str(ABS(il_domid+1)) )//&
6257	& " domain found" )
6258
6259	END FUNCTION grid_split_domain
6260	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
6261	SUBROUTINE grid_fill_small_dom(td_var, id_mask, id_minsize)
6262	!-------------------------------------------------------------------
6263	!> @brief This subroutine fill small closed sea with fill value.
6264	!>
6265	!> @details
6266	!> the minimum size (number of point) of closed sea to be kept could be
6267	!> sepcify with id_minsize.
6268	!> By default only the biggest sea is preserve.
6269	!>
6270	!> @author J.Paul
6271	!> @date November, 2013 - Initial Version
6272	!>
6273	!> @param[inout] td_var variable structure
6274	!> @param[in] id_mask domain mask (from grid_split_domain)
6275	!> @param[in] id_minsize minimum size of sea to be kept
6276	!-------------------------------------------------------------------
6277
6278	IMPLICIT NONE
6279
6280	! Argument
6281	TYPE(TVAR) , INTENT(INOUT) :: td_var
6282	INTEGER(i4), DIMENSION(:,:), INTENT(IN ) :: id_mask
6283	INTEGER(i4), INTENT(IN ), OPTIONAL :: id_minsize
6284
6285	! local variable
6286	INTEGER(i4) :: il_ndom
6287	INTEGER(i4) :: il_minsize
6288	INTEGER(i4), DIMENSION(2) :: il_shape
6289	INTEGER(i4), DIMENSION(:,:), ALLOCATABLE :: il_tmp
6290
6291	! loop indices
6292	INTEGER(i4) :: ji
6293	INTEGER(i4) :: jk
6294	INTEGER(i4) :: jl
6295	!----------------------------------------------------------------
6296
6297	il_shape(:)=SHAPE(id_mask(:,:))
6298	IF( ANY(il_shape(:) /= td_var%t_dim(1:2)%i_len) )THEN
6299	CALL logger_error("GRID FILL SMALL DOM: variable and mask "//&
6300	& "dimension differ")
6301	ELSE
6302
6303	il_ndom=MINVAL(id_mask(:,:))
6304
6305	ALLOCATE( il_tmp(il_shape(1),il_shape(2)) )
6306	il_tmp(:,:)=0
6307	DO ji=-1,il_ndom,-1
6308	WHERE( id_mask(:,:)==ji )
6309	il_tmp(:,:)=SUM(id_mask(:,:),id_mask(:,:)==ji)/ji
6310	END WHERE
6311	ENDDO
6312
6313	il_minsize=MAXVAL(il_tmp(:,:))
6314	IF( PRESENT(id_minsize) ) il_minsize=id_minsize
6315
6316	DO jl=1,td_var%t_dim(4)%i_len
6317	DO jk=1,td_var%t_dim(3)%i_len
6318	WHERE( il_tmp(:,:) < il_minsize )
6319	td_var%d_value(:,:,jk,jl)=td_var%d_fill
6320	END WHERE
6321	ENDDO
6322	ENDDO
6323
6324	DEALLOCATE( il_tmp )
6325
6326	ENDIF
6327
6328	END SUBROUTINE grid_fill_small_dom
6329	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
6330	SUBROUTINE grid_fill_small_msk(id_mask, id_minsize)
6331	!-------------------------------------------------------------------
6332	!> @brief This subroutine fill small domain inside bigger one.
6333	!>
6334	!> @details
6335	!> the minimum size (number of point) of domain sea to be kept
6336	!> is specified by id_minsize.
6337	!> smaller domain are included in the one they are embedded.
6338	!>
6339	!> @author J.Paul
6340	!> @date Ferbruay, 2015 - Initial Version
6341	!>
6342	!> @param[inout] id_mask domain mask (from grid_split_domain)
6343	!> @param[in] id_minsize minimum size of sea to be kept
6344	!-------------------------------------------------------------------
6345
6346	IMPLICIT NONE
6347
6348	! Argument
6349	INTEGER(i4), DIMENSION(:,:), INTENT(INOUT) :: id_mask
6350	INTEGER(i4), INTENT(IN ) :: id_minsize
6351
6352	! local variable
6353	INTEGER(i4) :: il_ndom
6354	INTEGER(i4) :: il_minsize
6355	INTEGER(i4) :: il_msk
6356
6357	INTEGER(i4) :: jim
6358	INTEGER(i4) :: jjm
6359	INTEGER(i4) :: jip
6360	INTEGER(i4) :: jjp
6361
6362	INTEGER(i4), DIMENSION(2) :: il_shape
6363	INTEGER(i4), DIMENSION(:,:), ALLOCATABLE :: il_tmp
6364
6365	! loop indices
6366	INTEGER(i4) :: ii
6367	INTEGER(i4) :: ij
6368
6369	INTEGER(i4) :: ji
6370	INTEGER(i4) :: jj
6371	!----------------------------------------------------------------
6372
6373	il_shape(:)=SHAPE(id_mask(:,:))
6374	il_ndom=MINVAL(id_mask(:,:))
6375
6376	ALLOCATE( il_tmp(il_shape(1),il_shape(2)) )
6377	il_tmp(:,:)=0
6378	DO ji=-1,il_ndom,-1
6379	WHERE( id_mask(:,:)==ji )
6380	il_tmp(:,:)=SUM(id_mask(:,:),id_mask(:,:)==ji)/ji
6381	END WHERE
6382	ENDDO
6383
6384	DO WHILE( id_minsize > MINVAL(il_tmp(:,:)) )
6385
6386	DO jj=1,il_shape(2)
6387	DO ji=1,il_shape(1)
6388
6389	IF( il_tmp(ji,jj) < id_minsize )THEN
6390	jim=MAX(1,ji-1) ; jip=MIN(il_shape(1),ji+1)
6391	jjm=MAX(1,jj-1) ; jjp=MIN(il_shape(2),jj+1)
6392
6393	il_msk=0
6394	DO ij=jjm,jjp
6395	DO ii=jim,jip
6396	IF( id_mask(ii,ij) /= id_mask(ji,jj) )THEN
6397	IF( il_msk == 0 )THEN
6398	il_msk=id_mask(ii,ij)
6399	ELSEIF( il_msk /= id_mask(ii,ij) )THEN
6400	CALL logger_error("GRID FILL SMALL MSK:"//&
6401	& " small domain not embedded in bigger one."//&
6402	& " point should be between two different"//&
6403	& " domain.")
6404	ENDIF
6405	ENDIF
6406	ENDDO
6407	ENDDO
6408	IF( il_msk /= 0 ) id_mask(ji,jj)=il_msk
6409
6410	ENDIF
6411
6412	ENDDO
6413	ENDDO
6414
6415
6416	il_tmp(:,:)=0
6417	DO ji=-1,il_ndom,-1
6418	WHERE( id_mask(:,:)==ji )
6419	il_tmp(:,:)=SUM(id_mask(:,:),id_mask(:,:)==ji)/ji
6420	END WHERE
6421	ENDDO
6422
6423	ENDDO
6424
6425	DEALLOCATE( il_tmp )
6426
6427	END SUBROUTINE grid_fill_small_msk
6428	!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
6429	END MODULE grid
6430

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

Download in other formats: