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.

create_boundary.F90 in branches/2016/dev_INGV_UKMO_2016/NEMOGCM/TOOLS/SIREN/src – NEMO

Context Navigation

source: branches/2016/dev_INGV_UKMO_2016/NEMOGCM/TOOLS/SIREN/src/create_boundary.F90 @ 7351

Last change on this file since 7351 was 7351, checked in by emanuelaclementi, 7 years ago
ticket #1805 step 3: /2016/dev_INGV_UKMO_2016 aligned to the trunk at revision 7161
File size: 62.2 KB

Line
1	!----------------------------------------------------------------------
2	! NEMO system team, System and Interface for oceanic RElocable Nesting
3	!----------------------------------------------------------------------
4	!
5	!
6	! PROGRAM: create_boundary
7	!
8	! DESCRIPTION:
9	!> @file
10	!> @brief
11	!> This program creates boundary files.
12	!>
13	!> @details
14	!> @section sec1 method
15	!> Variables are read from coarse grid standard output,
16	!> extracted or interpolated on fine grid.
17	!> Variables could also be manually written.<br/>
18	!> @note
19	!> method could be different for each variable.
20	!>
21	!> @section sec2 how to
22	!> to create boundaries files:<br/>
23	!> @code{.sh}
24	!> ./SIREN/bin/create_boundary create_boundary.nam
25	!> @endcode
26	!> <br/>
27	!> \image html boundary_NEATL36_70.png
28	!> \image latex boundary_NEATL36_70.png
29	!>
30	!> @note
31	!> you could find a template of the namelist in templates directory.
32	!>
33	!> create_boundary.nam contains 9 namelists:<br/>
34	!> - logger namelist (namlog)
35	!> - config namelist (namcfg)
36	!> - coarse grid namelist (namcrs)
37	!> - fine grid namelist (namfin)
38	!> - variable namelist (namvar)
39	!> - nesting namelist (namnst)
40	!> - boundary namelist (nambdy)
41	!> - vertical grid namelist (namzgr)
42	!> - output namelist (namout)
43	!>
44	!> * _logger namelist (namlog)_:<br/>
45	!> - cn_logfile : log filename
46	!> - cn_verbosity : verbosity ('trace','debug','info',
47	!> 'warning','error','fatal','none')
48	!> - in_maxerror : maximum number of error allowed
49	!>
50	!> * _config namelist (namcfg)_:<br/>
51	!> - cn_varcfg : variable configuration file
52	!> (see ./SIREN/cfg/variable.cfg)
53	!> - cn_dumcfg : useless (dummy) configuration file, for useless
54	!> dimension or variable (see ./SIREN/cfg/dummy.cfg).
55	!>
56	!> * _coarse grid namelist (namcrs)_:<br/>
57	!> - cn_coord0 : coordinate file
58	!> - in_perio0 : NEMO periodicity index (see Model Boundary Condition in
59	!> [NEMO documentation](http://www.nemo-ocean.eu/About-NEMO/Reference-manuals))
60	!>
61	!> * _fine grid namelist (namfin)_:<br/>
62	!> - cn_coord1 : coordinate file
63	!> - cn_bathy1 : bathymetry file
64	!> - in_perio1 : periodicity index
65	!>
66	!> * _vertical grid namelist (namzgr)_:<br/>
67	!> - dn_pp_to_be_computed :
68	!> - dn_ppsur :
69	!> - dn_ppa0 :
70	!> - dn_ppa1 :
71	!> - dn_ppa2 :
72	!> - dn_ppkth :
73	!> - dn_ppkth2 :
74	!> - dn_ppacr :
75	!> - dn_ppacr2 :
76	!> - dn_ppdzmin :
77	!> - dn_pphmax :
78	!> - in_nlevel : number of vertical level
79	!>
80	!> * _partial step namelist (namzps)_:<br/>
81	!> - dn_e3zps_min :
82	!> - dn_e3zps_rat :
83	!>
84	!> * _variable namelist (namvar)_:<br/>
85	!> - cn_varfile : list of variable, and associated file<br/>
86	!> cn_varfile is the path and filename of the file where find
87	!> variable.<br/>
88	!> @note
89	!> cn_varfile could be a matrix of value, if you want to filled
90	!> manually variable value.<br/>
91	!> the variable array of value is split into equal subdomain.<br/>
92	!> Each subdomain is filled with the corresponding value
93	!> of the matrix.<br/>
94	!> separators used to defined matrix are:
95	!> - ',' for line
96	!> - '/' for row
97	!> - '\' for level<br/>
98	!> Example:<br/>
99	!> 3,2,3/1,4,5 => @f$ \left( \begin{array}{ccc}
100	!> 3 & 2 & 3 \\\\
101	!> 1 & 4 & 5 \end{array} \right) @f$
102	!> @warning
103	!> the same matrix is used for all boundaries.
104	!>
105	!> Examples:
106	!> - 'votemper:gridT.nc', 'vozocrtx:gridU.nc'
107	!> - 'votemper:10\25', 'vozocrtx:gridU.nc'
108	!>
109	!> - cn_varinfo : list of variable and extra information about request(s)
110	!> to be used (separated by ',').<br/>
111	!> each elements of cn_varinfo is a string character.<br/>
112	!> it is composed of the variable name follow by ':',
113	!> then request(s) to be used on this variable.<br/>
114	!> request could be:
115	!> - int = interpolation method
116	!> - ext = extrapolation method
117	!> - flt = filter method
118	!> - min = minimum value
119	!> - max = maximum value
120	!> - unt = new units
121	!> - unf = unit scale factor (linked to new units)
122	!>
123	!> requests must be separated by ';'.<br/>
124	!> order of requests does not matter.
125	!>
126	!> informations about available method could be find in @ref interp,
127	!> @ref extrap and @ref filter.<br/>
128	!>
129	!> Example: 'votemper:int=linear;flt=hann;ext=dist_weight',
130	!> 'vosaline:int=cubic'
131	!> @note
132	!> If you do not specify a method which is required,
133	!> default one is apply.
134	!>
135	!> * _nesting namelist (namnst)_:<br/>
136	!> - in_rhoi : refinement factor in i-direction
137	!> - in_rhoj : refinement factor in j-direction
138	!>
139	!> * _boundary namelist (nambdy)_:<br/>
140	!> - ln_north : use north boundary
141	!> - ln_south : use south boundary
142	!> - ln_east : use east boundary
143	!> - ln_west : use west boundary
144	!> - cn_north : north boundary indices on fine grid
145	!> cn_north is a string character defining boundary
146	!> segmentation.<br/>
147	!> segments are separated by '\|'.<br/>
148	!> each segments of the boundary is composed of:
149	!> - indice of velocity (orthogonal to boundary .ie.
150	!> for north boundary, J-indice).
151	!> - indice of segment start (I-indice for north boundary)
152	!> - indice of segment end (I-indice for north boundary)<br/>
153	!> indices must be separated by ':' .<br/>
154	!> - optionally, boundary size could be added between '(' and ')'
155	!> in the definition of the first segment.
156	!> @note
157	!> boundary width is the same for all segments of one boundary.
158	!>
159	!> Examples:
160	!> - cn_north='index1,first1:last1(width)'
161	!> - cn_north='index1(width),first1:last1\|index2,first2:last2'
162	!> \image html boundary_50.png
163	!> \image latex boundary_50.png
164	!> - cn_south : south boundary indices on fine grid
165	!> - cn_east : east boundary indices on fine grid
166	!> - cn_west : west boundary indices on fine grid
167	!> - ln_oneseg : force to use only one segment for each boundary or not
168	!>
169	!> * _output namelist (namout)_:<br/>
170	!> - cn_fileout : fine grid boundary basename
171	!> (cardinal point and segment number will be automatically added)
172	!> - dn_dayofs : date offset in day (change only ouput file name)
173	!> - ln_extrap : extrapolate land point or not
174	!>
175	!> Examples:
176	!> - cn_fileout='boundary.nc'<br/>
177	!> if time_counter (16/07/2015 00h) is read on input file (see varfile),
178	!> west boundary will be named boundary_west_y2015m07d16
179	!> - dn_dayofs=-2.<br/>
180	!> if you use day offset you get boundary_west_y2015m07d14
181	!>
182	!>
183	!> @author J.Paul
184	! REVISION HISTORY:
185	!> @date November, 2013 - Initial Version
186	!> @date September, 2014
187	!> - add header for user
188	!> - take into account grid point to compue boundaries
189	!> - reorder output dimension for north and south boundaries
190	!> @date June, 2015
191	!> - extrapolate all land points, and add ln_extrap in namelist.
192	!> - allow to change unit.
193	!> @date July, 2015
194	!> - add namelist parameter to shift date of output file name.
195	!> @date September, 2015
196	!> - manage useless (dummy) variable, attributes, and dimension
197	!> - allow to run on multi processors with key_mpp_mpi
198	!> @date January, 2016
199	!> - same process use for variable extracted or interpolated from input file.
200	!>
201	!> @note Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
202	!----------------------------------------------------------------------
203	PROGRAM create_boundary
204
205	USE netcdf ! nf90 library
206	USE global ! global variable
207	USE phycst ! physical constant
208	USE kind ! F90 kind parameter
209	USE fct ! basic useful function
210	USE date ! date manager
211	USE att ! attribute manager
212	USE dim ! dimension manager
213	USE var ! variable manager
214	USE file ! file manager
215	USE multi ! multi file manager
216	USE boundary ! boundary manager
217	USE iom ! I/O manager
218	USE dom ! domain manager
219	USE grid ! grid manager
220	USE vgrid ! vertical grid manager
221	USE extrap ! extrapolation manager
222	USE interp ! interpolation manager
223	USE filter ! filter manager
224	USE mpp ! MPP manager
225	USE iom_mpp ! MPP I/O manager
226
227	IMPLICIT NONE
228
229	! local variable
230	INTEGER(i4) :: il_narg
231
232	#if defined key_mpp_mpi
233	! mpp variable
234	CHARACTER(LEN=lc), DIMENSION(:) , ALLOCATABLE :: cl_namelist
235	INTEGER(i4) :: ierror
236	INTEGER(i4) :: iproc
237	INTEGER(i4) :: nproc
238	INTEGER(i4) , DIMENSION(:) , ALLOCATABLE :: il_nprog
239
240	! loop indices
241	INTEGER(i4) :: jm
242	#else
243	CHARACTER(LEN=lc) :: cl_namelist
244	#endif
245	!-------------------------------------------------------------------
246	#if defined key_mpp_mpi
247	INCLUDE 'mpif.h'
248	#endif
249	!-------------------------------------------------------------------
250
251	il_narg=COMMAND_ARGUMENT_COUNT() !f03 intrinsec
252	#if ! defined key_mpp_mpi
253
254	IF( il_narg/=1 )THEN
255	PRINT *,"CREATE BOUNDARY: ERROR. need one namelist"
256	STOP
257	ELSE
258	CALL GET_COMMAND_ARGUMENT(1,cl_namelist) !f03 intrinsec
259	ENDIF
260
261	CALL create__boundary(cl_namelist)
262
263	#else
264
265	! Initialize MPI
266	CALL mpi_init(ierror)
267	CALL mpi_comm_rank(mpi_comm_world,iproc,ierror)
268	CALL mpi_comm_size(mpi_comm_world,nproc,ierror)
269
270	IF( il_narg==0 )THEN
271	PRINT *,"CREATE BOUNDARY: ERROR. need at least one namelist"
272	STOP
273	ELSE
274	ALLOCATE(cl_namelist(il_narg))
275	DO jm=1,il_narg
276	CALL GET_COMMAND_ARGUMENT(jm,cl_namelist(jm))
277	ENDDO
278	ENDIF
279
280	ALLOCATE(il_nprog(il_narg))
281	DO jm=1, il_narg
282	il_nprog(jm)= MOD(jm,nproc)
283	ENDDO
284
285	DO jm=1, il_narg
286	IF ( il_nprog(jm) .eq. iproc ) THEN
287	CALL create__boundary(cl_namelist(jm))
288	ENDIF
289	ENDDO
290
291	CALL mpi_finalize(ierror)
292
293	DEALLOCATE(cl_namelist)
294	DEALLOCATE(il_nprog)
295	#endif
296
297	CONTAINS
298	SUBROUTINE create__boundary(cd_namelist)
299	!-------------------------------------------------------------------
300	!> @brief
301	!> This subroutine create boundary files.
302	!>
303	!> @details
304	!>
305	!> @author J.Paul
306	!> @date January, 2016 - Initial Version
307	!>
308	!> @param[in] cd_namelist namelist file
309	!-------------------------------------------------------------------
310
311	USE logger ! log file manager
312
313	IMPLICIT NONE
314	! Argument
315	CHARACTER(LEN=lc), INTENT(IN) :: cd_namelist
316
317	! local variable
318	CHARACTER(LEN=lc) :: cl_date
319	CHARACTER(LEN=lc) :: cl_name
320	CHARACTER(LEN=lc) :: cl_bdyout
321	CHARACTER(LEN=lc) :: cl_data
322	CHARACTER(LEN=lc) :: cl_dimorder
323	CHARACTER(LEN=lc) :: cl_fmt
324
325	INTEGER(i4) :: il_status
326	INTEGER(i4) :: il_fileid
327	INTEGER(i4) :: il_imin0
328	INTEGER(i4) :: il_imax0
329	INTEGER(i4) :: il_jmin0
330	INTEGER(i4) :: il_jmax0
331	INTEGER(i4) :: il_shift
332	INTEGER(i4) , DIMENSION(ip_maxdim) :: il_rho
333	INTEGER(i4) , DIMENSION(2,2) :: il_offset
334	INTEGER(i4) , DIMENSION(2,2) :: il_ind
335
336	LOGICAL :: ll_exist
337
338	TYPE(TATT) :: tl_att
339
340	TYPE(TVAR) :: tl_depth
341	TYPE(TVAR) :: tl_time
342	TYPE(TVAR) :: tl_var1
343	TYPE(TVAR) :: tl_var0
344	TYPE(TVAR) :: tl_lon1
345	TYPE(TVAR) :: tl_lat1
346	TYPE(TVAR) :: tl_lvl1
347	TYPE(TVAR) , DIMENSION(:) , ALLOCATABLE :: tl_level
348	TYPE(TVAR) , DIMENSION(:,:,:) , ALLOCATABLE :: tl_seglvl1
349	TYPE(TVAR) , DIMENSION(:,:,:) , ALLOCATABLE :: tl_segvar1
350
351	TYPE(TDIM) , DIMENSION(ip_maxdim) :: tl_dim
352
353	TYPE(TDATE) :: tl_date
354
355	TYPE(TBDY) , DIMENSION(ip_ncard) :: tl_bdy
356
357	TYPE(TDOM) :: tl_dom0
358	TYPE(TDOM) :: tl_dom1
359	TYPE(TDOM) , DIMENSION(:,:,:) , ALLOCATABLE :: tl_segdom1
360
361	TYPE(TFILE) :: tl_fileout
362
363	TYPE(TMPP) :: tl_coord0
364	TYPE(TMPP) :: tl_coord1
365	TYPE(TMPP) :: tl_bathy1
366	TYPE(TMPP) :: tl_mpp
367
368	TYPE(TMULTI) :: tl_multi
369
370	! loop indices
371	INTEGER(i4) :: jvar
372	INTEGER(i4) :: jpoint
373	INTEGER(i4) :: ji
374	INTEGER(i4) :: jj
375	INTEGER(i4) :: jk
376	INTEGER(i4) :: jl
377
378	! namelist variable
379	! namlog
380	CHARACTER(LEN=lc) :: cn_logfile = 'create_boundary.log'
381	CHARACTER(LEN=lc) :: cn_verbosity = 'warning'
382	INTEGER(i4) :: in_maxerror = 5
383
384	! namcfg
385	CHARACTER(LEN=lc) :: cn_varcfg = 'variable.cfg'
386	CHARACTER(LEN=lc) :: cn_dumcfg = 'dummy.cfg'
387
388	! namcrs
389	CHARACTER(LEN=lc) :: cn_coord0 = ''
390	INTEGER(i4) :: in_perio0 = -1
391
392	! namfin
393	CHARACTER(LEN=lc) :: cn_coord1 = ''
394	CHARACTER(LEN=lc) :: cn_bathy1 = ''
395	INTEGER(i4) :: in_perio1 = -1
396
397	!namzgr
398	REAL(dp) :: dn_pp_to_be_computed = 0._dp
399	REAL(dp) :: dn_ppsur = -3958.951371276829_dp
400	REAL(dp) :: dn_ppa0 = 103.953009600000_dp
401	REAL(dp) :: dn_ppa1 = 2.415951269000_dp
402	REAL(dp) :: dn_ppa2 = 100.760928500000_dp
403	REAL(dp) :: dn_ppkth = 15.351013700000_dp
404	REAL(dp) :: dn_ppkth2 = 48.029893720000_dp
405	REAL(dp) :: dn_ppacr = 7.000000000000_dp
406	REAL(dp) :: dn_ppacr2 = 13.000000000000_dp
407	REAL(dp) :: dn_ppdzmin = 6._dp
408	REAL(dp) :: dn_pphmax = 5750._dp
409	INTEGER(i4) :: in_nlevel = 75
410
411	!namzps
412	REAL(dp) :: dn_e3zps_min = 25._dp
413	REAL(dp) :: dn_e3zps_rat = 0.2_dp
414
415	! namvar
416	CHARACTER(LEN=lc), DIMENSION(ip_maxvar) :: cn_varfile = ''
417	CHARACTER(LEN=lc), DIMENSION(ip_maxvar) :: cn_varinfo = ''
418
419	! namnst
420	INTEGER(i4) :: in_rhoi = 0
421	INTEGER(i4) :: in_rhoj = 0
422
423	! nambdy
424	LOGICAL :: ln_north = .TRUE.
425	LOGICAL :: ln_south = .TRUE.
426	LOGICAL :: ln_east = .TRUE.
427	LOGICAL :: ln_west = .TRUE.
428	LOGICAL :: ln_oneseg = .TRUE.
429	CHARACTER(LEN=lc) :: cn_north = ''
430	CHARACTER(LEN=lc) :: cn_south = ''
431	CHARACTER(LEN=lc) :: cn_east = ''
432	CHARACTER(LEN=lc) :: cn_west = ''
433
434	! namout
435	CHARACTER(LEN=lc) :: cn_fileout = 'boundary.nc'
436	REAL(dp) :: dn_dayofs = 0._dp
437	LOGICAL :: ln_extrap = .FALSE.
438	!-------------------------------------------------------------------
439
440	NAMELIST /namlog/ & !< logger namelist
441	& cn_logfile, & !< log file
442	& cn_verbosity, & !< log verbosity
443	& in_maxerror
444
445	NAMELIST /namcfg/ & !< config namelist
446	& cn_varcfg, & !< variable configuration file
447	& cn_dumcfg !< dummy configuration file
448
449	NAMELIST /namcrs/ & !< coarse grid namelist
450	& cn_coord0, & !< coordinate file
451	& in_perio0 !< periodicity index
452
453	NAMELIST /namfin/ & !< fine grid namelist
454	& cn_coord1, & !< coordinate file
455	& cn_bathy1, & !< bathymetry file
456	& in_perio1 !< periodicity index
457
458	NAMELIST /namzgr/ &
459	& dn_pp_to_be_computed, &
460	& dn_ppsur, &
461	& dn_ppa0, &
462	& dn_ppa1, &
463	& dn_ppa2, &
464	& dn_ppkth, &
465	& dn_ppkth2, &
466	& dn_ppacr, &
467	& dn_ppacr2, &
468	& dn_ppdzmin, &
469	& dn_pphmax, &
470	& in_nlevel !< number of vertical level
471
472	NAMELIST /namzps/ &
473	& dn_e3zps_min, &
474	& dn_e3zps_rat
475
476	NAMELIST /namvar/ & !< variable namelist
477	& cn_varfile, & !< list of variable and file where find it. (ex: 'votemper:GLORYS_gridT.nc' )
478	& cn_varinfo !< list of variable and method to apply on. (ex: 'votemper:linear','vosaline:cubic' )
479
480	NAMELIST /namnst/ & !< nesting namelist
481	& in_rhoi, & !< refinement factor in i-direction
482	& in_rhoj !< refinement factor in j-direction
483
484	NAMELIST /nambdy/ & !< boundary namelist
485	& ln_north, & !< use north boundary
486	& ln_south, & !< use south boundary
487	& ln_east , & !< use east boundary
488	& ln_west , & !< use west boundary
489	& cn_north, & !< north boundary indices on fine grid
490	& cn_south, & !< south boundary indices on fine grid
491	& cn_east , & !< east boundary indices on fine grid
492	& cn_west , & !< west boundary indices on fine grid
493	& ln_oneseg !< use only one segment for each boundary or not
494
495	NAMELIST /namout/ & !< output namelist
496	& cn_fileout, & !< fine grid boundary file basename
497	& dn_dayofs, & !< date offset in day (change only ouput file name)
498	& ln_extrap !< extrapolate or not
499	!-------------------------------------------------------------------
500
501	! read namelist
502	INQUIRE(FILE=TRIM(cd_namelist), EXIST=ll_exist)
503
504	IF( ll_exist )THEN
505
506	il_fileid=fct_getunit()
507
508	OPEN( il_fileid, FILE=TRIM(cd_namelist), &
509	& FORM='FORMATTED', &
510	& ACCESS='SEQUENTIAL', &
511	& STATUS='OLD', &
512	& ACTION='READ', &
513	& IOSTAT=il_status)
514	CALL fct_err(il_status)
515	IF( il_status /= 0 )THEN
516	PRINT *,"CREATE BOUNDARY: ERROR opening "//TRIM(cd_namelist)
517	STOP
518	ENDIF
519
520	READ( il_fileid, NML = namlog )
521	! define log file
522	CALL logger_open(TRIM(cn_logfile),TRIM(cn_verbosity),in_maxerror)
523	CALL logger_header()
524
525	READ( il_fileid, NML = namcfg )
526	! get variable extra information
527	CALL var_def_extra(TRIM(cn_varcfg))
528
529	! get dummy variable
530	CALL var_get_dummy(TRIM(cn_dumcfg))
531	! get dummy dimension
532	CALL dim_get_dummy(TRIM(cn_dumcfg))
533	! get dummy attribute
534	CALL att_get_dummy(TRIM(cn_dumcfg))
535
536	READ( il_fileid, NML = namcrs )
537	READ( il_fileid, NML = namfin )
538	READ( il_fileid, NML = namzgr )
539	READ( il_fileid, NML = namvar )
540	! add user change in extra information
541	CALL var_chg_extra(cn_varinfo)
542	! match variable with file
543	tl_multi=multi_init(cn_varfile)
544
545	READ( il_fileid, NML = namnst )
546	READ( il_fileid, NML = nambdy )
547	READ( il_fileid, NML = namout )
548
549	CLOSE( il_fileid, IOSTAT=il_status )
550	CALL fct_err(il_status)
551	IF( il_status /= 0 )THEN
552	CALL logger_error("CREATE BOUNDARY: ERROR closing "//TRIM(cd_namelist))
553	ENDIF
554
555	ELSE
556
557	PRINT *,"CREATE BOUNDARY: ERROR. can not find "//TRIM(cd_namelist)
558	STOP
559
560	ENDIF
561
562	CALL multi_print(tl_multi)
563	IF( tl_multi%i_nvar <= 0 )THEN
564	CALL logger_fatal("CREATE BOUNDARY: no variable to be used."//&
565	& " check namelist.")
566	ENDIF
567
568	! open files
569	IF( TRIM(cn_coord0) /= '' )THEN
570	tl_coord0=mpp_init( file_init(TRIM(cn_coord0)), id_perio=in_perio0)
571	CALL grid_get_info(tl_coord0)
572	ELSE
573	CALL logger_fatal("CREATE BOUNDARY: can not find coarse grid "//&
574	& "coordinate file. check namelist")
575	ENDIF
576
577	IF( TRIM(cn_coord1) /= '' )THEN
578	tl_coord1=mpp_init( file_init(TRIM(cn_coord1)), id_perio=in_perio1)
579	CALL grid_get_info(tl_coord1)
580	ELSE
581	CALL logger_fatal("CREATE BOUNDARY: can not find fine grid coordinate "//&
582	& "file. check namelist")
583	ENDIF
584
585	IF( TRIM(cn_bathy1) /= '' )THEN
586	tl_bathy1=mpp_init( file_init(TRIM(cn_bathy1)), id_perio=in_perio1)
587	CALL grid_get_info(tl_bathy1)
588	ELSE
589	CALL logger_fatal("CREATE BOUNDARY: can not find fine grid bathymetry "//&
590	& "file. check namelist")
591	ENDIF
592
593	! check
594	! check output file do not already exist
595	! WARNING: do not work when use time to create output file name
596	DO jk=1,ip_ncard
597	cl_bdyout=boundary_set_filename( TRIM(cn_fileout), &
598	& TRIM(cp_card(jk)), 1 )
599	INQUIRE(FILE=TRIM(cl_bdyout), EXIST=ll_exist)
600	IF( ll_exist )THEN
601	CALL logger_fatal("CREATE BOUNDARY: output file "//TRIM(cl_bdyout)//&
602	& " already exist.")
603	ENDIF
604
605	cl_bdyout=boundary_set_filename( TRIM(cn_fileout), &
606	& TRIM(cp_card(jk)) )
607	INQUIRE(FILE=TRIM(cl_bdyout), EXIST=ll_exist)
608	IF( ll_exist )THEN
609	CALL logger_fatal("CREATE BOUNDARY: output file "//TRIM(cl_bdyout)//&
610	& " already exist.")
611	ENDIF
612	ENDDO
613
614	! check namelist
615	! check refinement factor
616	il_rho(:)=1
617	IF( in_rhoi < 1 .OR. in_rhoj < 1 )THEN
618	CALL logger_error("CREATE BOUNDARY: invalid refinement factor."//&
619	& " check namelist "//TRIM(cd_namelist))
620	ELSE
621	il_rho(jp_I)=in_rhoi
622	il_rho(jp_J)=in_rhoj
623	ENDIF
624
625	!
626	! compute coarse grid indices around fine grid
627	il_ind(:,:)=grid_get_coarse_index(tl_coord0, tl_coord1, &
628	& id_rho=il_rho(:))
629
630	il_imin0=il_ind(1,1) ; il_imax0=il_ind(1,2)
631	il_jmin0=il_ind(2,1) ; il_jmax0=il_ind(2,2)
632
633	! check domain validity
634	CALL grid_check_dom(tl_coord0, il_imin0, il_imax0, il_jmin0, il_jmax0)
635
636	! check coordinate file
637	CALL grid_check_coincidence( tl_coord0, tl_coord1, &
638	& il_imin0, il_imax0, &
639	& il_jmin0, il_jmax0, &
640	& il_rho(:) )
641
642	! read or compute boundary
643	CALL mpp_get_contour(tl_bathy1)
644
645	CALL iom_mpp_open(tl_bathy1)
646
647	tl_var1=iom_mpp_read_var(tl_bathy1,'Bathymetry')
648
649	CALL iom_mpp_close(tl_bathy1)
650
651	! get boundaries indices
652	tl_bdy(:)=boundary_init(tl_var1, ln_north, ln_south, ln_east, ln_west, &
653	& cn_north, cn_south, cn_east, cn_west, &
654	& ln_oneseg )
655
656
657	CALL var_clean(tl_var1)
658
659	! compute level
660	ALLOCATE(tl_level(ip_npoint))
661	tl_level(:)=vgrid_get_level(tl_bathy1, cd_namelist )
662
663	! get coordinate for each segment of each boundary
664	ALLOCATE( tl_segdom1(ip_npoint,ip_maxseg,ip_ncard) )
665	ALLOCATE( tl_seglvl1(ip_npoint,ip_maxseg,ip_ncard) )
666
667	DO jl=1,ip_ncard
668	IF( tl_bdy(jl)%l_use )THEN
669	DO jk=1,tl_bdy(jl)%i_nseg
670
671	! get fine grid segment domain
672	tl_segdom1(:,jk,jl)=create_boundary_get_dom( tl_bathy1, &
673	& tl_bdy(jl), jk )
674
675	IF( .NOT. ln_extrap )THEN
676	! get fine grid level
677	tl_seglvl1(:,jk,jl)= &
678	& create_boundary_get_level( tl_level(:), &
679	& tl_segdom1(:,jk,jl))
680	ENDIF
681
682	! add extra band to fine grid domain (if possible)
683	! to avoid dimension of one and so be able to compute offset
684	DO jj=1,ip_npoint
685	CALL dom_add_extra(tl_segdom1(jj,jk,jl), &
686	& il_rho(jp_I), il_rho(jp_J))
687	ENDDO
688
689	ENDDO
690	ENDIF
691	ENDDO
692
693	! clean
694	CALL var_clean(tl_level(:))
695	DEALLOCATE(tl_level)
696
697	! clean bathy
698	CALL mpp_clean(tl_bathy1)
699
700	ALLOCATE( tl_segvar1(tl_multi%i_nvar,ip_maxseg,ip_ncard) )
701	! compute boundary for variable to be used (see namelist)
702	IF( .NOT. ASSOCIATED(tl_multi%t_mpp) )THEN
703	CALL logger_error("CREATE BOUNDARY: no file to work on. "//&
704	& "check cn_varfile in namelist.")
705	ELSE
706
707	jvar=0
708	! for each file
709	DO ji=1,tl_multi%i_nmpp
710
711	WRITE(cl_data,'(a,i2.2)') 'data-',jvar+1
712
713	IF( .NOT. ASSOCIATED(tl_multi%t_mpp(ji)%t_proc(1)%t_var) )THEN
714
715	CALL logger_error("CREATE BOUNDARY: no variable to work on for "//&
716	& "mpp "//TRIM(tl_multi%t_mpp(ji)%c_name)//&
717	& ". check cn_varfile in namelist.")
718
719	ELSEIF( TRIM(tl_multi%t_mpp(ji)%c_name) == TRIM(cl_data) )THEN
720	!- use input matrix to fill variable
721
722	WRITE(*,'(a)') "work on data"
723	! for each variable initialise from matrix
724	DO jj=1,tl_multi%t_mpp(ji)%t_proc(1)%i_nvar
725
726	jvar=jvar+1
727	WRITE(*,'(2x,a,a)') "work on variable "//&
728	& TRIM(tl_multi%t_mpp(ji)%t_proc(1)%t_var(jj)%c_name)
729
730	tl_var1=var_copy(tl_multi%t_mpp(ji)%t_proc(1)%t_var(jj))
731
732	SELECT CASE(TRIM(tl_var1%c_point))
733	CASE DEFAULT !'T'
734	jpoint=jp_T
735	CASE('U')
736	jpoint=jp_U
737	CASE('V')
738	jpoint=jp_V
739	CASE('F')
740	jpoint=jp_F
741	END SELECT
742
743	WRITE(*,'(4x,a,a)') 'work on '//TRIM(tl_var1%c_name)
744	DO jl=1,ip_ncard
745	IF( tl_bdy(jl)%l_use )THEN
746
747	DO jk=1,tl_bdy(jl)%i_nseg
748
749	! fill value with matrix data
750	tl_segvar1(jvar,jk,jl)=create_boundary_matrix( &
751	& tl_var1, &
752	& tl_segdom1(jpoint,jk,jl), &
753	& in_nlevel )
754
755	!del extra
756	CALL dom_del_extra( tl_segvar1(jvar,jk,jl), &
757	& tl_segdom1(jpoint,jk,jl) )
758
759	ENDDO
760
761	ENDIF
762	ENDDO
763
764	! clean
765	CALL var_clean(tl_var1)
766
767	ENDDO
768
769	!- end of use input matrix to fill variable
770	ELSE
771	!- use mpp file to fill variable
772
773	WRITE(*,'(a)') "work on file "//TRIM(tl_multi%t_mpp(ji)%c_name)
774	!
775	tl_mpp=mpp_init(file_init(TRIM(tl_multi%t_mpp(ji)%t_proc(1)%c_name)))
776	CALL grid_get_info(tl_mpp)
777
778	DO jl=1,ip_ncard
779	IF( tl_bdy(jl)%l_use )THEN
780
781	WRITE(*,'(2x,a,a)') 'work on '//TRIM(tl_bdy(jl)%c_card)//&
782	& ' boundary'
783	DO jk=1,tl_bdy(jl)%i_nseg
784
785	! for each variable of this file
786	DO jj=1,tl_multi%t_mpp(ji)%t_proc(1)%i_nvar
787
788	WRITE(*,'(4x,a,a)') "work on variable "//&
789	& TRIM(tl_multi%t_mpp(ji)%t_proc(1)%t_var(jj)%c_name)
790
791	tl_var0=var_copy(tl_multi%t_mpp(ji)%t_proc(1)%t_var(jj))
792
793	! open mpp file
794	CALL iom_mpp_open(tl_mpp)
795
796	! get or check depth value
797	CALL create_boundary_check_depth( tl_var0, tl_mpp, &
798	& in_nlevel, tl_depth )
799
800	! get or check time value
801	CALL create_boundary_check_time( tl_var0, tl_mpp, &
802	& tl_time )
803
804	! close mpp file
805	CALL iom_mpp_close(tl_mpp)
806
807	! open mpp file on domain
808	SELECT CASE(TRIM(tl_var0%c_point))
809	CASE DEFAULT !'T'
810	jpoint=jp_T
811	CASE('U')
812	jpoint=jp_U
813	CASE('V')
814	jpoint=jp_V
815	CASE('F')
816	jpoint=jp_F
817	END SELECT
818
819	tl_dom1=dom_copy(tl_segdom1(jpoint,jk,jl))
820
821	CALL create_boundary_get_coord( tl_coord1, tl_dom1, &
822	& tl_var0%c_point, &
823	& tl_lon1, tl_lat1 )
824
825	! get coarse grid indices of this segment
826	il_ind(:,:)=grid_get_coarse_index(tl_coord0, &
827	& tl_lon1, tl_lat1, &
828	& id_rho=il_rho(:) )
829
830	IF( ANY(il_ind(:,:)==0) )THEN
831	CALL logger_error("CREATE BOUNDARY: error "//&
832	& "computing coarse grid indices")
833	ELSE
834	il_imin0=il_ind(1,1)
835	il_imax0=il_ind(1,2)
836
837	il_jmin0=il_ind(2,1)
838	il_jmax0=il_ind(2,2)
839	ENDIF
840
841	il_offset(:,:)= grid_get_fine_offset( &
842	& tl_coord0, &
843	& il_imin0, il_jmin0,&
844	& il_imax0, il_jmax0,&
845	& tl_lon1%d_value(:,:,1,1),&
846	& tl_lat1%d_value(:,:,1,1),&
847	& il_rho(:),&
848	& TRIM(tl_var0%c_point) )
849
850	! compute coarse grid segment domain
851	tl_dom0=dom_init( tl_coord0, &
852	& il_imin0, il_imax0,&
853	& il_jmin0, il_jmax0 )
854
855	! add extra band (if possible) to compute interpolation
856	CALL dom_add_extra(tl_dom0)
857
858	! open mpp files
859	CALL iom_dom_open(tl_mpp, tl_dom0)
860
861	cl_name=tl_var0%c_name
862	! read variable value on domain
863	tl_segvar1(jvar+jj,jk,jl)= &
864	& iom_dom_read_var(tl_mpp, TRIM(cl_name), tl_dom0)
865
866	IF( ANY(il_rho(:)/=1) )THEN
867	WRITE(*,'(4x,a,a)') "interp variable "//TRIM(cl_name)
868	! work on variable
869	CALL create_boundary_interp( &
870	& tl_segvar1(jvar+jj,jk,jl),&
871	& il_rho(:), il_offset(:,:) )
872	ENDIF
873
874	! remove extraband added to domain
875	CALL dom_del_extra( tl_segvar1(jvar+jj,jk,jl), &
876	& tl_dom0, il_rho(:) )
877
878	! del extra point on fine grid
879	CALL dom_del_extra( tl_segvar1(jvar+jj,jk,jl), &
880	& tl_dom1 )
881	! clean extra point information on coarse grid domain
882	CALL dom_clean_extra( tl_dom0 )
883
884	! add attribute to variable
885	tl_att=att_init('src_file',&
886	& TRIM(fct_basename(tl_mpp%c_name)))
887	CALL var_move_att(tl_segvar1(jvar+jj,jk,jl), &
888	& tl_att)
889
890	!
891	tl_att=att_init('src_i_indices',&
892	& (/tl_dom0%i_imin, tl_dom0%i_imax/))
893	CALL var_move_att(tl_segvar1(jvar+jj,jk,jl), &
894	& tl_att)
895
896	tl_att=att_init('src_j_indices', &
897	& (/tl_dom0%i_jmin, tl_dom0%i_jmax/))
898	CALL var_move_att(tl_segvar1(jvar+jj,jk,jl), &
899	& tl_att)
900
901	IF( ANY(il_rho(:)/=1) )THEN
902	tl_att=att_init("refinment_factor", &
903	& (/il_rho(jp_I),il_rho(jp_J)/))
904	CALL var_move_att(tl_segvar1(jvar+jj,jk,jl), &
905	& tl_att)
906	ENDIF
907
908	! clean structure
909	CALL att_clean(tl_att)
910
911	! clean
912	CALL dom_clean(tl_dom0)
913	CALL dom_clean(tl_dom1)
914
915	! close mpp files
916	CALL iom_dom_close(tl_mpp)
917
918	! clean structure
919	CALL var_clean(tl_lon1)
920	CALL var_clean(tl_lat1)
921	CALL var_clean(tl_lvl1)
922
923	ENDDO ! jj
924
925	! clean
926	CALL var_clean(tl_var0)
927
928	ENDDO ! jk
929
930	ENDIF
931	ENDDO ! jl
932
933	jvar=jvar+tl_multi%t_mpp(ji)%t_proc(1)%i_nvar
934
935	! clean
936	CALL mpp_clean(tl_mpp)
937
938	!- end of use file to fill variable
939	ENDIF
940	ENDDO ! ji
941	ENDIF
942
943	IF( jvar /= tl_multi%i_nvar )THEN
944	CALL logger_error("CREATE BOUNDARY: it seems some variable "//&
945	& "can not be read")
946	ENDIF
947
948	! write file for each segment of each boundary
949	DO jl=1,ip_ncard
950	IF( tl_bdy(jl)%l_use )THEN
951
952	DO jk=1,tl_bdy(jl)%i_nseg
953	!-
954	CALL create_boundary_get_coord( tl_coord1, tl_segdom1(jp_T,jk,jl),&
955	& 'T', tl_lon1, tl_lat1 )
956
957	! force to use nav_lon, nav_lat as variable name
958	tl_lon1%c_name='nav_lon'
959	tl_lat1%c_name='nav_lat'
960
961	! del extra point on fine grid
962	CALL dom_del_extra( tl_lon1, tl_segdom1(jp_T,jk,jl) )
963	CALL dom_del_extra( tl_lat1, tl_segdom1(jp_T,jk,jl) )
964
965	! clean
966	DO jpoint=1,ip_npoint
967	CALL dom_clean(tl_segdom1(jpoint,jk,jl))
968	ENDDO
969
970	! swap array
971	CALL boundary_swap(tl_lon1, tl_bdy(jl))
972	CALL boundary_swap(tl_lat1, tl_bdy(jl))
973	DO jvar=1,tl_multi%i_nvar
974
975	! use additional request
976	! change unit and apply factor
977	CALL var_chg_unit(tl_segvar1(jvar,jk,jl))
978
979	! forced min and max value
980	CALL var_limit_value(tl_segvar1(jvar,jk,jl))
981
982	! filter
983	CALL filter_fill_value(tl_segvar1(jvar,jk,jl))
984
985	IF( .NOT. ln_extrap )THEN
986	! use mask
987	SELECT CASE(TRIM(tl_segvar1(jvar,jk,jl)%c_point))
988	CASE DEFAULT !'T'
989	jpoint=jp_T
990	CASE('U')
991	jpoint=jp_U
992	CASE('V')
993	jpoint=jp_V
994	CASE('F')
995	jpoint=jp_F
996	END SELECT
997
998	CALL create_boundary_use_mask(tl_segvar1(jvar,jk,jl), &
999	& tl_seglvl1(jpoint,jk,jl))
1000	ENDIF
1001
1002	! swap dimension order
1003	CALL boundary_swap(tl_segvar1(jvar,jk,jl), tl_bdy(jl))
1004
1005	ENDDO
1006
1007	! create file
1008	! create file structure
1009	! set file namearray of level variable structure
1010	IF( tl_bdy(jl)%i_nseg > 1 )THEN
1011	IF( ASSOCIATED(tl_time%d_value) )THEN
1012	cl_fmt="('y',i0.4,'m',i0.2,'d',i0.2)"
1013	tl_date=var_to_date(tl_time)
1014	tl_date=tl_date+dn_dayofs
1015	cl_date=date_print( tl_date, cl_fmt )
1016
1017	cl_bdyout=boundary_set_filename( TRIM(cn_fileout), &
1018	& TRIM(tl_bdy(jl)%c_card), jk,&
1019	& cd_date=TRIM(cl_date) )
1020	ELSE
1021	cl_bdyout=boundary_set_filename( TRIM(cn_fileout), &
1022	& TRIM(tl_bdy(jl)%c_card), jk )
1023	ENDIF
1024	ELSE
1025	IF( ASSOCIATED(tl_time%d_value) )THEN
1026	cl_fmt="('y',i0.4,'m',i0.2,'d',i0.2)"
1027	tl_date=var_to_date(tl_time)
1028	tl_date=tl_date+dn_dayofs
1029	cl_date=date_print( tl_date, cl_fmt )
1030
1031	cl_bdyout=boundary_set_filename( TRIM(cn_fileout), &
1032	& TRIM(tl_bdy(jl)%c_card), &
1033	& cd_date=TRIM(cl_date) )
1034	ELSE
1035	cl_bdyout=boundary_set_filename( TRIM(cn_fileout), &
1036	& TRIM(tl_bdy(jl)%c_card) )
1037	ENDIF
1038	ENDIF
1039	!
1040	tl_fileout=file_init(TRIM(cl_bdyout),id_perio=in_perio1)
1041
1042	! add dimension
1043	tl_dim(:)=var_max_dim(tl_segvar1(:,jk,jl))
1044
1045	SELECT CASE(TRIM(tl_bdy(jl)%c_card))
1046	CASE DEFAULT ! 'north','south'
1047	cl_dimorder='xyzt'
1048	CASE('east','west')
1049	cl_dimorder='yxzt'
1050	END SELECT
1051
1052	DO ji=1,ip_maxdim
1053	IF( tl_dim(ji)%l_use ) CALL file_add_dim(tl_fileout, tl_dim(ji))
1054	ENDDO
1055
1056	! add variables
1057	IF( ALL( tl_dim(1:2)%l_use ) )THEN
1058	! add longitude
1059	CALL file_add_var(tl_fileout, tl_lon1)
1060	CALL var_clean(tl_lon1)
1061
1062	! add latitude
1063	CALL file_add_var(tl_fileout, tl_lat1)
1064	CALL var_clean(tl_lat1)
1065	ENDIF
1066
1067
1068
1069	IF( tl_dim(3)%l_use )THEN
1070	IF( ASSOCIATED(tl_depth%d_value) )THEN
1071	! add depth
1072	CALL file_add_var(tl_fileout, tl_depth)
1073	ENDIF
1074	ENDIF
1075
1076	IF( tl_dim(4)%l_use )THEN
1077	IF( ASSOCIATED(tl_time%d_value) )THEN
1078	! add time
1079	CALL file_add_var(tl_fileout, tl_time)
1080	ENDIF
1081	ENDIF
1082
1083	! add other variable
1084	DO jvar=tl_multi%i_nvar,1,-1
1085	CALL file_add_var(tl_fileout, tl_segvar1(jvar,jk,jl))
1086	CALL var_clean(tl_segvar1(jvar,jk,jl))
1087	ENDDO
1088
1089	! add some attribute
1090	tl_att=att_init("Created_by","SIREN create_boundary")
1091	CALL file_add_att(tl_fileout, tl_att)
1092
1093	cl_date=date_print(date_now())
1094	tl_att=att_init("Creation_date",cl_date)
1095	CALL file_add_att(tl_fileout, tl_att)
1096
1097	! add shift on north and east boundary
1098	! boundary compute on T point but express on U or V point
1099	SELECT CASE(TRIM(tl_bdy(jl)%c_card))
1100	CASE DEFAULT ! 'south','west'
1101	il_shift=0
1102	CASE('north','east')
1103	il_shift=1
1104	END SELECT
1105
1106	! add indice of velocity row or column
1107	tl_att=att_init('bdy_ind',tl_bdy(jl)%t_seg(jk)%i_index-il_shift)
1108	CALL file_move_att(tl_fileout, tl_att)
1109
1110	! add width of the relaxation zone
1111	tl_att=att_init('bdy_width',tl_bdy(jl)%t_seg(jk)%i_width)
1112	CALL file_move_att(tl_fileout, tl_att)
1113
1114	! add indice of segment start
1115	tl_att=att_init('bdy_deb',tl_bdy(jl)%t_seg(jk)%i_first)
1116	CALL file_move_att(tl_fileout, tl_att)
1117
1118	! add indice of segment end
1119	tl_att=att_init('bdy_end',tl_bdy(jl)%t_seg(jk)%i_last)
1120	CALL file_move_att(tl_fileout, tl_att)
1121
1122	! clean
1123	CALL att_clean(tl_att)
1124
1125	! create file
1126	CALL iom_create(tl_fileout)
1127
1128	! write file
1129	CALL iom_write_file(tl_fileout, cl_dimorder)
1130
1131	! close file
1132	CALL iom_close(tl_fileout)
1133	CALL file_clean(tl_fileout)
1134
1135	ENDDO ! jk
1136
1137	ENDIF
1138	! clean
1139	CALL boundary_clean(tl_bdy(jl))
1140	ENDDO !jl
1141
1142	! clean
1143	IF( ASSOCIATED(tl_depth%d_value) ) CALL var_clean(tl_depth)
1144	IF( ASSOCIATED(tl_time%d_value) ) CALL var_clean(tl_time)
1145	DEALLOCATE( tl_segdom1 )
1146	DEALLOCATE( tl_segvar1 )
1147	CALL var_clean(tl_seglvl1(:,:,:))
1148	DEALLOCATE( tl_seglvl1 )
1149
1150
1151	CALL mpp_clean(tl_coord1)
1152	CALL mpp_clean(tl_coord0)
1153
1154	CALL multi_clean(tl_multi)
1155
1156	! close log file
1157	CALL logger_footer()
1158	CALL logger_close()
1159	CALL logger_clean()
1160
1161	END SUBROUTINE create__boundary
1162	!-------------------------------------------------------------------
1163	!> @brief
1164	!> This subroutine compute boundary domain for each grid point (T,U,V,F)
1165	!>
1166	!> @author J.Paul
1167	!> @date November, 2013 - Initial Version
1168	!> @date September, 2014
1169	!> - take into account grid point to compute boundary indices
1170	!>
1171	!> @param[in] td_bathy1 file structure
1172	!> @param[in] td_bdy boundary structure
1173	!> @param[in] id_seg segment indice
1174	!> @return array of domain structure
1175	!-------------------------------------------------------------------
1176	FUNCTION create_boundary_get_dom( td_bathy1, td_bdy, id_seg )
1177
1178	IMPLICIT NONE
1179
1180	! Argument
1181	TYPE(TMPP) , INTENT(IN ) :: td_bathy1
1182	TYPE(TBDY) , INTENT(IN ) :: td_bdy
1183	INTEGER(i4), INTENT(IN ) :: id_seg
1184
1185	! function
1186	TYPE(TDOM), DIMENSION(ip_npoint) :: create_boundary_get_dom
1187
1188	! local variable
1189	INTEGER(i4) :: il_imin1
1190	INTEGER(i4) :: il_imax1
1191	INTEGER(i4) :: il_jmin1
1192	INTEGER(i4) :: il_jmax1
1193
1194	INTEGER(i4) :: il_imin
1195	INTEGER(i4) :: il_imax
1196	INTEGER(i4) :: il_jmin
1197	INTEGER(i4) :: il_jmax
1198
1199	INTEGER(i4), DIMENSION(ip_npoint) :: il_ishift
1200	INTEGER(i4), DIMENSION(ip_npoint) :: il_jshift
1201
1202	! loop indices
1203	INTEGER(i4) :: ji
1204	INTEGER(i4) :: jk
1205	!----------------------------------------------------------------
1206	! init
1207	jk=id_seg
1208
1209	il_ishift(:)=0
1210	il_jshift(:)=0
1211
1212	! get boundary definition
1213	SELECT CASE(TRIM(td_bdy%c_card))
1214	CASE('north')
1215
1216	il_imin1=td_bdy%t_seg(jk)%i_first
1217	il_imax1=td_bdy%t_seg(jk)%i_last
1218	il_jmin1=td_bdy%t_seg(jk)%i_index-(td_bdy%t_seg(jk)%i_width-1)
1219	il_jmax1=td_bdy%t_seg(jk)%i_index
1220
1221	il_jshift(jp_V)=-1
1222	il_jshift(jp_F)=-1
1223
1224	CASE('south')
1225
1226	il_imin1=td_bdy%t_seg(jk)%i_first
1227	il_imax1=td_bdy%t_seg(jk)%i_last
1228	il_jmin1=td_bdy%t_seg(jk)%i_index
1229	il_jmax1=td_bdy%t_seg(jk)%i_index+(td_bdy%t_seg(jk)%i_width-1)
1230
1231	CASE('east')
1232
1233	il_imin1=td_bdy%t_seg(jk)%i_index-(td_bdy%t_seg(jk)%i_width-1)
1234	il_imax1=td_bdy%t_seg(jk)%i_index
1235	il_jmin1=td_bdy%t_seg(jk)%i_first
1236	il_jmax1=td_bdy%t_seg(jk)%i_last
1237
1238	il_ishift(jp_U)=-1
1239	il_ishift(jp_F)=-1
1240
1241	CASE('west')
1242
1243	il_imin1=td_bdy%t_seg(jk)%i_index
1244	il_imax1=td_bdy%t_seg(jk)%i_index+(td_bdy%t_seg(jk)%i_width-1)
1245	il_jmin1=td_bdy%t_seg(jk)%i_first
1246	il_jmax1=td_bdy%t_seg(jk)%i_last
1247
1248	END SELECT
1249
1250	!-read fine grid domain
1251	DO ji=1,ip_npoint
1252
1253	! shift domain
1254	il_imin=il_imin1+il_ishift(ji)
1255	il_imax=il_imax1+il_ishift(ji)
1256
1257	il_jmin=il_jmin1+il_jshift(ji)
1258	il_jmax=il_jmax1+il_jshift(ji)
1259
1260	! compute domain
1261	create_boundary_get_dom(ji)=dom_init( td_bathy1, &
1262	& il_imin, il_imax,&
1263	& il_jmin, il_jmax,&
1264	& TRIM(td_bdy%c_card) )
1265
1266	ENDDO
1267
1268	END FUNCTION create_boundary_get_dom
1269	!-------------------------------------------------------------------
1270	!> @brief
1271	!> This subroutine get coordinates over boundary domain
1272	!>
1273	!> @author J.Paul
1274	!> @date November, 2013 - Initial Version
1275	!> @date September, 2014
1276	!> - take into account grid point
1277	!>
1278	!> @param[in] td_coord1 coordinates file structure
1279	!> @param[in] td_dom1 boundary domain structure
1280	!> @param[in] cd_point grid point
1281	!> @param[out] td_lon1 longitude variable structure
1282	!> @param[out] td_lat1 latitude variable structure
1283	!-------------------------------------------------------------------
1284	SUBROUTINE create_boundary_get_coord( td_coord1, td_dom1, cd_point, &
1285	& td_lon1, td_lat1 )
1286
1287	IMPLICIT NONE
1288	! Argument
1289	TYPE(TMPP) , INTENT(IN ) :: td_coord1
1290	TYPE(TDOM) , INTENT(IN ) :: td_dom1
1291	CHARACTER(LEN=*), INTENT(IN ) :: cd_point
1292	TYPE(TVAR) , INTENT( OUT) :: td_lon1
1293	TYPE(TVAR) , INTENT( OUT) :: td_lat1
1294
1295	! local variable
1296	TYPE(TMPP) :: tl_coord1
1297
1298	CHARACTER(LEN=lc) :: cl_name
1299	! loop indices
1300	!----------------------------------------------------------------
1301	!read variables on domain (ugly way to do it, have to work on it)
1302	! init mpp structure
1303	tl_coord1=mpp_copy(td_coord1)
1304
1305	! open mpp files
1306	CALL iom_dom_open(tl_coord1, td_dom1)
1307
1308	! read variable value on domain
1309	WRITE(cl_name,*) 'longitude_'//TRIM(cd_point)
1310	td_lon1=iom_dom_read_var( tl_coord1, TRIM(cl_name), td_dom1)
1311	WRITE(cl_name,*) 'latitude_'//TRIM(cd_point)
1312	td_lat1=iom_dom_read_var( tl_coord1, TRIM(cl_name), td_dom1)
1313
1314	! close mpp files
1315	CALL iom_dom_close(tl_coord1)
1316
1317	! clean structure
1318	CALL mpp_clean(tl_coord1)
1319
1320	END SUBROUTINE create_boundary_get_coord
1321	!-------------------------------------------------------------------
1322	!> @brief
1323	!> This subroutine interpolate variable on boundary
1324	!>
1325	!> @details
1326	!>
1327	!> @author J.Paul
1328	!> @date November, 2013 - Initial Version
1329	!>
1330	!> @param[inout] td_var variable structure
1331	!> @param[in] id_rho array of refinment factor
1332	!> @param[in] id_offset array of offset between fine and coarse grid
1333	!> @param[in] id_iext i-direction size of extra bands (default=im_minext)
1334	!> @param[in] id_jext j-direction size of extra bands (default=im_minext)
1335	!-------------------------------------------------------------------
1336	SUBROUTINE create_boundary_interp( td_var, &
1337	& id_rho, &
1338	& id_offset, &
1339	& id_iext, id_jext )
1340
1341	IMPLICIT NONE
1342
1343	! Argument
1344	TYPE(TVAR) , INTENT(INOUT) :: td_var
1345	INTEGER(I4), DIMENSION(:) , INTENT(IN ) :: id_rho
1346	INTEGER(i4), DIMENSION(:,:), INTENT(IN ) :: id_offset
1347
1348	INTEGER(i4), INTENT(IN ), OPTIONAL :: id_iext
1349	INTEGER(i4), INTENT(IN ), OPTIONAL :: id_jext
1350
1351
1352	! local variable
1353	INTEGER(i4) :: il_iext
1354	INTEGER(i4) :: il_jext
1355	! loop indices
1356	!----------------------------------------------------------------
1357
1358	!WARNING: at least two extrabands are required for cubic interpolation
1359	il_iext=2
1360	IF( PRESENT(id_iext) ) il_iext=id_iext
1361
1362	il_jext=2
1363	IF( PRESENT(id_jext) ) il_jext=id_jext
1364
1365	IF( il_iext < 2 .AND. td_var%c_interp(1) == 'cubic' )THEN
1366	CALL logger_warn("CREATE BOUNDARY INTERP: at least extrapolation "//&
1367	& "on two points are required with cubic interpolation ")
1368	il_iext=2
1369	ENDIF
1370
1371	IF( il_jext < 2 .AND. td_var%c_interp(1) == 'cubic' )THEN
1372	CALL logger_warn("CREATE BOUNDARY INTERP: at least extrapolation "//&
1373	& "on two points are required with cubic interpolation ")
1374	il_jext=2
1375	ENDIF
1376
1377	! work on variable
1378	! add extraband
1379	CALL extrap_add_extrabands(td_var, il_iext, il_jext)
1380
1381	! extrapolate variable
1382	CALL extrap_fill_value( td_var )
1383
1384	! interpolate Bathymetry
1385	CALL interp_fill_value( td_var, id_rho(:), &
1386	& id_offset=id_offset(:,:) )
1387
1388	! remove extraband
1389	CALL extrap_del_extrabands(td_var, il_iext*id_rho(jp_I), &
1390	& il_jext*id_rho(jp_J))
1391
1392	END SUBROUTINE create_boundary_interp
1393	!-------------------------------------------------------------------
1394	!> @brief
1395	!> This function create variable, filled with matrix value
1396	!>
1397	!> @details
1398	!> A variable is create with the same name that the input variable,
1399	!> and with dimension of the coordinate file.
1400	!> Then the variable array of value is split into equal subdomain.
1401	!> Each subdomain is fill with the associated value of the matrix.
1402	!>
1403	!> @author J.Paul
1404	!> @date November, 2013 - Initial Version
1405	!>
1406	!> @param[in] td_var variable structure
1407	!> @param[in] td_dom domain structure
1408	!> @param[in] id_nlevel number of levels
1409	!> @return variable structure
1410	!-------------------------------------------------------------------
1411	FUNCTION create_boundary_matrix(td_var, td_dom, id_nlevel)
1412	IMPLICIT NONE
1413	! Argument
1414	TYPE(TVAR) , INTENT(IN) :: td_var
1415	TYPE(TDOM) , INTENT(IN) :: td_dom
1416	INTEGER(i4), INTENT(IN) :: id_nlevel
1417
1418	! function
1419	TYPE(TVAR) :: create_boundary_matrix
1420
1421	! local variable
1422	INTEGER(i4) , DIMENSION(3) :: il_dim
1423	INTEGER(i4) , DIMENSION(3) :: il_size
1424	INTEGER(i4) , DIMENSION(3) :: il_rest
1425
1426	INTEGER(i4) , DIMENSION(:) , ALLOCATABLE :: il_ishape
1427	INTEGER(i4) , DIMENSION(:) , ALLOCATABLE :: il_jshape
1428	INTEGER(i4) , DIMENSION(:) , ALLOCATABLE :: il_kshape
1429
1430	REAL(dp) , DIMENSION(:,:,:,:), ALLOCATABLE :: dl_value
1431
1432	TYPE(TDIM) , DIMENSION(ip_maxdim) :: tl_dim
1433
1434	! loop indices
1435	INTEGER(i4) :: ji
1436	INTEGER(i4) :: jj
1437	INTEGER(i4) :: jk
1438	!----------------------------------------------------------------
1439
1440	! write value on grid
1441	! get matrix dimension
1442	il_dim(:)=td_var%t_dim(1:3)%i_len
1443
1444	tl_dim(jp_I:jp_J)=dim_copy(td_dom%t_dim(jp_I:jp_J))
1445	tl_dim(jp_K)%i_len=id_nlevel
1446
1447	! split output domain in N subdomain depending of matrix dimension
1448	il_size(:) = tl_dim(1:3)%i_len / il_dim(:)
1449	il_rest(:) = MOD(tl_dim(1:3)%i_len, il_dim(:))
1450
1451	ALLOCATE( il_ishape(il_dim(1)+1) )
1452	il_ishape(:)=0
1453	DO ji=2,il_dim(1)+1
1454	il_ishape(ji)=il_ishape(ji-1)+il_size(1)
1455	ENDDO
1456	! add rest to last cell
1457	il_ishape(il_dim(1)+1)=il_ishape(il_dim(1)+1)+il_rest(1)
1458
1459	ALLOCATE( il_jshape(il_dim(2)+1) )
1460	il_jshape(:)=0
1461	DO jj=2,il_dim(2)+1
1462	il_jshape(jj)=il_jshape(jj-1)+il_size(2)
1463	ENDDO
1464	! add rest to last cell
1465	il_jshape(il_dim(2)+1)=il_jshape(il_dim(2)+1)+il_rest(2)
1466
1467	ALLOCATE( il_kshape(il_dim(3)+1) )
1468	il_kshape(:)=0
1469	DO jk=2,il_dim(3)+1
1470	il_kshape(jk)=il_kshape(jk-1)+il_size(3)
1471	ENDDO
1472	! add rest to last cell
1473	il_kshape(il_dim(3)+1)=il_kshape(il_dim(3)+1)+il_rest(3)
1474
1475	! write ouput array of value
1476	ALLOCATE(dl_value( tl_dim(1)%i_len, &
1477	& tl_dim(2)%i_len, &
1478	& tl_dim(3)%i_len, &
1479	& tl_dim(4)%i_len) )
1480
1481	dl_value(:,:,:,:)=0
1482
1483	DO jk=2,il_dim(3)+1
1484	DO jj=2,il_dim(2)+1
1485	DO ji=2,il_dim(1)+1
1486
1487	dl_value( 1+il_ishape(ji-1):il_ishape(ji), &
1488	& 1+il_jshape(jj-1):il_jshape(jj), &
1489	& 1+il_kshape(jk-1):il_kshape(jk), &
1490	& 1 ) = td_var%d_value(ji-1,jj-1,jk-1,1)
1491
1492	ENDDO
1493	ENDDO
1494	ENDDO
1495
1496	! initialise variable with value
1497	create_boundary_matrix=var_init(TRIM(td_var%c_name),dl_value(:,:,:,:))
1498
1499	DEALLOCATE(dl_value)
1500
1501	END FUNCTION create_boundary_matrix
1502	!-------------------------------------------------------------------
1503	!> @brief
1504	!> This subroutine use mask to filled land point with _FillValue
1505	!>
1506	!> @details
1507	!>
1508	!> @author J.Paul
1509	!> @date November, 2013 - Initial Version
1510	!>
1511	!> @param[inout] td_var variable structure
1512	!> @param[in] td_mask mask variable structure
1513	!-------------------------------------------------------------------
1514	SUBROUTINE create_boundary_use_mask( td_var, td_mask )
1515
1516	IMPLICIT NONE
1517
1518	! Argument
1519	TYPE(TVAR), INTENT(INOUT) :: td_var
1520	TYPE(TVAR), INTENT(IN ) :: td_mask
1521
1522	! local variable
1523	INTEGER(i4), DIMENSION(:,:), ALLOCATABLE :: il_mask
1524
1525	! loop indices
1526	INTEGER(i4) :: jk
1527	INTEGER(i4) :: jl
1528	!----------------------------------------------------------------
1529
1530	IF( ANY(td_var%t_dim(1:2)%i_len /= &
1531	& td_mask%t_dim(1:2)%i_len) )THEN
1532	CALL logger_debug(" mask dimension ( "//&
1533	& TRIM(fct_str(td_mask%t_dim(1)%i_len))//","//&
1534	& TRIM(fct_str(td_mask%t_dim(2)%i_len))//")" )
1535	CALL logger_debug(" variable dimension ( "//&
1536	& TRIM(fct_str(td_var%t_dim(1)%i_len))//","//&
1537	& TRIM(fct_str(td_var%t_dim(2)%i_len))//")" )
1538	CALL logger_fatal("CREATE BOUNDARY USE MASK: mask and "//&
1539	& "variable dimension differ." )
1540	ENDIF
1541
1542	ALLOCATE( il_mask(td_var%t_dim(1)%i_len, &
1543	& td_var%t_dim(2)%i_len) )
1544
1545	il_mask(:,:)=INT(td_mask%d_value(:,:,1,1))
1546
1547	DO jl=1,td_var%t_dim(4)%i_len
1548	DO jk=1,td_var%t_dim(3)%i_len
1549	WHERE( il_mask(:,:) < jk ) td_var%d_value(:,:,jk,jl)=td_var%d_fill
1550	ENDDO
1551	ENDDO
1552
1553	DEALLOCATE( il_mask )
1554
1555	END SUBROUTINE create_boundary_use_mask
1556	!-------------------------------------------------------------------
1557	!> @brief
1558	!> This function extract level over domain on each grid point, and return
1559	!> array of variable structure
1560	!>
1561	!> @author J.Paul
1562	!> @date November, 2013 - Initial Version
1563	!>
1564	!> @param[in] td_level array of level variable structure
1565	!> @param[in] td_dom array of domain structure
1566	!> @return array of variable structure
1567	!-------------------------------------------------------------------
1568	FUNCTION create_boundary_get_level(td_level, td_dom)
1569	IMPLICIT NONE
1570	! Argument
1571	TYPE(TVAR), DIMENSION(:), INTENT(IN) :: td_level
1572	TYPE(TDOM), DIMENSION(:), INTENT(IN) :: td_dom
1573
1574	! function
1575	TYPE(TVAR), DIMENSION(ip_npoint) :: create_boundary_get_level
1576
1577	! local variable
1578	TYPE(TVAR), DIMENSION(ip_npoint) :: tl_var
1579
1580	! loop indices
1581	INTEGER(i4) :: ji
1582	!----------------------------------------------------------------
1583
1584	IF( SIZE(td_level(:)) /= ip_npoint .OR. &
1585	& SIZE(td_dom(:)) /= ip_npoint )THEN
1586	CALL logger_error("CREATE BDY GET LEVEL: invalid dimension. "//&
1587	& "check input array of level and domain.")
1588	ELSE
1589
1590	DO ji=1,ip_npoint
1591
1592	tl_var(ji)=var_copy(td_level(ji))
1593
1594	IF( ASSOCIATED(tl_var(ji)%d_value) ) DEALLOCATE(tl_var(ji)%d_value)
1595
1596	tl_var(ji)%t_dim(1)%i_len=td_dom(ji)%t_dim(1)%i_len
1597	tl_var(ji)%t_dim(2)%i_len=td_dom(ji)%t_dim(2)%i_len
1598	ALLOCATE(tl_var(ji)%d_value(tl_var(ji)%t_dim(1)%i_len, &
1599	& tl_var(ji)%t_dim(2)%i_len, &
1600	& tl_var(ji)%t_dim(3)%i_len, &
1601	& tl_var(ji)%t_dim(4)%i_len) )
1602
1603	tl_var(ji)%d_value(:,:,:,:) = &
1604	& td_level(ji)%d_value( td_dom(ji)%i_imin:td_dom(ji)%i_imax, &
1605	& td_dom(ji)%i_jmin:td_dom(ji)%i_jmax, :, : )
1606
1607	ENDDO
1608	! save result
1609	create_boundary_get_level(:)=var_copy(tl_var(:))
1610
1611	! clean
1612	CALL var_clean(tl_var(:))
1613
1614	ENDIF
1615	END FUNCTION create_boundary_get_level
1616	!-------------------------------------------------------------------
1617	!> @brief
1618	!> This subroutine check if variable need depth dimension,
1619	!> get depth variable value in an open mpp structure
1620	!> and check if agree with already input depth variable.
1621	!>
1622	!> @details
1623	!>
1624	!> @author J.Paul
1625	!> @date November, 2014 - Initial Version
1626	!> @date January, 2016
1627	!> - check if variable need/use depth dimension
1628	!>
1629	!> @param[in] td_var variable structure
1630	!> @param[in] td_mpp mpp structure
1631	!> @param[in] id_nlevel mpp structure
1632	!> @param[inout] td_depth depth variable structure
1633	!-------------------------------------------------------------------
1634	SUBROUTINE create_boundary_check_depth( td_var, td_mpp, id_nlevel, td_depth )
1635
1636	IMPLICIT NONE
1637
1638	! Argument
1639	TYPE(TVAR) , INTENT(IN ) :: td_var
1640	TYPE(TMPP) , INTENT(IN ) :: td_mpp
1641	INTEGER(i4), INTENT(IN ) :: id_nlevel
1642	TYPE(TVAR) , INTENT(INOUT) :: td_depth
1643
1644	! local variable
1645	INTEGER(i4) :: il_varid
1646	TYPE(TVAR) :: tl_depth
1647	! loop indices
1648	!----------------------------------------------------------------
1649
1650	IF( td_var%t_dim(jp_K)%l_use .AND. &
1651	& ( TRIM(td_var%c_axis) == '' .OR. &
1652	& INDEX(TRIM(td_var%c_axis),'Z') /= 0 )&
1653	& )THEN
1654
1655	! check vertical dimension
1656	IF( td_mpp%t_dim(jp_K)%l_use )THEN
1657	IF( td_mpp%t_dim(jp_K)%i_len /= id_nlevel .AND. &
1658	& td_mpp%t_dim(jp_K)%i_len /= 1 )THEN
1659	CALL logger_error("CREATE BOUNDARY: dimension in file "//&
1660	& TRIM(td_mpp%c_name)//" not agree with namelist in_nlevel ")
1661	ENDIF
1662	ENDIF
1663
1664	! get or check depth value
1665	IF( td_mpp%t_proc(1)%i_depthid /= 0 )THEN
1666
1667	il_varid=td_mpp%t_proc(1)%i_depthid
1668	IF( ASSOCIATED(td_depth%d_value) )THEN
1669
1670	tl_depth=iom_mpp_read_var(td_mpp, il_varid)
1671
1672	IF( ANY( td_depth%d_value(:,:,:,:) /= &
1673	& tl_depth%d_value(:,:,:,:) ) )THEN
1674
1675	CALL logger_error("CREATE BOUNDARY: depth value "//&
1676	& "for variable "//TRIM(td_var%c_name)//&
1677	& "from "//TRIM(td_mpp%c_name)//" not conform "//&
1678	& " to those from former file(s).")
1679
1680	ENDIF
1681	CALL var_clean(tl_depth)
1682
1683	ELSE
1684	td_depth=iom_mpp_read_var(td_mpp,il_varid)
1685	ENDIF
1686
1687	ENDIF
1688	ELSE
1689	CALL logger_debug("CREATE BOUNDARY: no depth dimension use"//&
1690	& " for variable "//TRIM(td_var%c_name))
1691	ENDIF
1692
1693	END SUBROUTINE create_boundary_check_depth
1694	!-------------------------------------------------------------------
1695	!> @brief
1696	!> This subroutine check if variable need time dimension,
1697	!> get date and time in an open mpp structure
1698	!> and check if agree with date and time already read.
1699	!>
1700	!> @details
1701	!>
1702	!> @author J.Paul
1703	!> @date November, 2014 - Initial Version
1704	!> @date January, 2016
1705	!> - check if variable need/use time dimension
1706	!>
1707	!> @param[in] td_var variable structure
1708	!> @param[in] td_mpp mpp structure
1709	!> @param[inout] td_time time variable structure
1710	!-------------------------------------------------------------------
1711	SUBROUTINE create_boundary_check_time( td_var, td_mpp, td_time )
1712
1713	IMPLICIT NONE
1714
1715	! Argument
1716	TYPE(TVAR), INTENT(IN ) :: td_var
1717	TYPE(TMPP), INTENT(IN ) :: td_mpp
1718	TYPE(TVAR), INTENT(INOUT) :: td_time
1719
1720	! local variable
1721	INTEGER(i4) :: il_varid
1722	TYPE(TVAR) :: tl_time
1723
1724	TYPE(TDATE) :: tl_date1
1725	TYPE(TDATE) :: tl_date2
1726	! loop indices
1727	!----------------------------------------------------------------
1728	IF( td_var%t_dim(jp_L)%l_use .AND. &
1729	& ( TRIM(td_var%c_axis) == '' .OR. &
1730	& INDEX(TRIM(td_var%c_axis),'T') /= 0 )&
1731	& )THEN
1732
1733	! get or check depth value
1734	IF( td_mpp%t_proc(1)%i_timeid /= 0 )THEN
1735
1736	il_varid=td_mpp%t_proc(1)%i_timeid
1737	IF( ASSOCIATED(td_time%d_value) )THEN
1738
1739	tl_time=iom_mpp_read_var(td_mpp, il_varid)
1740
1741	tl_date1=var_to_date(td_time)
1742	tl_date2=var_to_date(tl_time)
1743	IF( tl_date1 - tl_date2 /= 0 )THEN
1744
1745	CALL logger_warn("CREATE BOUNDARY: date from "//&
1746	& TRIM(td_mpp%c_name)//" not conform "//&
1747	& " to those from former file(s).")
1748
1749	ENDIF
1750	CALL var_clean(tl_time)
1751
1752	ELSE
1753	td_time=iom_mpp_read_var(td_mpp,il_varid)
1754	ENDIF
1755
1756	ENDIF
1757
1758	ELSE
1759	CALL logger_debug("CREATE BOUNDARY: no time dimension use"//&
1760	& " for variable "//TRIM(td_var%c_name))
1761	ENDIF
1762
1763	END SUBROUTINE create_boundary_check_time
1764	END PROGRAM create_boundary

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

Download in other formats: