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_bathy.f90 in NEMO/trunk/tools/SIREN/src – NEMO

Context Navigation

source: NEMO/trunk/tools/SIREN/src/create_bathy.f90 @ 9598

Last change on this file since 9598 was 9598, checked in by nicolasmartin, 6 years ago
Reorganisation plan for NEMO repository: changes to make compilation succeed with new structure Juste one issue left with AGRIF_NORDIC with AGRIF preprocessing Standardisation of routines header with version 4.0 and year 2018 Fix for some broken symlinks
File size: 40.8 KB

Line
1	!----------------------------------------------------------------------
2	! NEMO system team, System and Interface for oceanic RElocable Nesting
3	!----------------------------------------------------------------------
4	!
5	! PROGRAM: create_bathy
6	!
7	! DESCRIPTION:
8	!> @file
9	!> @brief
10	!> This program creates fine grid bathymetry file.
11	!>
12	!> @details
13	!> @section sec1 method
14	!> Bathymetry could be extracted from fine grid Bathymetry file, interpolated
15	!> from coarse grid Bathymetry file, or manually written.
16	!>
17	!> @section sec2 how to
18	!> to create fine grid bathymetry file:<br/>
19	!> @code{.sh}
20	!> ./SIREN/bin/create_bathy create_bathy.nam
21	!> @endcode
22	!> <br/>
23	!> \image html bathy_40.png
24	!> <center>\image latex bathy_30.png
25	!> </center>
26	!>
27	!> @note
28	!> you could find a template of the namelist in templates directory.
29	!>
30	!> create_bathy.nam contains 7 namelists:<br/>
31	!> - logger namelist (namlog)
32	!> - config namelist (namcfg)
33	!> - coarse grid namelist (namcrs)
34	!> - fine grid namelist (namfin)
35	!> - variable namelist (namvar)
36	!> - nesting namelist (namnst)
37	!> - output namelist (namout)
38	!>
39	!> * _logger namelist (namlog)_:<br/>
40	!> - cn_logfile : log filename
41	!> - cn_verbosity : verbosity ('trace','debug','info',
42	!> 'warning','error','fatal','none')
43	!> - in_maxerror : maximum number of error allowed
44	!>
45	!> * _config namelist (namcfg)_:<br/>
46	!> - cn_varcfg : variable configuration file
47	!> (see ./SIREN/cfg/variable.cfg)
48	!> - cn_dimcfg : dimension configuration file. defines dimension allowed
49	!> (see ./SIREN/cfg/dimension.cfg).
50	!> - cn_dumcfg : useless (dummy) configuration file, for useless
51	!> dimension or variable (see ./SIREN/cfg/dummy.cfg).
52	!>
53	!> * _coarse grid namelist (namcrs)_:<br/>
54	!> - cn_coord0 : coordinate file
55	!> - in_perio0 : NEMO periodicity index (see Model Boundary Condition in
56	!> [NEMO documentation](http://www.nemo-ocean.eu/About-NEMO/Reference-manuals))
57	!>
58	!> * _fine grid namelist (namfin)_:<br/>
59	!> - cn_coord1 : coordinate file
60	!> - in_perio1 : periodicity index
61	!> - ln_fillclosed : fill closed sea or not (default is .TRUE.)
62	!>
63	!> * _variable namelist (namvar)_:<br/>
64	!> - cn_varfile : list of variable, and corresponding file.<br/>
65	!> cn_varfile is the path and filename of the file where find
66	!> variable.
67	!> @note
68	!> cn_varfile could be a matrix of value, if you want to filled
69	!> manually variable value.<br/>
70	!> the variable array of value is split into equal subdomain.<br/>
71	!> Each subdomain is filled with the corresponding value
72	!> of the matrix.<br/>
73	!> separators used to defined matrix are:
74	!> - ',' for line
75	!> - '/' for row
76	!> Example:<br/>
77	!> 3,2,3/1,4,5 => @f$ \left( \begin{array}{ccc}
78	!> 3 & 2 & 3 \\
79	!> 1 & 4 & 5 \end{array} \right) @f$
80	!>
81	!> Examples:
82	!> - 'Bathymetry:gridT.nc'
83	!> - 'Bathymetry:5000,5000,5000/5000,3000,5000/5000,5000,5000'
84	!>
85	!> - cn_varinfo : list of variable and extra information about request(s)
86	!> to be used.<br/>
87	!> each elements of cn_varinfo is a string character
88	!> (separated by ',').<br/>
89	!> it is composed of the variable name follow by ':',
90	!> then request(s) to be used on this variable.<br/>
91	!> request could be:
92	!> - int = interpolation method
93	!> - ext = extrapolation method
94	!> - flt = filter method
95	!> - min = minimum value
96	!> - max = maximum value
97	!> - unt = new units
98	!> - unf = unit scale factor (linked to new units)
99	!>
100	!> requests must be separated by ';'.<br/>
101	!> order of requests does not matter.<br/>
102	!>
103	!> informations about available method could be find in @ref interp,
104	!> @ref extrap and @ref filter modules.<br/>
105	!> Example: 'Bathymetry: flt=2*hamming(2,3); min=0'
106	!> @note
107	!> If you do not specify a method which is required,
108	!> default one is apply.
109	!> @warning
110	!> variable name must be __Bathymetry__ here.
111	!>
112	!> * _nesting namelist (namnst)_:<br/>
113	!> - in_rhoi : refinement factor in i-direction
114	!> - in_rhoj : refinement factor in j-direction
115	!> @note
116	!> coarse grid indices will be deduced from fine grid
117	!> coordinate file.
118	!>
119	!> * _output namelist (namout)_:<br/>
120	!> - cn_fileout : output bathymetry file
121	!>
122	!> @author J.Paul
123	! REVISION HISTORY:
124	!> @date November, 2013 - Initial Version
125	!> @date Sepember, 2014
126	!> - add header for user
127	!> - Bug fix, compute offset depending of grid point
128	!> @date June, 2015
129	!> - extrapolate all land points.
130	!> - allow to change unit.
131	!> @date September, 2015
132	!> - manage useless (dummy) variable, attributes, and dimension
133	!> @date January,2016
134	!> - add create_bathy_check_depth as in create_boundary
135	!> - add create_bathy_check_time as in create_boundary
136	!> @date February, 2016
137	!> - do not closed sea for east-west cyclic domain
138	!> @date October, 2016
139	!> - dimension to be used select from configuration file
140	!
141	!> @todo
142	!> - check tl_multi is not empty
143	!>
144	!> @note Software governed by the CeCILL licence (./LICENSE)
145	!----------------------------------------------------------------------
146	PROGRAM create_bathy
147
148	USE global ! global variable
149	USE kind ! F90 kind parameter
150	USE logger ! log file manager
151	USE fct ! basic useful function
152	USE date ! date manager
153	USE att ! attribute manager
154	USE dim ! dimension manager
155	USE var ! variable manager
156	USE file ! file manager
157	USE multi ! multi file manager
158	USE iom ! I/O manager
159	USE grid ! grid manager
160	USE extrap ! extrapolation manager
161	USE interp ! interpolation manager
162	USE filter ! filter manager
163	USE mpp ! MPP manager
164	USE dom ! domain manager
165	USE iom_mpp ! MPP I/O manager
166	USE iom_dom ! DOM I/O manager
167
168	IMPLICIT NONE
169
170	! local variable
171	CHARACTER(LEN=lc) :: cl_namelist
172	CHARACTER(LEN=lc) :: cl_date
173	CHARACTER(LEN=lc) :: cl_data
174
175	INTEGER(i4) :: il_narg
176	INTEGER(i4) :: il_status
177	INTEGER(i4) :: il_fileid
178	INTEGER(i4) :: il_attid
179	INTEGER(i4) :: il_imin0
180	INTEGER(i4) :: il_imax0
181	INTEGER(i4) :: il_jmin0
182	INTEGER(i4) :: il_jmax0
183	INTEGER(i4) , DIMENSION(ip_maxdim) :: il_rho
184	INTEGER(i4) , DIMENSION(2,2) :: il_offset
185	INTEGER(i4) , DIMENSION(2,2) :: il_ind
186	INTEGER(i4) , DIMENSION(:,:) , ALLOCATABLE :: il_mask
187
188	LOGICAL :: ll_exist
189	LOGICAL :: ll_fillclosed
190
191	TYPE(TMPP) :: tl_coord0
192	TYPE(TMPP) :: tl_coord1
193	TYPE(TMPP) :: tl_mpp
194	TYPE(TFILE) :: tl_fileout
195
196	TYPE(TATT) :: tl_att
197
198	TYPE(TVAR) :: tl_lon
199	TYPE(TVAR) :: tl_lat
200	TYPE(TVAR) :: tl_depth
201	TYPE(TVAR) :: tl_time
202
203	TYPE(TVAR) :: tl_tmp
204	TYPE(TVAR) , DIMENSION(:), ALLOCATABLE :: tl_var
205
206	TYPE(TDIM) , DIMENSION(ip_maxdim) :: tl_dim
207
208	TYPE(TMULTI) :: tl_multi
209
210	REAL(dp) :: dl_minbat
211
212	! loop indices
213	INTEGER(i4) :: ji
214	INTEGER(i4) :: jj
215	INTEGER(i4) :: jk
216
217	! namelist variable
218	! namlog
219	CHARACTER(LEN=lc) :: cn_logfile = 'create_bathy.log'
220	CHARACTER(LEN=lc) :: cn_verbosity = 'warning'
221	INTEGER(i4) :: in_maxerror = 5
222
223	! namcfg
224	CHARACTER(LEN=lc) :: cn_varcfg = './cfg/variable.cfg'
225	CHARACTER(LEN=lc) :: cn_dimcfg = './cfg/dimension.cfg'
226	CHARACTER(LEN=lc) :: cn_dumcfg = './cfg/dummy.cfg'
227
228	! namcrs
229	CHARACTER(LEN=lc) :: cn_coord0 = ''
230	INTEGER(i4) :: in_perio0 = -1
231
232	! namfin
233	CHARACTER(LEN=lc) :: cn_coord1 = ''
234	INTEGER(i4) :: in_perio1 = -1
235	LOGICAL :: ln_fillclosed = .TRUE.
236
237	! namvar
238	CHARACTER(LEN=lc), DIMENSION(ip_maxvar) :: cn_varfile = ''
239	CHARACTER(LEN=lc), DIMENSION(ip_maxvar) :: cn_varinfo = ''
240
241	! namnst
242	INTEGER(i4) :: in_rhoi = 1
243	INTEGER(i4) :: in_rhoj = 1
244
245	! namout
246	CHARACTER(LEN=lc) :: cn_fileout = 'bathy_fine.nc'
247	!-------------------------------------------------------------------
248
249	NAMELIST /namlog/ & !< logger namelist
250	& cn_logfile, & !< log file
251	& cn_verbosity, & !< log verbosity
252	& in_maxerror !< logger maximum error
253
254	NAMELIST /namcfg/ & !< configuration namelist
255	& cn_varcfg, & !< variable configuration file
256	& cn_dimcfg, & !< dimension configuration file
257	& cn_dumcfg !< dummy configuration file
258
259	NAMELIST /namcrs/ & !< coarse grid namelist
260	& cn_coord0, & !< coordinate file
261	& in_perio0 !< periodicity index
262
263	NAMELIST /namfin/ & !< fine grid namelist
264	& cn_coord1, & !< coordinate file
265	& in_perio1, & !< periodicity index
266	& ln_fillclosed !< fill closed sea
267
268	NAMELIST /namvar/ & !< variable namelist
269	& cn_varfile, & !< list of variable file
270	& cn_varinfo !< list of variable and interpolation method to be used. (ex: 'votemper:linear','vosaline:cubic' )
271
272	NAMELIST /namnst/ & !< nesting namelist
273	& in_rhoi, & !< refinement factor in i-direction
274	& in_rhoj !< refinement factor in j-direction
275
276	NAMELIST /namout/ & !< output namlist
277	& cn_fileout !< fine grid bathymetry file
278	!-------------------------------------------------------------------
279
280	! namelist
281	! get namelist
282	il_narg=COMMAND_ARGUMENT_COUNT() !f03 intrinsec
283	IF( il_narg/=1 )THEN
284	PRINT *,"ERROR in create_bathy: need a namelist"
285	STOP
286	ELSE
287	CALL GET_COMMAND_ARGUMENT(1,cl_namelist) !f03 intrinsec
288	ENDIF
289
290	! read namelist
291	INQUIRE(FILE=TRIM(cl_namelist), EXIST=ll_exist)
292	IF( ll_exist )THEN
293
294	il_fileid=fct_getunit()
295
296	OPEN( il_fileid, FILE=TRIM(cl_namelist), &
297	& FORM='FORMATTED', &
298	& ACCESS='SEQUENTIAL', &
299	& STATUS='OLD', &
300	& ACTION='READ', &
301	& IOSTAT=il_status)
302	CALL fct_err(il_status)
303	IF( il_status /= 0 )THEN
304	PRINT *,"ERROR in create_bathy: error opening "//TRIM(cl_namelist)
305	STOP
306	ENDIF
307
308	READ( il_fileid, NML = namlog )
309	! define log file
310	CALL logger_open(TRIM(cn_logfile),TRIM(cn_verbosity),in_maxerror)
311	CALL logger_header()
312
313	READ( il_fileid, NML = namcfg )
314	! get variable extra information
315	CALL var_def_extra(TRIM(cn_varcfg))
316
317	! get dimension allowed
318	CALL dim_def_extra(TRIM(cn_dimcfg))
319
320	! get dummy variable
321	CALL var_get_dummy(TRIM(cn_dumcfg))
322	! get dummy dimension
323	CALL dim_get_dummy(TRIM(cn_dumcfg))
324	! get dummy attribute
325	CALL att_get_dummy(TRIM(cn_dumcfg))
326
327	READ( il_fileid, NML = namcrs )
328	READ( il_fileid, NML = namfin )
329	READ( il_fileid, NML = namvar )
330	! add user change in extra information
331	CALL var_chg_extra( cn_varinfo )
332	! match variable with file
333	tl_multi=multi_init(cn_varfile)
334
335	READ( il_fileid, NML = namnst )
336	READ( il_fileid, NML = namout )
337
338	CLOSE( il_fileid, IOSTAT=il_status )
339	CALL fct_err(il_status)
340	IF( il_status /= 0 )THEN
341	CALL logger_error("CREATE BATHY: closing "//TRIM(cl_namelist))
342	ENDIF
343
344	ELSE
345
346	PRINT *,"ERROR in create_bathy: can't find "//TRIM(cl_namelist)
347	STOP
348
349	ENDIF
350
351	CALL multi_print(tl_multi)
352
353	! open files
354	IF( cn_coord0 /= '' )THEN
355	tl_coord0=mpp_init( file_init(TRIM(cn_coord0)), id_perio=in_perio0)
356	CALL grid_get_info(tl_coord0)
357	ELSE
358	CALL logger_fatal("CREATE BATHY: no coarse grid coordinate found. "//&
359	& "check namelist")
360	ENDIF
361
362	IF( TRIM(cn_coord1) /= '' )THEN
363	tl_coord1=mpp_init( file_init(TRIM(cn_coord1)),id_perio=in_perio1)
364	CALL grid_get_info(tl_coord1)
365	ELSE
366	CALL logger_fatal("CREATE BATHY: no fine grid coordinate found. "//&
367	& "check namelist")
368	ENDIF
369
370	! do not closed sea for east-west cyclic domain
371	ll_fillclosed=ln_fillclosed
372	IF( tl_coord1%i_perio == 1 ) ll_fillclosed=.FALSE.
373
374	! check
375	! check output file do not already exist
376	INQUIRE(FILE=TRIM(cn_fileout), EXIST=ll_exist)
377	IF( ll_exist )THEN
378	CALL logger_fatal("CREATE BATHY: output file "//TRIM(cn_fileout)//&
379	& " already exist.")
380	ENDIF
381
382	! check namelist
383	! check refinement factor
384	il_rho(:)=1
385	IF( in_rhoi < 1 .OR. in_rhoj < 1 )THEN
386	CALL logger_error("CREATE BATHY: invalid refinement factor."//&
387	& " check namelist "//TRIM(cl_namelist))
388	ELSE
389	il_rho(jp_I)=in_rhoi
390	il_rho(jp_J)=in_rhoj
391	ENDIF
392
393	! check domain indices
394	! compute coarse grid indices around fine grid
395	il_ind(:,:)=grid_get_coarse_index( tl_coord0, tl_coord1, &
396	& id_rho=il_rho(:) )
397
398	il_imin0=il_ind(jp_I,1) ; il_imax0=il_ind(jp_I,2)
399	il_jmin0=il_ind(jp_J,1) ; il_jmax0=il_ind(jp_J,2)
400
401	! check domain validity
402	CALL grid_check_dom(tl_coord0, il_imin0, il_imax0, il_jmin0, il_jmax0)
403
404	! check coincidence between coarse and fine grid
405	CALL grid_check_coincidence( tl_coord0, tl_coord1, &
406	& il_imin0, il_imax0, &
407	& il_jmin0, il_jmax0, &
408	& il_rho(:) )
409
410	IF( .NOT. ASSOCIATED(tl_multi%t_mpp) )THEN
411	CALL logger_error("CREATE BATHY: no mpp file to work on. "//&
412	& "check cn_varfile in namelist.")
413	ELSE
414
415	ALLOCATE( tl_var( tl_multi%i_nvar ) )
416	jk=0
417	DO ji=1,tl_multi%i_nmpp
418
419	WRITE(cl_data,'(a,i2.2)') 'data-',jk+1
420	IF( .NOT. ASSOCIATED(tl_multi%t_mpp(ji)%t_proc(1)%t_var) )THEN
421
422	CALL logger_fatal("CREATE BATHY: no variable to work on for "//&
423	& "mpp file"//TRIM(tl_multi%t_mpp(ji)%c_name)//&
424	& ". check cn_varfile in namelist.")
425
426	ELSEIF( TRIM(tl_multi%t_mpp(ji)%c_name) == TRIM(cl_data) )THEN
427
428	!- use input matrix to initialise variable
429	DO jj=1,tl_multi%t_mpp(ji)%t_proc(1)%i_nvar
430	jk=jk+1
431	tl_tmp=var_copy(tl_multi%t_mpp(ji)%t_proc(1)%t_var(jj))
432
433	tl_var(jk)=create_bathy_matrix(tl_tmp, tl_coord1)
434	ENDDO
435	! clean
436	CALL var_clean(tl_tmp)
437
438	ELSE
439
440	tl_mpp=mpp_init( file_init(TRIM(tl_multi%t_mpp(ji)%c_name)) )
441	CALL grid_get_info(tl_mpp)
442
443	! open mpp file
444	CALL iom_mpp_open(tl_mpp)
445
446	! get or check depth value
447	CALL create_bathy_check_depth( tl_mpp, tl_depth )
448
449	! get or check time value
450	CALL create_bathy_check_time( tl_mpp, tl_time )
451
452	! close mpp file
453	CALL iom_mpp_close(tl_mpp)
454
455	IF( ANY(tl_mpp%t_dim(1:2)%i_len /= tl_coord0%t_dim(1:2)%i_len).OR.&
456	& ALL(il_rho(:)==1) )THEN
457	!- extract bathymetry from fine grid bathymetry
458	DO jj=1,tl_multi%t_mpp(ji)%t_proc(1)%i_nvar
459	jk=jk+1
460	tl_tmp=var_copy(tl_multi%t_mpp(ji)%t_proc(1)%t_var(jj))
461
462	tl_var(jk)=create_bathy_extract( tl_tmp, tl_mpp, &
463	& tl_coord1 )
464	ENDDO
465	! clean
466	CALL var_clean(tl_tmp)
467	ELSE
468	!- get bathymetry from coarse grid bathymetry
469	DO jj=1,tl_multi%t_mpp(ji)%t_proc(1)%i_nvar
470	jk=jk+1
471	tl_tmp=var_copy(tl_multi%t_mpp(ji)%t_proc(1)%t_var(jj))
472
473	il_offset(:,:)= grid_get_fine_offset(tl_coord0, &
474	& il_imin0, il_jmin0, &
475	& il_imax0, il_jmax0, &
476	& tl_coord1, &
477	& il_rho(:), &
478	& TRIM(tl_tmp%c_point) )
479
480	tl_var(jk)=create_bathy_get_var( tl_tmp, tl_mpp, &
481	& il_imin0, il_jmin0, &
482	& il_imax0, il_jmax0, &
483	& il_offset(:,:), &
484	& il_rho(:) )
485	ENDDO
486	! clean
487	CALL var_clean(tl_tmp)
488	ENDIF
489
490	! clean structure
491	CALL mpp_clean(tl_mpp)
492
493	ENDIF
494	ENDDO
495	ENDIF
496
497	! use additional request
498	DO jk=1,tl_multi%i_nvar
499
500	! change unit and apply factor
501	CALL var_chg_unit(tl_var(jk))
502
503	! forced min and max value
504	CALL var_limit_value(tl_var(jk))
505
506	! fill closed sea
507	IF( ll_fillclosed )THEN
508	ALLOCATE( il_mask(tl_var(jk)%t_dim(1)%i_len, &
509	& tl_var(jk)%t_dim(2)%i_len) )
510
511	! split domain in N sea subdomain
512	il_mask(:,:)=grid_split_domain(tl_var(jk))
513
514	! fill smallest domain
515	CALL grid_fill_small_dom( tl_var(jk), il_mask(:,:) )
516
517	DEALLOCATE( il_mask )
518	ENDIF
519
520	! filter
521	CALL filter_fill_value(tl_var(jk))
522
523	! check bathymetry
524	dl_minbat=MINVAL(tl_var(jk)%d_value(:,:,:,:))
525	IF( TRIM(tl_var(jk)%c_stdname) == 'bathymetry' .AND. &
526	& dl_minbat <= 0._dp )THEN
527	CALL logger_debug("CREATE BATHY: min value "//TRIM(fct_str(dl_minbat)))
528	CALL logger_fatal("CREATE BATHY: Bathymetry has value <= 0")
529	ENDIF
530
531	ENDDO
532
533	! create file
534	tl_fileout=file_init(TRIM(cn_fileout))
535
536	! add dimension
537	tl_dim(:)=var_max_dim(tl_var(:))
538
539	DO ji=1,ip_maxdim
540	IF( tl_dim(ji)%l_use ) CALL file_add_dim(tl_fileout, tl_dim(ji))
541	ENDDO
542
543	! add variables
544	IF( ALL( tl_dim(1:2)%l_use ) )THEN
545
546	! open mpp files
547	CALL iom_mpp_open(tl_coord1)
548
549	! add longitude
550	tl_lon=iom_mpp_read_var(tl_coord1,'longitude')
551	CALL file_add_var(tl_fileout, tl_lon)
552	CALL var_clean(tl_lon)
553
554	! add latitude
555	tl_lat=iom_mpp_read_var(tl_coord1,'latitude')
556	CALL file_add_var(tl_fileout, tl_lat)
557	CALL var_clean(tl_lat)
558
559	! close mpp files
560	CALL iom_mpp_close(tl_coord1)
561
562	ENDIF
563
564	IF( tl_dim(3)%l_use )THEN
565	! add depth
566	CALL file_add_var(tl_fileout, tl_depth)
567	CALL var_clean(tl_depth)
568	ENDIF
569
570	IF( tl_dim(4)%l_use )THEN
571	! add time
572	CALL file_add_var(tl_fileout, tl_time)
573	CALL var_clean(tl_time)
574	ENDIF
575
576	! add other variables
577	DO jk=tl_multi%i_nvar,1,-1
578	CALL file_add_var(tl_fileout, tl_var(jk))
579	CALL var_clean(tl_var(jk))
580	ENDDO
581	DEALLOCATE(tl_var)
582
583	! add some attribute
584	tl_att=att_init("Created_by","SIREN create_bathy")
585	CALL file_add_att(tl_fileout, tl_att)
586
587	cl_date=date_print(date_now())
588	tl_att=att_init("Creation_date",cl_date)
589	CALL file_add_att(tl_fileout, tl_att)
590
591	! add attribute periodicity
592	il_attid=0
593	IF( ASSOCIATED(tl_fileout%t_att) )THEN
594	il_attid=att_get_index(tl_fileout%t_att(:),'periodicity')
595	ENDIF
596	IF( tl_coord1%i_perio >= 0 .AND. il_attid == 0 )THEN
597	tl_att=att_init('periodicity',tl_coord1%i_perio)
598	CALL file_add_att(tl_fileout,tl_att)
599	ENDIF
600
601	! add attribute east west overlap
602	il_attid=0
603	IF( ASSOCIATED(tl_fileout%t_att) )THEN
604	il_attid=att_get_index(tl_fileout%t_att(:),'ew_overlap')
605	ENDIF
606	IF( tl_coord1%i_ew >= 0 .AND. il_attid == 0 )THEN
607	tl_att=att_init('ew_overlap',tl_coord1%i_ew)
608	CALL file_add_att(tl_fileout,tl_att)
609	ENDIF
610
611	! create file
612	CALL iom_create(tl_fileout)
613
614	! write file
615	CALL iom_write_file(tl_fileout)
616
617	! close file
618	CALL iom_close(tl_fileout)
619
620	! clean
621	CALL att_clean(tl_att)
622
623	CALL file_clean(tl_fileout)
624	CALL mpp_clean(tl_coord1)
625	CALL mpp_clean(tl_coord0)
626	CALL var_clean_extra()
627
628	! close log file
629	CALL logger_footer()
630	CALL logger_close()
631
632	CONTAINS
633	!-------------------------------------------------------------------
634	!> @brief
635	!> This function create variable, filled with matrix value
636	!>
637	!> @details
638	!> A variable is create with the same name that the input variable,
639	!> and with dimension of the coordinate file.<br/>
640	!> Then the variable array of value is split into equal subdomain.
641	!> Each subdomain is filled with the corresponding value of the matrix.
642	!>
643	!> @author J.Paul
644	!> @date November, 2013 - Initial Version
645	!>
646	!> @param[in] td_var variable structure
647	!> @param[in] td_coord coordinate file structure
648	!> @return variable structure
649	!-------------------------------------------------------------------
650	FUNCTION create_bathy_matrix(td_var, td_coord)
651	IMPLICIT NONE
652	! Argument
653	TYPE(TVAR), INTENT(IN) :: td_var
654	TYPE(TMPP), INTENT(IN) :: td_coord
655
656	! function
657	TYPE(TVAR) :: create_bathy_matrix
658
659	! local variable
660	INTEGER(i4) , DIMENSION(2,2) :: il_xghost
661	INTEGER(i4) , DIMENSION(3) :: il_dim
662	INTEGER(i4) , DIMENSION(3) :: il_size
663	INTEGER(i4) , DIMENSION(3) :: il_rest
664
665	INTEGER(i4) , DIMENSION(:) , ALLOCATABLE :: il_ishape
666	INTEGER(i4) , DIMENSION(:) , ALLOCATABLE :: il_jshape
667	INTEGER(i4) , DIMENSION(:) , ALLOCATABLE :: il_kshape
668
669	REAL(dp) , DIMENSION(:,:,:,:), ALLOCATABLE :: dl_value
670
671	TYPE(TVAR) :: tl_lon
672	TYPE(TDIM) , DIMENSION(ip_maxdim) :: tl_dim
673
674	TYPE(TMPP) :: tl_coord
675
676	! loop indices
677	INTEGER(i4) :: ji
678	INTEGER(i4) :: jj
679	INTEGER(i4) :: jk
680	!----------------------------------------------------------------
681
682	! copy structure
683	tl_coord=mpp_copy(td_coord)
684
685	! use only edge processor
686	CALL mpp_get_contour(tl_coord)
687
688	! open useful processor
689	CALL iom_mpp_open(tl_coord)
690
691	! read output grid
692	tl_lon=iom_mpp_read_var(tl_coord,'longitude')
693
694	! look for ghost cell
695	il_xghost(:,:)=grid_get_ghost( tl_coord )
696
697	! close processor
698	CALL iom_mpp_close(tl_coord)
699	! clean
700	CALL mpp_clean(tl_coord)
701
702	! remove ghost cell
703	CALL grid_del_ghost(tl_lon, il_xghost(:,:))
704
705	! write value on grid
706	! get matrix dimension
707	il_dim(:)=td_var%t_dim(1:3)%i_len
708	! output dimension
709	tl_dim(:)=dim_copy(tl_lon%t_dim(:))
710	! clean
711	CALL var_clean(tl_lon)
712
713	! split output domain in N subdomain depending of matrix dimension
714	il_size(:) = tl_dim(1:3)%i_len / il_dim(:)
715	il_rest(:) = MOD(tl_dim(1:3)%i_len, il_dim(:))
716
717	ALLOCATE( il_ishape(il_dim(1)+1) )
718	il_ishape(:)=0
719	DO ji=2,il_dim(1)+1
720	il_ishape(ji)=il_ishape(ji-1)+il_size(1)
721	ENDDO
722	! add rest to last cell
723	il_ishape(il_dim(1)+1)=il_ishape(il_dim(1)+1)+il_rest(1)
724
725	ALLOCATE( il_jshape(il_dim(2)+1) )
726	il_jshape(:)=0
727	DO jj=2,il_dim(2)+1
728	il_jshape(jj)=il_jshape(jj-1)+il_size(2)
729	ENDDO
730	! add rest to last cell
731	il_jshape(il_dim(2)+1)=il_jshape(il_dim(2)+1)+il_rest(2)
732
733	ALLOCATE( il_kshape(il_dim(3)+1) )
734	il_kshape(:)=0
735	DO jk=2,il_dim(3)+1
736	il_kshape(jk)=il_kshape(jk-1)+il_size(3)
737	ENDDO
738	! add rest to last cell
739	il_kshape(il_dim(3)+1)=il_kshape(il_dim(3)+1)+il_rest(3)
740
741	! write ouput array of value
742	ALLOCATE(dl_value( tl_dim(1)%i_len, &
743	& tl_dim(2)%i_len, &
744	& tl_dim(3)%i_len, &
745	& tl_dim(4)%i_len) )
746
747	dl_value(:,:,:,:)=0
748
749	DO jk=2,il_dim(3)+1
750	DO jj=2,il_dim(2)+1
751	DO ji=2,il_dim(1)+1
752
753	dl_value( 1+il_ishape(ji-1):il_ishape(ji), &
754	& 1+il_jshape(jj-1):il_jshape(jj), &
755	& 1+il_kshape(jk-1):il_kshape(jk), &
756	& 1 ) = td_var%d_value(ji-1,jj-1,jk-1,1)
757
758	ENDDO
759	ENDDO
760	ENDDO
761
762	! initialise variable with value
763	create_bathy_matrix=var_init(TRIM(td_var%c_name),dl_value(:,:,:,:))
764
765	DEALLOCATE(dl_value)
766
767	! add ghost cell
768	CALL grid_add_ghost(create_bathy_matrix, il_xghost(:,:))
769
770	! clean
771	CALL dim_clean(tl_dim(:))
772
773	END FUNCTION create_bathy_matrix
774	!-------------------------------------------------------------------
775	!> @brief
776	!> This function extract variable from file over coordinate domain and
777	!> return variable structure
778	!>
779	!> @author J.Paul
780	!> @date November, 2013 - Initial Version
781	!>
782	!> @param[in] td_var variable structure
783	!> @param[in] td_mpp mpp file structure
784	!> @param[in] td_coord coordinate file structure
785	!> @return variable structure
786	!-------------------------------------------------------------------
787	FUNCTION create_bathy_extract(td_var, td_mpp, &
788	& td_coord)
789	IMPLICIT NONE
790	! Argument
791	TYPE(TVAR), INTENT(IN) :: td_var
792	TYPE(TMPP), INTENT(IN) :: td_mpp
793	TYPE(TMPP), INTENT(IN) :: td_coord
794
795	! function
796	TYPE(TVAR) :: create_bathy_extract
797
798	! local variable
799	INTEGER(i4), DIMENSION(2,2) :: il_ind
800
801	INTEGER(i4) :: il_imin
802	INTEGER(i4) :: il_jmin
803	INTEGER(i4) :: il_imax
804	INTEGER(i4) :: il_jmax
805
806	TYPE(TMPP) :: tl_mpp
807
808	TYPE(TATT) :: tl_att
809
810	TYPE(TVAR) :: tl_var
811
812	TYPE(TDOM) :: tl_dom
813	! loop indices
814	!----------------------------------------------------------------
815
816	IF( .NOT. ASSOCIATED(td_mpp%t_proc) )THEN
817	CALL logger_error("CREATE BATHY EXTRACT: no processor associated "//&
818	& "to mpp "//TRIM(td_mpp%c_name))
819	ELSE
820
821	!init
822	tl_mpp=mpp_copy(td_mpp)
823
824	! compute file grid indices around coord grid
825	il_ind(:,:)=grid_get_coarse_index(tl_mpp, td_coord )
826
827	il_imin=il_ind(1,1) ; il_imax=il_ind(1,2)
828	il_jmin=il_ind(2,1) ; il_jmax=il_ind(2,2)
829
830	! check grid coincidence
831	CALL grid_check_coincidence( tl_mpp, td_coord, &
832	& il_imin, il_imax, &
833	& il_jmin, il_jmax, &
834	& (/1, 1, 1/) )
835
836	! compute domain
837	tl_dom=dom_init(tl_mpp, &
838	& il_imin, il_imax, &
839	& il_jmin, il_jmax)
840
841	! open mpp files over domain
842	CALL iom_dom_open(tl_mpp, tl_dom)
843
844	! read variable on domain
845	tl_var=iom_dom_read_var(tl_mpp,TRIM(td_var%c_name),tl_dom)
846
847	! close mpp file
848	CALL iom_dom_close(tl_mpp)
849
850	! add ghost cell
851	CALL grid_add_ghost(tl_var,tl_dom%i_ghost(:,:))
852
853	! check result
854	IF( ANY( tl_var%t_dim(:)%l_use .AND. &
855	& tl_var%t_dim(:)%i_len /= td_coord%t_dim(:)%i_len) )THEN
856	CALL logger_debug("CREATE BATHY EXTRACT: "//&
857	& "dimensoin of variable "//TRIM(td_var%c_name)//" "//&
858	& TRIM(fct_str(tl_var%t_dim(1)%i_len))//","//&
859	& TRIM(fct_str(tl_var%t_dim(2)%i_len))//","//&
860	& TRIM(fct_str(tl_var%t_dim(3)%i_len))//","//&
861	& TRIM(fct_str(tl_var%t_dim(4)%i_len)) )
862	CALL logger_debug("CREATE BATHY EXTRACT: "//&
863	& "dimensoin of coordinate file "//&
864	& TRIM(fct_str(td_coord%t_dim(1)%i_len))//","//&
865	& TRIM(fct_str(td_coord%t_dim(2)%i_len))//","//&
866	& TRIM(fct_str(td_coord%t_dim(3)%i_len))//","//&
867	& TRIM(fct_str(td_coord%t_dim(4)%i_len)) )
868	CALL logger_fatal("CREATE BATHY EXTRACT: "//&
869	& "dimensoin of extracted "//&
870	& "variable and coordinate file dimension differ")
871	ENDIF
872
873	! add attribute to variable
874	tl_att=att_init('src_file',TRIM(fct_basename(tl_mpp%c_name)))
875	CALL var_move_att(tl_var, tl_att)
876
877	tl_att=att_init('src_i_indices',(/tl_dom%i_imin, tl_dom%i_imax/))
878	CALL var_move_att(tl_var, tl_att)
879
880	tl_att=att_init('src_j_indices',(/tl_dom%i_jmin, tl_dom%i_jmax/))
881	CALL var_move_att(tl_var, tl_att)
882
883	! save result
884	create_bathy_extract=var_copy(tl_var)
885
886	! clean structure
887	CALL att_clean(tl_att)
888	CALL var_clean(tl_var)
889	CALL mpp_clean(tl_mpp)
890	ENDIF
891
892	END FUNCTION create_bathy_extract
893	!-------------------------------------------------------------------
894	!> @brief
895	!> This function get coarse grid variable, interpolate variable, and return
896	!> variable structure over fine grid
897	!>
898	!> @author J.Paul
899	!> @date November, 2013 - Initial Version
900	!>
901	!> @param[in] td_var variable structure
902	!> @param[in] td_mpp mpp file structure
903	!> @param[in] id_imin i-direction lower left corner indice
904	!> @param[in] id_imax i-direction upper right corner indice
905	!> @param[in] id_jmin j-direction lower left corner indice
906	!> @param[in] id_jmax j-direction upper right corner indice
907	!> @param[in] id_offset offset between fine grid and coarse grid
908	!> @param[in] id_rho array of refinement factor
909	!> @return variable structure
910	!-------------------------------------------------------------------
911	FUNCTION create_bathy_get_var(td_var, td_mpp, &
912	& id_imin, id_jmin, &
913	& id_imax, id_jmax, &
914	& id_offset, &
915	& id_rho )
916	IMPLICIT NONE
917	! Argument
918	TYPE(TVAR) , INTENT(IN) :: td_var
919	TYPE(TMPP) , INTENT(IN) :: td_mpp
920	INTEGER(i4) , INTENT(IN) :: id_imin
921	INTEGER(i4) , INTENT(IN) :: id_imax
922	INTEGER(i4) , INTENT(IN) :: id_jmin
923	INTEGER(i4) , INTENT(IN) :: id_jmax
924	INTEGER(i4), DIMENSION(:,:), INTENT(IN) :: id_offset
925	INTEGER(i4), DIMENSION(:) , INTENT(IN) :: id_rho
926
927	! function
928	TYPE(TVAR) :: create_bathy_get_var
929
930	! local variable
931	TYPE(TMPP) :: tl_mpp
932	TYPE(TATT) :: tl_att
933	TYPE(TVAR) :: tl_var
934	TYPE(TDOM) :: tl_dom
935
936	INTEGER(i4) :: il_size
937	INTEGER(i4), DIMENSION(:), ALLOCATABLE :: il_rho
938
939	! loop indices
940	!----------------------------------------------------------------
941	IF( ANY(SHAPE(id_offset(:,:)) /= 2) )THEN
942	CALL logger_error("CREATE BATHY GET VAR: invalid dimension of "//&
943	& "offset array")
944	ENDIF
945
946	! copy structure
947	tl_mpp=mpp_copy(td_mpp)
948
949	!- compute domain
950	tl_dom=dom_init(tl_mpp, &
951	& id_imin, id_imax, &
952	& id_jmin, id_jmax)
953
954	!- add extra band (if possible) to compute interpolation
955	CALL dom_add_extra(tl_dom)
956
957	!- open mpp files over domain
958	CALL iom_dom_open(tl_mpp, tl_dom)
959
960	!- read variable value on domain
961	tl_var=iom_dom_read_var(tl_mpp,TRIM(td_var%c_name),tl_dom)
962
963	!- close mpp files
964	CALL iom_dom_close(tl_mpp)
965
966	il_size=SIZE(id_rho(:))
967	ALLOCATE( il_rho(il_size) )
968	il_rho(:)=id_rho(:)
969
970	!- interpolate variable
971	CALL create_bathy_interp(tl_var, il_rho(:), id_offset(:,:))
972
973	!- remove extraband added to domain
974	CALL dom_del_extra( tl_var, tl_dom, il_rho(:) )
975
976	CALL dom_clean_extra( tl_dom )
977
978	!- add ghost cell
979	CALL grid_add_ghost(tl_var,tl_dom%i_ghost(:,:))
980
981	!- add attribute to variable
982	tl_att=att_init('src_file',TRIM(fct_basename(tl_mpp%c_name)))
983	CALL var_move_att(tl_var, tl_att)
984
985	tl_att=att_init('src_i_indices',(/tl_dom%i_imin, tl_dom%i_imax/))
986	CALL var_move_att(tl_var, tl_att)
987
988	tl_att=att_init('src_j_indices',(/tl_dom%i_jmin, tl_dom%i_jmax/))
989	CALL var_move_att(tl_var, tl_att)
990
991	IF( .NOT. ALL(id_rho(:)==1) )THEN
992	tl_att=att_init("refinment_factor",(/id_rho(jp_I),id_rho(jp_J)/))
993	CALL var_move_att(tl_var, tl_att)
994	ENDIF
995
996	DEALLOCATE( il_rho )
997
998	!- save result
999	create_bathy_get_var=var_copy(tl_var)
1000
1001	!- clean structure
1002	CALL att_clean(tl_att)
1003	CALL var_clean(tl_var)
1004	CALL mpp_clean(tl_mpp)
1005
1006	END FUNCTION create_bathy_get_var
1007	!-------------------------------------------------------------------
1008	!> @brief
1009	!> This subroutine interpolate variable
1010	!>
1011	!> @author J.Paul
1012	!> @date November, 2013 - Initial Version
1013	!>
1014	!> @param[inout] td_var variable structure
1015	!> @param[in] id_rho array of refinment factor
1016	!> @param[in] id_offset array of offset between fine and coarse grid
1017	!> @param[in] id_iext i-direction size of extra bands (default=im_minext)
1018	!> @param[in] id_jext j-direction size of extra bands (default=im_minext)
1019	!-------------------------------------------------------------------
1020	SUBROUTINE create_bathy_interp( td_var, &
1021	& id_rho, &
1022	& id_offset, &
1023	& id_iext, id_jext)
1024
1025	IMPLICIT NONE
1026
1027	! Argument
1028	TYPE(TVAR) , INTENT(INOUT) :: td_var
1029	INTEGER(i4), DIMENSION(:) , INTENT(IN ) :: id_rho
1030	INTEGER(i4), DIMENSION(:,:), INTENT(IN ) :: id_offset
1031	INTEGER(i4), INTENT(IN ), OPTIONAL :: id_iext
1032	INTEGER(i4), INTENT(IN ), OPTIONAL :: id_jext
1033
1034	! local variable
1035	TYPE(TVAR) :: tl_mask
1036
1037	INTEGER(i1), DIMENSION(:,:,:,:), ALLOCATABLE :: bl_mask
1038
1039	INTEGER(i4) :: il_iext
1040	INTEGER(i4) :: il_jext
1041
1042	! loop indices
1043	!----------------------------------------------------------------
1044
1045	!WARNING: two extrabands are required for cubic interpolation
1046	il_iext=3
1047	IF( PRESENT(id_iext) ) il_iext=id_iext
1048
1049	il_jext=3
1050	IF( PRESENT(id_jext) ) il_jext=id_jext
1051
1052	IF( il_iext < 2 .AND. td_var%c_interp(1) == 'cubic' )THEN
1053	CALL logger_warn("CREATE BATHY INTERP: at least extrapolation "//&
1054	& "on two points are required with cubic interpolation ")
1055	il_iext=2
1056	ENDIF
1057
1058	IF( il_jext < 2 .AND. td_var%c_interp(1) == 'cubic' )THEN
1059	CALL logger_warn("CREATE BATHY INTERP: at least extrapolation "//&
1060	& "on two points are required with cubic interpolation ")
1061	il_jext=2
1062	ENDIF
1063
1064	! work on mask
1065	! create mask
1066	ALLOCATE(bl_mask(td_var%t_dim(1)%i_len, &
1067	& td_var%t_dim(2)%i_len, &
1068	& td_var%t_dim(3)%i_len, &
1069	& td_var%t_dim(4)%i_len) )
1070
1071	bl_mask(:,:,:,:)=1
1072	WHERE(td_var%d_value(:,:,:,:)==td_var%d_fill) bl_mask(:,:,:,:)=0
1073
1074	SELECT CASE(TRIM(td_var%c_point))
1075	CASE DEFAULT ! 'T'
1076	tl_mask=var_init('tmask', bl_mask(:,:,:,:), td_dim=td_var%t_dim(:), &
1077	& id_ew=td_var%i_ew )
1078	CASE('U','V','F')
1079	CALL logger_fatal("CREATE BATHY INTERP: can not computed "//&
1080	& "interpolation on "//TRIM(td_var%c_point)//&
1081	& " grid point (variable "//TRIM(td_var%c_name)//&
1082	& "). check namelist.")
1083	END SELECT
1084
1085	DEALLOCATE(bl_mask)
1086
1087	! interpolate mask
1088	CALL interp_fill_value( tl_mask, id_rho(:), &
1089	& id_offset=id_offset(:,:) )
1090
1091	! work on variable
1092	! add extraband
1093	CALL extrap_add_extrabands(td_var, il_iext, il_jext)
1094
1095	! extrapolate variable
1096	CALL extrap_fill_value( td_var )
1097
1098	! interpolate Bathymetry
1099	CALL interp_fill_value( td_var, id_rho(:), &
1100	& id_offset=id_offset(:,:) )
1101
1102	! remove extraband
1103	CALL extrap_del_extrabands(td_var, il_iextid_rho(jp_I), il_jextid_rho(jp_J))
1104
1105	! keep original mask
1106	WHERE( tl_mask%d_value(:,:,:,:) == 0 )
1107	td_var%d_value(:,:,:,:)=td_var%d_fill
1108	END WHERE
1109
1110	! clean variable structure
1111	CALL var_clean(tl_mask)
1112
1113	END SUBROUTINE create_bathy_interp
1114	!-------------------------------------------------------------------
1115	!> @brief
1116	!> This subroutine get depth variable value in an open mpp structure
1117	!> and check if agree with already input depth variable.
1118	!>
1119	!> @details
1120	!>
1121	!> @author J.Paul
1122	!> @date January, 2016 - Initial Version
1123	!>
1124	!> @param[in] td_mpp mpp structure
1125	!> @param[inout] td_depth depth variable structure
1126	!-------------------------------------------------------------------
1127	SUBROUTINE create_bathy_check_depth( td_mpp, td_depth )
1128
1129	IMPLICIT NONE
1130
1131	! Argument
1132	TYPE(TMPP) , INTENT(IN ) :: td_mpp
1133	TYPE(TVAR) , INTENT(INOUT) :: td_depth
1134
1135	! local variable
1136	INTEGER(i4) :: il_varid
1137	TYPE(TVAR) :: tl_depth
1138	! loop indices
1139	!----------------------------------------------------------------
1140
1141	! get or check depth value
1142	IF( td_mpp%t_proc(1)%i_depthid /= 0 )THEN
1143
1144	il_varid=td_mpp%t_proc(1)%i_depthid
1145	IF( ASSOCIATED(td_depth%d_value) )THEN
1146
1147	tl_depth=iom_mpp_read_var(td_mpp, il_varid)
1148
1149	IF( ANY( td_depth%d_value(:,:,:,:) /= &
1150	& tl_depth%d_value(:,:,:,:) ) )THEN
1151
1152	CALL logger_warn("CREATE BATHY: depth value from "//&
1153	& TRIM(td_mpp%c_name)//" not conform "//&
1154	& " to those from former file(s).")
1155
1156	ENDIF
1157	CALL var_clean(tl_depth)
1158
1159	ELSE
1160	td_depth=iom_mpp_read_var(td_mpp,il_varid)
1161	ENDIF
1162
1163	ENDIF
1164
1165	END SUBROUTINE create_bathy_check_depth
1166	!-------------------------------------------------------------------
1167	!> @brief
1168	!> This subroutine get date and time in an open mpp structure
1169	!> and check if agree with date and time already read.
1170	!>
1171	!> @details
1172	!>
1173	!> @author J.Paul
1174	!> @date January, 2016 - Initial Version
1175	!>
1176	!> @param[in] td_mpp mpp structure
1177	!> @param[inout] td_time time variable structure
1178	!-------------------------------------------------------------------
1179	SUBROUTINE create_bathy_check_time( td_mpp, td_time )
1180
1181	IMPLICIT NONE
1182
1183	! Argument
1184	TYPE(TMPP), INTENT(IN ) :: td_mpp
1185	TYPE(TVAR), INTENT(INOUT) :: td_time
1186
1187	! local variable
1188	INTEGER(i4) :: il_varid
1189	TYPE(TVAR) :: tl_time
1190
1191	TYPE(TDATE) :: tl_date1
1192	TYPE(TDATE) :: tl_date2
1193	! loop indices
1194	!----------------------------------------------------------------
1195
1196	! get or check depth value
1197	IF( td_mpp%t_proc(1)%i_timeid /= 0 )THEN
1198
1199	il_varid=td_mpp%t_proc(1)%i_timeid
1200	IF( ASSOCIATED(td_time%d_value) )THEN
1201
1202	tl_time=iom_mpp_read_var(td_mpp, il_varid)
1203
1204	tl_date1=var_to_date(td_time)
1205	tl_date2=var_to_date(tl_time)
1206	IF( tl_date1 - tl_date2 /= 0 )THEN
1207
1208	CALL logger_warn("CREATE BATHY: date from "//&
1209	& TRIM(td_mpp%c_name)//" not conform "//&
1210	& " to those from former file(s).")
1211
1212	ENDIF
1213	CALL var_clean(tl_time)
1214
1215	ELSE
1216	td_time=iom_mpp_read_var(td_mpp,il_varid)
1217	ENDIF
1218
1219	ENDIF
1220
1221	END SUBROUTINE create_bathy_check_time
1222	END PROGRAM create_bathy

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

Download in other formats: