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 branches/2015/nemo_v3_6_STABLE/NEMOGCM/TOOLS/SIREN/src – NEMO

Context Navigation

source: branches/2015/nemo_v3_6_STABLE/NEMOGCM/TOOLS/SIREN/src/create_bathy.f90 @ 6392

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

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

Download in other formats: