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/UKMO/r5518_INGV1_WAVE-coupling/NEMOGCM/TOOLS/SIREN/src – NEMO

Context Navigation

source: branches/UKMO/r5518_INGV1_WAVE-coupling/NEMOGCM/TOOLS/SIREN/src/create_bathy.f90 @ 7152

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

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

Download in other formats: