Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Blame
Revision Log

mod_psmile_io.F90 @ 1677

Last change on this file since 1677 was 1677, checked in by aclsce, 12 years ago
Imported oasis3 (tag ipslcm5a) from cvs server to svn server (igcmg project).
File size: 75.1 KB

Line
1	#undef __DEBUG
2
3	!#define __VERBOSE
4	#define __READ_BY_LOOKUP
5	!----------------------------------------------------------------------
6	! BOP
7	!
8	! !MODULE: mod_psmile_io
9	! !REMARKS: Initially programed by Reiner Vogelsang, SGI GmbH, and integrated
10	! into OASIS 3.0 and tested by Arnaud Caubel, CERFACS
11	! Bugs should be reported to reiner@sgi.com or Sophie.Valcke@cerfacs.fr.
12	! !REVISION HISTORY:
13	! 2003.01.31 Reiner Vogelsang Finalized initial version
14	! 2003.02.28 Arnaud Caubel Tests and correction of the READ_BY_LOOKUP
15	! 2003.04.24 Reiner Vogelsang Corrections with regard to CF naming convention.
16	! 2003.04.28 Reiner Vogelsang Added ProTex documentation header and
17	! static CVS versions string.
18	! 2004.03.08 Sophie Valcke Added the writting of bounds,lons and lats.
19	! 2004.04.07 Reiner Vogelsang Correction: bounds are written with
20	! number of corners as leading dimension
21	! 2004.19.06 Reiner Vogelsang Considering cyclic boundaries for longitudes.
22	! which lead to differences between corners
23	! and bounds by 360 degrees
24	! 2004.26.08 Reiner Vogelsang Correction of allocation scheme for arrays
25	! hosting the longitudes.
26	!
27	! !PUBLIC MEMBER FUNCTIONS:
28	!
29	! subroutine psmile_io_init_comp(id_error)
30	! This subroutine finds a local communicator for the MPI tasks
31	! of a model component participating in coupling.
32	!
33	! subroutine psmile_def_domains(id_error)
34	! This subroutine initializes the MPP_IO domains.
35	!
36	! subroutine psmile_def_files(id_error)
37	! This subroutine open files acccording to the kewords
38	!
39	! subroutine psmile_def_metadata(id_error)
40	! This subroutine writes CF metadata or finds a variable according
41	! to CF metadata information.
42	!
43	! subroutine psmile_close_files(id_error)
44	! This subroutine closes all opened psmile_io files
45	!
46	! subroutine psmile_io_cleanup(id_error)
47	! This subroutine deallocates all PSMILe I/O related data
48	! and exits MPP_IO
49	!
50	! subroutine psmile_read_8(id_port_id,rd_field,id_newtime)
51	! This subroutine reads REAL(KIND=4) data.
52	!
53	! subroutine psmile_read_4(id_port_id,rd_field,id_newtime)
54	! This subroutine reads REAL(KIND=8) data.
55	!
56	! subroutine psmile_write_4(id_port_id,rd_field,id_newtime)
57	! This subroutine writes REAL(KIND=4) data.
58	!
59	! subroutine psmile_write_8(id_port_id,rd_field,id_newtime)
60	! This subroutine writes REAL(KIND=8) data.
61	!
62	! subroutine indexi(n,arr,indx)
63	! This subroutine sorts by indexing.
64	!
65	! subroutine combine_with_date(cd_in,cd_mode,id_initial_date,cd_on)
66	! This subroutines appends the date in ISO format to a strings.
67	!
68
69	module mod_psmile_io
70
71	! !USES:
72	use mod_kinds_model
73	use mpp_mod
74	use mpp_io_mod
75	use mpp_domains_mod
76	! !PUBLIC TYPES:
77	!Highest unit number permittable by the OS.
78	integer(kind=ip_intwp_p),parameter::PSMILE_MAX_UNIT=299
79	!Highest unit number reserved by the OS.
80	integer(kind=ip_intwp_p),parameter::PSMILE_MAX_RESERVED_UNIT=104
81	!Lowest unit number reserved by the OS.
82	integer(kind=ip_intwp_p),parameter::PSMILE_MIN_RESERVED_UNIT=100
83	integer(kind=ip_intwp_p),parameter::PSMILE_MAX_CF_TABLE_ENTRIES=1024
84	integer(kind=ip_intwp_p),allocatable::nbr_corners(:)
85
86	integer(kind=ip_intwp_p)::ig_mpp_io_comm,ig_color,ig_mpp_io_comm_size,ig_mpp_io_comm_rank
87	integer(kind=ip_intwp_p)::ig_stackmax=1310720,ig_stackmaxd=655360
88	integer(kind=ip_intwp_p)::ig_layout(2)
89	!Global arrays to specify IO domains
90	type(domain2D),allocatable::domain(:)
91	type(domain1D),allocatable::xdom(:),ydom(:)
92	logical,allocatable::lg_is_odd_distribution(:)
93	!Global arrarys for opening files
94	integer(kind=ip_intwp_p),allocatable::ig_units(:)
95	integer(kind=ip_intwp_p),allocatable::ig_action(:)
96	integer(kind=ip_intwp_p),allocatable::ig_form(:)
97	integer(kind=ip_intwp_p),allocatable::ig_threading(:)
98	integer(kind=ip_intwp_p),allocatable::ig_fileset(:)
99	character(len=128),allocatable::cg_my_input_file(:)
100	character(len=128),allocatable::cg_output_file(:)
101
102
103	!Global arrays for netcdf metadata informations
104	integer(kind=ip_intwp_p),allocatable::ig_var_index(:)
105	integer(kind=ip_intwp_p),allocatable::ig_ntimes(:)
106	!RV: 22.01.2003
107	TYPE time_stamps
108	! real(kind=ip_double_p),pointer::rl_times(:)
109	!AC:31/01
110	real(kind=ip_double_p),dimension(:),pointer :: rl_times
111	END TYPE time_stamps
112	! !PUBLIC DATA MEMBERS:
113	type(time_stamps),allocatable :: dg_times(:)
114
115	type(fieldtype),allocatable::ig_vars(:),field(:),latij(:),lonij(:)
116	type(fieldtype),allocatable::crnlatij(:), crnlonij(:)
117	type(axistype),allocatable::x_axis(:),y_axis(:),z_axis(:),t_axis(:)
118	type(axistype),allocatable::c_axis(:)
119	! These arrays are supposed to be part of a different module!
120	! However, no information is available in Oasis 3.0. So, let them live here
121	! for a while.
122	!RV
123	character(len=64),allocatable::cports_lgname(:)
124	character(len=64),allocatable::cports_units(:)
125	character(len=64),allocatable::cfldlab(:),cfunitslab(:)
126	!RV
127	! !DESCRIPTION:
128	! This module contains the global data objects which hold informations for
129	! performing I/O on coupled fields and are the glue between routines defined
130	! within this file, mod_psmile_io.F90.
131	! !EOP
132	!-------------------------------------------------------------------------------
133	!$Id: mod_psmile_io.F90,v 1.2 2009/01/15 15:10:14 adm Exp $
134	!-------------------------------------------------------------------------------
135	end module mod_psmile_io
136
137
138
139	subroutine psmile_io_init_comp(id_error)
140	!--------------------------------------------------------------------------
141
142	!Routine must be called in prism_init_comp.
143	!The assumption is that each process of a model calls prism_init_comp.
144	use mod_kinds_model
145	use mod_prism_proto
146	use mod_comprism_proto
147	use mod_psmile_io
148
149	implicit none
150
151	#include <mpif.h>
152
153
154
155	!Arguments:
156	!Error return code
157	integer(kind=ip_intwp_p),intent(out)::id_error
158
159	!Local Variables:
160	!Don't touch the following character variable!
161	character(len=80),save::cl_version_string= &
162	'$Id: mod_psmile_io.F90,v 1.2 2009/01/15 15:10:14 adm Exp $'
163	!
164	!Communicator which covers the whole model.
165	integer(kind=ip_intwp_p)::il_local_comm
166
167	integer(kind=ip_intwp_p)::il_i,il_key
168	integer(kind=ip_intwp_p),allocatable::il_cpl_ranks(:)
169	integer(kind=ip_intwp_p)::ig_local_comm_grp,il_ncplprocs,ig_local_comm_size,il_mynum
170
171	integer(kind=ip_intwp_p)::il_file_unit,il_max_entry_id,il_no_of_entries
172	integer(kind=ip_intwp_p)::il_pos
173	character(len=64) :: cl_cfname,cl_cfunit
174	logical :: ll_exist
175
176
177	#ifdef __VERBOSE
178	write(nulprt,*) 'Start - - psmile_io_init_comp'
179	#endif
180	id_error=CLIM_OK
181	il_local_comm=ig_local_comm
182
183	!Create a communicator which covers all processes of a model invoLved with
184	!coupling.
185	il_ncplprocs=kbcplproc(ig_mynummod)
186	allocate(il_cpl_ranks(1:il_ncplprocs))
187
188	do il_i=0,il_ncplprocs-1
189	il_cpl_ranks(il_i+1)=il_i
190	enddo
191
192	il_key=1
193	call MPI_Comm_rank(il_local_comm,il_mynum,id_error)
194
195	!
196	!Create a color by finding the cpl ranks in the model communicator.
197	!
198	if(ANY(il_cpl_ranks.eq.il_mynum)) then
199	ig_color=1000*ig_mynummod
200	else
201	ig_color=1
202	endif
203	#ifdef __VERBOSE
204	write(nulprt,*)'psmile_io_init_comp: rank ',il_mynum,' color ',ig_color
205	#endif
206
207	!
208	!Do the splitting.
209	!
210	#ifdef __VERBOSE
211	write(nulprt,*) 'psmile_io_init_comp: Communicator splitting'
212	#endif
213	call MPI_Comm_split(il_local_comm,ig_color,il_key,ig_mpp_io_comm,id_error)
214	call MPI_Comm_rank(ig_mpp_io_comm,ig_mpp_io_comm_rank,id_error)
215	call MPI_Comm_size(ig_mpp_io_comm,ig_mpp_io_comm_size,id_error)
216
217	!
218	!Initialize FMS mpp_init,mpp_io,mpp_domains
219	!
220	if (ig_color .eq. 1000*ig_mynummod ) then
221	#ifdef __VERBOSE
222	write(nulprt,*) 'psmile_io_init_comp: MPP_INIT : Model',ig_mynummod &
223	,ig_mpp_io_comm_rank
224	#endif
225	call mpp_init(mpp_comm=ig_mpp_io_comm,logfile=trim(cnaprt))
226	call mpp_domains_set_stack_size(ig_stackmaxd)
227	!rv
228	!rv Oasis 3 uses for the tracefiles the I/O units >= 1000. Since
229	!rv mpp_io stores informations regarding to file specs in an array of
230	!rv structs that array has to be allocated apropiately.
231	!rv maxresunit specifies that the top maxresunits are excluded from the list
232	!rv of the files to be closed on mpp_io_exit.
233	!rv
234	call mpp_io_init(maxunit=nulprt,maxresunit=25)
235	call mpp_set_stack_size(ig_stackmax)
236	endif
237
238	!After this point I assume that I/O routines are called by processes
239	!involved with coupling only.
240
241	!RV At this point I am reading the CF name table definitions.
242
243	inquire(file='cf_name_table.txt',exist=ll_exist)
244	if(ll_exist) then
245	#ifdef __VERBOSE
246	write(nulprt,*) 'psmile_io_init_comp: Reading CF name table!'
247	call flush(nulprt)
248	#endif
249
250	il_file_unit=nulprt+1
251	open(file='cf_name_table.txt' &
252	,unit=il_file_unit &
253	,form='formatted' &
254	,status='old')
255
256	read(unit=il_file_unit,fmt=*,iostat=id_error)
257	read(unit=il_file_unit,fmt=*,iostat=id_error) &
258	il_max_entry_id, il_no_of_entries
259
260	if(id_error.ne.0) then
261	write(nulprt,*) &
262	'psmile_io_init_comp:cf_name_table.txt:' &
263	,' Reading of first record failed!'
264	call flush(nulprt)
265	call MPI_Abort(mpi_comm,0,mpi_error)
266	endif
267
268	if(il_max_entry_id.gt.0.and. &
269	il_max_entry_id.le.PSMILE_MAX_CF_TABLE_ENTRIES) then
270
271	allocate(cfldlab(1:il_max_entry_id),STAT=id_error)
272	if(id_error.ne.0) then
273	write(nulprt,*) &
274	'psmile_io_init_comp: Allocation of cfldlab failed!'
275	call flush(nulprt)
276	call MPI_Abort(mpi_comm,0,mpi_error)
277	endif
278
279	allocate(cfunitslab(1:il_max_entry_id),STAT=id_error)
280	if(id_error.ne.0) then
281	write(nulprt,*) &
282	'psmile_io_init_comp: Allocation of cfunitslab failed!'
283	call flush(nulprt)
284	call MPI_Abort(mpi_comm,0,mpi_error)
285	endif
286
287	else
288
289	write(nulprt,*) &
290	'psmile_io_init_comp:cf_name_table.txt:' &
291	,'Your max. numlab is out of range !'
292	call flush(nulprt)
293	call MPI_Abort(mpi_comm,0,mpi_error)
294
295	endif
296
297
298	read(unit=il_file_unit,fmt=*,iostat=id_error)
299	do il_i=1,il_no_of_entries
300
301	read(unit=il_file_unit,fmt=*,iostat=id_error) &
302	il_pos,cl_cfname,cl_cfunit
303
304	if(id_error.eq.0) then
305	if(il_pos.le.il_max_entry_id) then
306	cfldlab(il_pos)=trim(cl_cfname)
307	cfunitslab(il_pos)=trim(cl_cfunit)
308	else
309	write(nulprt,*) &
310	'psmile_io_init_comp:cf_name_table.txt:' &
311	,'Record ',il_i,': numlab =',il_pos,' out of range!'
312	call flush(nulprt)
313	call MPI_Abort(mpi_comm,0,mpi_error)
314	endif
315	else
316	write(nulprt,*) &
317	'psmile_io_init_comp:cf_name_table.txt:' &
318	,'Reading record ',il_i,' failed!'
319	call flush(nulprt)
320	call MPI_Abort(mpi_comm,0,mpi_error)
321	endif
322
323	enddo
324
325	endif
326
327
328
329	!RV
330	deallocate(il_cpl_ranks)
331	#ifdef __VERBOSE
332	write(nulprt,*) 'End - - psmile_io_init_comp'
333	#endif
334
335	return
336	end
337
338	subroutine psmile_def_domains(id_error)
339	!--------------------------------------------------------------------------
340	!
341	! Routine defines the IO domain partitioning
342	! Called at the end of prism_enddef.
343	!
344	use mod_kinds_model
345	use mod_prism_proto
346	use mod_comprism_proto
347	use mod_psmile_io
348	use mod_psmile_io_interfaces, only:indexi
349	implicit none
350
351	#include <mpif.h>
352
353	integer(kind=ip_intwp_p),intent(out)::id_error
354
355	integer(kind=ip_intwp_p)::il_port,il_i,il_j,il_ldx
356	integer(kind=ip_intwp_p),allocatable::il_gextents(:,:,:)
357	integer(kind=ip_intwp_p),allocatable::il_offsets(:,:,:)
358	integer(kind=ip_intwp_p),allocatable::il_extentsx(:)
359	integer(kind=ip_intwp_p),allocatable::il_extentsy(:)
360	integer(kind=ip_intwp_p),allocatable::il_offsetx(:)
361	integer(kind=ip_intwp_p),allocatable::il_offsety(:)
362	logical,allocatable::ll_mask(:,:)
363	integer(kind=ip_intwp_p),allocatable::il_pelist(:)
364	integer(kind=ip_intwp_p)::il_ngx,il_ngy
365	integer(kind=ip_intwp_p)::il_recvcount=1,il_sum
366
367	!Assumption: Processes beloging to the same model have the same
368	! strategy per port. Each process involved in the
369	! coupling has the same number of ports
370
371	#ifdef __VERBOSE
372	write(nulprt,*)'Start - - psmile_def_domains'
373	#endif
374
375	id_error=CLIM_ok
376
377	!Global arrays
378	allocate(domain(1:nports))
379	allocate(xdom(1:nports))
380	allocate(ydom(1:nports))
381	allocate(lg_is_odd_distribution(1:nports))
382
383	!Local arrays
384	!il_gextents(:,1,:) carries extents in x, il_gextents(:,2,:) extents in y
385	allocate(il_gextents(0:ig_mpp_io_comm_size-1,1:2,1:nports))
386	allocate(il_offsets(1:ig_mpp_io_comm_size,1:2,1:nports))
387	allocate(il_pelist(0:ig_mpp_io_comm_size-1))
388	allocate(il_extentsx(0:ig_mpp_io_comm_size-1))
389	allocate(il_extentsy(0:ig_mpp_io_comm_size-1))
390	allocate(il_offsetx(0:ig_mpp_io_comm_size-1))
391	allocate(il_offsety(0:ig_mpp_io_comm_size-1))
392	allocate(ll_mask(0:ig_mpp_io_comm_size-1,0:ig_mpp_io_comm_size-1))
393
394	if(ig_color/ig_mynummod.ne.1000) then
395	WRITE(nulprt,*)'psmile_def_domains: Called by proc which is not' &
396	,' involved with coupling!'
397	call flush(nulprt)
398	call MPI_ABORT(mpi_comm,0,mpi_err)
399	endif
400
401	do il_port=1,nports
402	#ifdef __VERBOSE
403	WRITE(nulprt,*)'ig_def_state ',ig_def_state(il_port)
404	#endif
405	IF (ig_def_state(il_port) .eq. ip_ignout .or. &
406	ig_def_state(il_port) .eq. ip_expout .or. &
407	ig_def_state(il_port) .eq. ip_output .or. &
408	ig_def_state(il_port) .eq. ip_input ) THEN
409	#ifdef __VERBOSE
410	WRITE(nulprt,*)'psmile_def_domains: port: ',il_port
411	call flush(nulprt)
412	#endif
413
414	if( mydist(CLIM_Strategy,il_port) .eq. CLIM_Serial) then
415
416	#ifdef __VERBOSE
417	WRITE(nulprt,*)'psmile_def_domains: port: ',il_port,' CLIM_Serial'
418	#endif
419
420	if(ig_mpp_io_comm_size.eq.1)then
421	#ifdef __VERBOSE
422	WRITE(nulprt,*)'psmile_def_domains: port: ',il_port &
423	,'(/1,mydist( CLIM_Segments+2,il_port),1,1/)' &
424	,1,mydist(CLIM_Segments+2,il_port),1,1
425	! call flush(nulprt)
426	#endif
427	ig_layout(1:2)=1
428	call mpp_define_domains((/1,mydist( CLIM_Segments+2,il_port),1,1/)&
429	, ig_layout,domain(il_port))
430	call mpp_get_domain_components( domain(il_port), xdom(il_port),&
431	ydom(il_port) )
432	else
433	write(nulprt, * ) 'psmile_def_domains: no. of procs for I/O ',&
434	'> 1 and CLIM_Serial does not make sense! '
435	call MPI_ABORT(mpi_comm,0,mpi_err)
436
437	endif
438
439	else if( mydist(CLIM_Strategy,il_port) .eq. CLIM_Apple) then
440
441	#ifdef __VERBOSE
442	WRITE(nulprt,*)'psmile_def_domains: port: ',il_port,' CLIM_Apple',mydist(:,il_port)
443	#endif
444
445	!Assuming 1D decomposition. At the this point this is the only information
446	!which I have considering CLIM_Apple!
447
448	ig_layout(1)=ig_mpp_io_comm_size
449	ig_layout(2)=1
450
451	!Build a list of the local extents which exists on the procs involved with
452	!coupling
453	call MPI_Allgather(mydist(CLIM_Length+1,il_port),1,MPI_INTEGER,&
454	il_gextents(0,1,il_port),il_recvcount,MPI_INTEGER,&
455	ig_mpp_io_comm,id_error)
456	call MPI_Allgather(mydist(CLIM_Offset+1,il_port),1,MPI_INTEGER,&
457	il_offsets(1,1,il_port),il_recvcount,MPI_INTEGER,&
458	ig_mpp_io_comm,id_error)
459
460	!Sort the extents and the PE-list according to the offsets
461
462	call indexi(ig_mpp_io_comm_size &
463	,il_offsets(1:ig_mpp_io_comm_size,1,il_port),il_pelist)
464
465	il_pelist(:)=il_pelist(:)-1
466	il_extentsx(0:ig_mpp_io_comm_size-1)=&
467	il_gextents(il_pelist(:),1,il_port)
468
469	!The max global index for the current 1D decomposition.
470	il_sum=sum(il_extentsx(0:ig_mpp_io_comm_size-1))
471
472	#ifdef __VERBOSE
473	WRITE(nulprt,*)'psmile_def_domains: port: ',il_port &
474	,'(/1,il_sum,1,1/),pelist',1,il_sum,1,1,il_pelist &
475	,'extents',il_extentsx
476	#endif
477	!Define a pseudo 2D domain
478	il_extentsy(0)=1
479	call mpp_define_domains((/1,il_sum,1,1/),(/ig_mpp_io_comm_size,1/),&
480	domain(il_port),&
481	pelist=il_pelist,&
482	xextent=il_extentsx,&
483	yextent=il_extentsy(0:0))
484	call mpp_get_domain_components( domain(il_port), xdom(il_port),&
485	ydom(il_port) )
486
487	else if( mydist(CLIM_Strategy,il_port) .eq. CLIM_Box) then
488
489	#ifdef __VERBOSE
490	WRITE(nulprt,*)'psmile_def_domains: port: ',il_port,' CLIM_Box'
491	! WRITE(nulprt,*)'psmile_def_domains: port: ',il_port,' CLIM_Box', &
492	! mydist(:,il_port),'::',myport(:,il_port)
493	#endif
494	!Build a list of the local extents which exists on the procs involved with
495	!coupling
496	!Block sizes in x
497	call MPI_Allgather(mydist(CLIM_SizeX+1,il_port),1,MPI_INTEGER,&
498	il_gextents(0,1,il_port),il_recvcount,MPI_INTEGER,&
499	ig_mpp_io_comm,id_error)
500	!Block sizes in y
501	call MPI_Allgather(mydist(CLIM_Segments,il_port),1,MPI_INTEGER,&
502	il_gextents(0,2,il_port),il_recvcount,MPI_INTEGER,&
503	ig_mpp_io_comm,id_error)
504	!Offsets of the blocks
505	call MPI_Allgather(mydist(CLIM_Offset+1,il_port),1,MPI_INTEGER,&
506	il_offsets(1,1,il_port),il_recvcount,MPI_INTEGER,&
507	ig_mpp_io_comm,id_error)
508
509	call indexi(ig_mpp_io_comm_size &
510	,il_offsets(1:ig_mpp_io_comm_size,1,il_port),&
511	il_pelist)
512	il_pelist(:)=il_pelist(:)-1
513
514	il_extentsx(0:ig_mpp_io_comm_size-1)=&
515	il_gextents(il_pelist(:),1,il_port)
516	il_extentsy(0:ig_mpp_io_comm_size-1)=&
517	il_gextents(il_pelist(:),2,il_port)
518	#ifdef __VERBOSE
519	! write(nulprt,*)'psmile_def_domains:il_extentsx',il_extentsx
520	! write(nulprt,*)'psmile_def_domains:il_extentsy',il_extentsy
521	#endif
522	! At this point I have to think about the domain decomposition
523	! 1st assumption : The leading dimension is a global parameter for
524	! the model. So, I don't do further checks on that.
525
526	il_ldx=mydist(CLIM_Segments+3,il_port)-mydist(CLIM_Offset+1,il_port)
527	#ifdef __VERBOSE
528	! write(nulprt,*)'psmile_def_domains: LDX = ',il_ldx
529	#endif
530
531
532	!Count the number of offsets which are less than LDX to get the partitioning
533	!in x-direction
534	ig_layout(1)=&
535	count(mask=il_offsets(:,1,il_port).lt.il_ldx)
536
537	if(ig_layout(1).gt.0) then
538	ig_layout(2)=ig_mpp_io_comm_size/ig_layout(1)
539	if(ig_layout(1)*ig_layout(2).ne.ig_mpp_io_comm_size) then
540	lg_is_odd_distribution(il_port)=.true.
541	else
542	!check within a column y for different extents of y
543	lg_is_odd_distribution(il_port)=.false.
544	do il_j=0,ig_layout(2)-1
545	do il_i=1,ig_layout(1)-1
546	lg_is_odd_distribution(il_port)= &
547	il_extentsy(ig_layout(1)*il_j).ne.&
548	il_extentsy(il_i+ig_layout(1)*il_j)
549	#ifdef __VERBOSE
550	! write(nulprt,*)'i j extentsy ',il_i,il_j, &
551	! il_extentsy(ig_layout(1)*il_j), &
552	! il_extentsy(il_i+ig_layout(1)*il_j)
553	#endif
554	enddo
555	enddo
556
557	!check within a column x for different extents of x
558	do il_i=0,ig_layout(2)-1
559	do il_j=1,ig_layout(1)-1
560	lg_is_odd_distribution(il_port)= &
561	il_extentsx(il_i).ne.il_extentsx(il_i+ig_layout(1)*il_j)
562	#ifdef __VERBOSE
563	! write(nulprt,*)'i j extentsy ',il_i,il_j, &
564	! il_extentsx(il_i), &
565	! il_extentsx(il_i+ig_layout(1)*il_j)
566	#endif
567	enddo
568	enddo
569
570	endif
571
572	!Okay, fine! There are global divisors of the x- and y-axis!
573	if(.not.lg_is_odd_distribution(il_port)) then
574
575
576	il_ngx= il_ldx
577	il_ngy=il_offsets(il_pelist(ig_mpp_io_comm_size-1)+1,1,il_port)&
578	/il_ngx + il_extentsy(ig_mpp_io_comm_size-1)
579
580	#ifdef __VERBOSE
581	WRITE(nulprt,*)'psmile_def_domains: port: ',il_port &
582	,'(/1,il_ngx,1,il_ngy/),pelist' &
583	,1,il_ngx,1,il_ngy,il_pelist
584	#endif
585
586	call mpp_define_domains((/1,il_ngx,1,il_ngy/) &
587	,ig_layout,domain(il_port) &
588	,xextent=il_extentsx(0:ig_mpp_io_comm_size-1:ig_layout(2)) &
589	,yextent=il_extentsy(0:ig_mpp_io_comm_size-1:ig_layout(1)) &
590	,pelist=il_pelist)
591	call mpp_get_domain_components( domain(il_port), xdom(il_port),&
592	ydom(il_port) )
593
594	else
595
596	!Okay, at this point it seems that have we have a block-structured
597	!distribution. Need to calculate the offsets of the blocks in x- and y-
598	!direction. These offsets migth be provided by prism_set_offset.
599	ll_mask(:,:)=.false.
600	do il_i=0,ig_mpp_io_comm_size-1
601	ll_mask(il_pelist(il_i),il_pelist(il_i))=.true.
602	enddo
603	il_ngx=il_ldx
604	il_offsetx(:)=il_offsets(il_pelist(0:ig_mpp_io_comm_size-1)+1,1,il_port) &
605	-il_ngx* &
606	(il_offsets(il_pelist(0:ig_mpp_io_comm_size-1)+1,1,il_port) &
607	/il_ngx)
608	il_offsety(:)=il_offsets(il_pelist(0:ig_mpp_io_comm_size-1)+1,1,il_port) &
609	/il_ngx
610	il_ngy=il_offsety(ig_mpp_io_comm_size-1) &
611	+il_extentsy(ig_mpp_io_comm_size-1)
612	ig_layout(1)=ig_mpp_io_comm_size
613	ig_layout(2)=ig_mpp_io_comm_size
614
615	#ifdef __VERBOSE
616	WRITE(nulprt,*)'psmile_def_domains: port: ',il_port &
617	,' :fall back: (/1,il_ngx,1,il_ngy/),pelist' &
618	,1,il_ngx,1,il_ngy,il_pelist,il_offsetx(:) &
619	,il_offsety(:)
620	#endif
621
622	call mpp_define_domains((/1,il_ngx,1,il_ngy/) &
623	,ig_layout,domain(il_port) &
624	,xextent=il_extentsx &
625	,yextent=il_extentsy &
626	,pelist=il_pelist,maskmap=ll_mask &
627	,offsetx=il_offsetx,offsety=il_offsety)
628	call mpp_get_domain_components( domain(il_port), xdom(il_port),&
629	ydom(il_port) )
630	endif
631
632	else
633
634	write(nulprt, * ) 'psmile_def_domains: ig_layout(1) =0! '
635	call MPI_ABORT(mpi_comm,0,mpi_err)
636	endif
637
638	else if( mydist(CLIM_Strategy,il_port) .eq. CLIM_Orange) then
639
640	write(nulprt, * ) 'psmile_def_domains: Strategy CLIM_orange',&
641	'not yet implemented for I/O!'
642	call MPI_ABORT(mpi_comm,0,mpi_err)
643
644	endif
645	endif
646	enddo
647
648	!Cleanup
649	deallocate(il_pelist)
650	deallocate(il_extentsx)
651	deallocate(il_extentsy)
652	deallocate(il_gextents)
653	deallocate(il_offsets)
654	deallocate(il_offsetx)
655	deallocate(il_offsety)
656	deallocate(ll_mask)
657
658	#ifdef __VERBOSE
659	WRITE(nulprt,*)'End - - psmile_def_domains:'
660	call flush(nulprt)
661	#endif
662
663
664	end subroutine psmile_def_domains
665
666	subroutine psmile_def_files(id_error)
667	!--------------------------------------------------------------------------
668	!
669	! This subroutine open files acccording to the kewords
670	! OUTPUT: write by model issueing prism_put
671	! INPUT: read by issueing prism_get
672	! EXPOUT: write by model issueing prism_put and prims_get after receive.
673	! Usage for debugging data exchange model->Oasis->model
674	! IGNOUT: write by model issueing prism_put and prims_get after receive.
675	! Usage for debugging data exchange model->model
676	!
677	! Called at the end of prism_enddef.
678	!
679	use mod_kinds_model
680	use mod_prism_proto
681	use mod_comprism_proto
682	use mod_psmile_io
683	use mod_psmile_io_interfaces,only:combine_with_date
684	use mod_psmile_date_and_time
685
686	implicit none
687
688	#include <mpif.h>
689
690	integer(kind=ip_intwp_p),intent(out)::id_error
691	!Local variables
692	integer(kind=ip_intwp_p)::il_i,il_port,il_unit
693	!Test integer(kind=ip_intwp_p),allocatable::il_units(:)
694	logical:: ll_opened,ll_exist
695
696	#ifdef __VERBOSE
697	write(nulprt,*)'Start - - psmile_def_files'
698	call flush(nulprt)
699	#endif
700
701	if(ig_color/ig_mynummod.ne.1000) then
702	WRITE(nulprt,*)'psmile_def_files: Called by proc which is not' &
703	,' involved with coupling!'
704	call flush(nulprt)
705	call MPI_ABORT(mpi_comm,0,mpi_err)
706	endif
707
708	if((.not.allocated(domain))) then
709	WRITE(nulprt,*)'psmile_def_files:' &
710	,' Not called after psmile_def_domains!'
711	call flush(nulprt)
712	call MPI_ABORT(mpi_comm,0,mpi_err)
713	endif
714
715
716	id_error=CLIM_Ok
717
718	!Find unused units.
719	!Problem is here that mpi codes usually are scattering/collecting from/on one PE
720	! which read/write. The rank of that PE is not known to PSMILE.
721	!Algorithm: Let each PE involved with coupling do a search beginning from
722	!PSMILE_MAX_UNIT. Check if the result is the same on each couple task.
723	!Currently I don't know if the check is absolutely necessary.
724	!Moreover, I assume that MPI is doing private I/O for each task, means
725	!each opened file is handled individually per task by the OS .
726	!I made my experience with UNICOS where file handling where global per default
727	!running MPI in autotasking mode. In that case the user had to specify a
728	!certain FFIO flag to get the normal behaviour.
729
730	!Test allocate(il_units(1:ig_mpp_io_comm_size))
731	allocate(ig_units(1:nports))
732	!RV:22.01.2003
733	ig_units=-1
734
735	il_unit=PSMILE_MAX_UNIT
736	do il_port=1,nports
737	IF (ig_def_state(il_port) .eq. ip_ignout .or. &
738	ig_def_state(il_port) .eq. ip_expout .or. &
739	ig_def_state(il_port) .eq. ip_output .or. &
740	ig_def_state(il_port) .eq. ip_input ) THEN
741	do
742	inquire(unit=il_unit,exist=ll_exist,opened=ll_opened)
743	if( ll_exist .and. (.not.ll_opened)) exit
744	il_unit=il_unit-1
745	if(il_unit.eq.PSMILE_MAX_RESERVED_UNIT) &
746	il_unit=PSMILE_MIN_RESERVED_UNIT-1
747	if(il_unit.eq.7) then
748	write(nulprt,*)'psmile_def_files: Could not find unit!'
749	call MPI_ABORT(mpi_comm,0,mpi_err)
750	endif
751	enddo
752	ig_units(il_port)=il_unit
753
754	#ifdef __VERBOSE
755	write(nulprt,*)'psmile_def_files: Found FORTRAN I/O unit ', il_unit
756	call flush(nulprt)
757	#endif
758
759	il_unit=il_unit-1
760	endif
761	enddo
762
763	!Opening files!
764	!I am opening the files such that read / write is performed on one file.
765	!One can easily change that to distributed i/o changing the fileset attribute.
766
767	allocate(ig_action(1:nports))
768	allocate(ig_form(1:nports))
769	allocate(ig_threading(1:nports))
770	allocate(ig_fileset(1:nports))
771	allocate(cg_my_input_file(1:nports))
772	allocate(cg_output_file(1:nports))
773
774	! call_mpp_sync()
775	do il_port=1,nports
776	if(ig_def_state(il_port).eq.ip_input) then
777
778	ig_action(il_port)=MPP_RDONLY
779	ig_form(il_port)=MPP_NETCDF
780	ig_threading(il_port)=MPP_MULTI
781	ig_fileset(il_port)=MPP_SINGLE
782
783	!I create here the name of the input file according to the convention
784	!cports(il_port)_in.<year>-<month>-<day>T<hh>:<mm>:ss.nc
785	!Arnaud, if you like you prefer approach of reading the name from the namcouple
786	!initialize cg_my_input_file accordingly and comment out the call of
787	!combine_with_date.
788
789	! call combine_with_date(cports(il_port),'in',il_my_initial_date &
790	! ,cg_my_input_file(il_port))
791
792	call mpp_open( ig_units(il_port),trim(cg_def_inpfile(il_port)) &
793	, action=ig_action(il_port) &
794	, form=ig_form(il_port) &
795	, threading=ig_threading(il_port) &
796	, fileset=ig_fileset(il_port))
797
798	elseif(ig_def_state(il_port).eq.ip_output) then
799
800	ig_action(il_port)=MPP_OVERWR
801	ig_form(il_port)=MPP_NETCDF
802	ig_threading(il_port)=MPP_SINGLE
803
804	call combine_with_date(cports(il_port),'out',ig_inidate &
805	,cg_output_file(il_port))
806
807	call mpp_open( ig_units(il_port), trim(cg_output_file(il_port)) &
808	, action=ig_action(il_port) &
809	, form=ig_form(il_port) &
810	, threading=ig_threading(il_port))
811
812
813	elseif(ig_def_state(il_port).eq.ip_expout) then
814
815	ig_action(il_port)=MPP_OVERWR
816	ig_form(il_port)=MPP_NETCDF
817	ig_threading(il_port)=MPP_SINGLE
818
819	call combine_with_date(cports(il_port),'out',ig_inidate &
820	,cg_output_file(il_port))
821
822	call mpp_open( ig_units(il_port),trim(cg_output_file(il_port)) &
823	, action=ig_action(il_port) &
824	, form=ig_form(il_port) &
825	, threading=ig_threading(il_port))
826
827	elseif(ig_def_state(il_port).eq.ip_ignout) then
828
829	ig_action(il_port)=MPP_OVERWR
830	ig_form(il_port)=MPP_NETCDF
831	ig_threading(il_port)=MPP_SINGLE
832
833	call combine_with_date(cg_ignout_field(il_port),'out',ig_inidate &
834	,cg_output_file(il_port))
835
836	call mpp_open( ig_units(il_port), trim(cg_output_file(il_port)) &
837	, action=ig_action(il_port) &
838	, form=ig_form(il_port) &
839	, threading=ig_threading(il_port))
840	endif
841	enddo
842
843
844	#ifdef __VERBOSE
845	write(nulprt,*)'End - - psmile_def_files'
846	call flush(nulprt)
847	#endif
848	end subroutine psmile_def_files
849
850	subroutine psmile_def_metadata(id_error)
851	!--------------------------------------------------------------------------
852	!
853	!This routine defines the metadata headers of netcdf files case or writing.
854	!In case or reading it defines the index to lookup a certain variable in
855	!a netcdf file
856	! Called at the end of prism_enddef.
857	use mod_kinds_model
858	use mod_psmile_io_interfaces,only:combine_with_date
859	use mod_prism_proto
860	use mod_comprism_proto
861	use mod_psmile_io
862	implicit none
863
864	#include <mpif.h>
865
866	integer(kind=ip_intwp_p),intent(out)::id_error
867	!Local
868	character(len=64)::cl_varname
869	integer(kind=ip_intwp_p)::il_port,il_i
870	integer(kind=ip_intwp_p)::il_xbegin,il_xend,il_ybegin,il_yend
871	integer(kind=ip_intwp_p)::il_ndim,il_nvar,il_natt,il_ntime
872	integer(kind=ip_intwp_p)::il_crn,il_ii,il_jj
873
874	type(fieldtype),allocatable::il_vars(:)
875	!
876	! Dimensions of the lons and lats contained in the OASIS 3 grid file
877	integer(kind=ip_intwp_p)::lens_lat(4),lens_lon(4)
878	REAL(kind=ip_realwp_p), ALLOCATABLE :: lon(:,:),lat(:,:)
879	REAL(kind=ip_realwp_p), ALLOCATABLE :: loncrn(:,:,:),latcrn(:,:,:)
880	LOGICAL :: lg_exist
881	CHARACTER(len=8):: cg_clim_cgrdnamnc, clstrg
882	INTEGER(kind=ip_intwp_p):: nc_grdid, il_varid, il_x, il_y
883	REAL(kind=ip_realwp_p), ALLOCATABLE :: rl_lons(:,:), rl_lats(:,:)
884	REAL(kind=ip_realwp_p), ALLOCATABLE :: rl_lonscrn(:,:,:), rl_latscrn(:,:,:)
885	REAL(kind=ip_realwp_p):: difflon
886	REAL(kind=ip_realwp_p):: threesixty
887	LOGICAL :: llylat,llxlon, ll_notfound
888
889	#ifdef __VERBOSE
890	write(nulprt,*)'Start -- psmile_def_metada'
891	call flush(nulprt)
892	#endif
893
894	if(ig_color/ig_mynummod.ne.1000) then
895	WRITE(nulprt,*)'psmile_def_metdata: Called by proc which is not' &
896	,' involved with coupling!'
897	call flush(nulprt)
898	call MPI_ABORT(mpi_comm,0,mpi_err)
899	endif
900
901	if((.not.allocated(ig_units)).and.(.not.allocated(domain))) then
902	WRITE(nulprt,*)'psmile_def_metdata:' &
903	,' Not called after psmile_def_domains' &
904	,' and psmile_def_files!'
905	call flush(nulprt)
906	call MPI_ABORT(mpi_comm,0,mpi_err)
907	endif
908
909
910	id_error=CLIM_Ok
911	threesixty=360.0
912
913	allocate(ig_var_index(nports))
914	allocate(ig_vars(1:nports))
915	allocate(ig_ntimes(1:nports))
916
917	!RV:22.01.2003
918	#ifdef __READ_BY_LOOKUP
919	allocate(dg_times(1:nports))
920	!AC:31/01
921	DO il_i = 1,nports
922	nullify(dg_times(il_i)%rl_times)
923	ENDDO
924	#endif
925
926	allocate(x_axis(1:nports))
927	allocate(y_axis(1:nports))
928	allocate(c_axis(1:nports))
929	#ifdef __VERBOSE
930	write(nulprt,*)'psmile_def_metada: Number of corners: ',ig_noc
931	call flush(nulprt)
932	#endif
933	allocate(nbr_corners(1:nports))
934	nbr_corners=ig_noc
935	!Was not allocated. Oasis sends 2D arrays. In case of bundles or 3d arrays
936	!we have to declare z_axis + plus one for bundles.
937	!Left here as a reminder!
938	! allocate(z_axis(1:nports))
939	allocate(t_axis(1:nports))
940	allocate(field(1:nports))
941	allocate(lonij(1:nports))
942	allocate(latij(1:nports))
943	allocate(crnlonij(1:nports))
944	allocate(crnlatij(1:nports))
945
946	!The following files are allocated here since no informations
947	!are given in Oasis regarding to long names and units. I am filling them with
948	!dummy informations. I put the data structure into mod_psmile_io.
949	allocate(cports_lgname(1:nports))
950	allocate(cports_units(1:nports))
951	cports_lgname='BLANK'
952	cports_units='1'
953
954
955	DO il_port=1,nports
956	IF(ig_def_state(il_port).EQ.ip_input) THEN
957
958	!Get information with regard to no. of dimensions, variables, attributes
959	!and time stamps.
960	call mpp_get_info(ig_units(il_port)&
961	,il_ndim,il_nvar,il_natt,il_ntime)
962
963	ig_ntimes(il_port)=il_ntime
964
965	allocate(il_vars(1:il_nvar))
966	!RV,16.12.2002:Get the all ids of the fields within the file.
967	call mpp_get_fields (ig_units(il_port),il_vars(:))
968
969	!Loop over variables and find the index of cports(il_port).
970	ig_var_index(il_port)=-1
971	do il_i=1,il_nvar
972
973	call mpp_get_atts(il_vars(il_i),name=cl_varname)
974	if(trim(cl_varname).eq.trim(cports(il_port)) ) then
975	ig_var_index(il_port)=il_i
976	ig_vars(il_port)=il_vars(il_i)
977	endif
978
979	enddo
980
981	!Check: Did I find the variable? If not, the user has given the wrong
982	! input file.
983	if(ig_var_index(il_port).lt.0) then
984
985	write(nulprt, * ) &
986	'psmile_def_metadata: Variable name given in namcouple '&
987	,cports(il_port) &
988	,'could not be found in file ',cg_my_input_file(il_port) &
989	,' port number: ',il_port
990	call MPI_ABORT(mpi_comm,0,mpi_err)
991
992	endif
993
994	deallocate(il_vars)
995
996	elseif((ig_def_state(il_port).eq.ip_output).or. &
997	(ig_def_state(il_port).eq.ip_expout).or. &
998	(ig_def_state(il_port).eq.ip_ignout))then
999
1000	!Define the global axis X and Y.
1001	!I have no informations on the units of the axis. I assume that
1002	!the axis are indices in order to perform a mapping like lat(i,j).
1003	call mpp_get_global_domain(domain(il_port) &
1004	,xbegin=il_xbegin,xend=il_xend &
1005	,ybegin=il_ybegin,yend=il_yend)
1006	call mpp_write_meta(ig_units(il_port) &
1007	,'Conventions' &
1008	,cval='CF-1.0')
1009	#ifdef __VERBOSE
1010	write(nulprt,*)'psmile_def_metadata: x-axis ',il_xbegin,il_xend
1011	#endif
1012
1013	call mpp_write_meta(ig_units(il_port),x_axis(il_port) &
1014	,'X','1','global index space for x' &
1015	,domain=xdom(il_port) &
1016	,data=(/(il_i-1.,il_i=il_xbegin,il_xend)/) )
1017
1018	#ifdef __VERBOSE
1019	write(nulprt,*)'psmile_def_metadata: y-axis ',il_ybegin,il_yend
1020	call flush(nulprt)
1021	#endif
1022
1023	call mpp_write_meta(ig_units(il_port),y_axis(il_port) &
1024	,'Y','1','global index space for y' &
1025	,domain=ydom(il_port) &
1026	,data=(/(il_i-1.,il_i=il_ybegin,il_yend)/) )
1027	#ifdef __VERBOSE
1028	write(nulprt,*)'psmile_def_metadata: c-axis '
1029	#endif
1030	call mpp_write_meta(ig_units(il_port),c_axis(il_port) &
1031	,'C','1','global index space for c' &
1032	,data=(/(il_i-1.,il_i=1,nbr_corners(il_port))/) )
1033
1034
1035	!Define the unlimited time axis
1036	call combine_with_date('second','since',ig_inidate &
1037	,cl_varname)
1038
1039	call mpp_write_meta( ig_units(il_port), t_axis(il_port) &
1040	, 'time', trim(cl_varname), 'Time' )
1041	!Define the real coordinates , here lattitudes and longitudes.
1042	! I have no information with regard to latitude and longitude
1043	!RV:17.04.2003 Wright this information anyway to pass the CF check.
1044	call mpp_write_meta( ig_units(il_port), latij(il_port) &
1045	, (/x_axis(il_port),y_axis(il_port)/) &
1046	, 'lat', 'degrees_north' &
1047	, 'latitude')
1048	call mpp_write_meta( ig_units(il_port), lonij(il_port) &
1049	, (/x_axis(il_port),y_axis(il_port)/) &
1050	, 'lon', 'degrees_east' &
1051	, 'longitude')
1052	!rv,sgi Number of corners has to be the leading dimension according the
1053	!rv,sgi CF convention
1054	call mpp_write_meta( ig_units(il_port), crnlatij(il_port) &
1055	, (/c_axis(il_port),x_axis(il_port),y_axis(il_port)/) &
1056	, 'bounds_lat', 'degrees_north' &
1057	, 'corner latitudes')
1058	call mpp_write_meta( ig_units(il_port), crnlonij(il_port) &
1059	, (/c_axis(il_port),x_axis(il_port),y_axis(il_port)/) &
1060	, 'bounds_lon', 'degrees_east' &
1061	, 'corner longitudes')
1062	!Define that variable cports(il_port) of unit cports_units(il_port) and
1063	!longname cports_lgname(il_port) will be written as field field(il_port)
1064	!according to the three axis x_axis(il_port),y_axis(il_port) and
1065	!t_axis(il_port). The precision is double (pack=1).
1066	!RV
1067	! Initialize cports_lgname and cports_units according to the CF name table
1068	! exploiting the field label.
1069	!
1070	IF(allocated(cfldlab)) then
1071	IF(ig_def_numlab(il_port).gt.0.and. &
1072	ig_def_numlab(il_port).le.PSMILE_MAX_CF_TABLE_ENTRIES) then
1073	cports_lgname(il_port)=cfldlab(ig_def_numlab(il_port))
1074	cports_units(il_port)=cfunitslab(ig_def_numlab(il_port))
1075	ENDIF
1076	ENDIF
1077	!RV
1078	IF (ig_def_state(il_port).EQ.ip_ignout) THEN
1079	call mpp_write_meta( ig_units(il_port), field(il_port) &
1080	, (/x_axis(il_port),y_axis(il_port) &
1081	,t_axis(il_port)/) &
1082	,trim(cg_ignout_field(il_port)) &
1083	,trim(cports_units(il_port)) &
1084	, trim(cports_lgname(il_port)), pack=1 )
1085	#ifdef __DEBUG
1086	write(nulprt,*)'psmile_def_metadata 1:',ig_def_numlab(il_port) &
1087	,' ',trim(cg_ignout_field(il_port)) &
1088	,' ',trim(cports_lgname(il_port)) &
1089	,' ',trim(cports_units(il_port))
1090	call flush(nulprt)
1091	#endif
1092	ELSE
1093	call mpp_write_meta( ig_units(il_port), field(il_port) &
1094	, (/x_axis(il_port),y_axis(il_port) &
1095	,t_axis(il_port)/) &
1096	,trim(cports(il_port)) &
1097	,trim(cports_units(il_port)) &
1098	, trim(cports_lgname(il_port)), pack=1 )
1099	#ifdef __DEBUG
1100	write(nulprt,*)'psmile_def_metadata 2:',ig_def_numlab(il_port) &
1101	,' ', trim(cports(il_port)) &
1102	,' ',trim(cports_lgname(il_port)) &
1103	,' ',trim(cports_units(il_port))
1104	call flush(nulprt)
1105	#endif
1106	ENDIF
1107
1108	!
1109	! Define an additional attribute bounds for lon and lats.
1110	!
1111	call mpp_write_meta(ig_units(il_port) &
1112	,mpp_get_id(latij(il_port)) &
1113	, 'bounds',cval='bounds_lat')
1114	call mpp_write_meta(ig_units(il_port) &
1115	,mpp_get_id(lonij(il_port)) &
1116	, 'bounds',cval='bounds_lon')
1117	!Define an additional attribute that field(il_port) lives on the coordinates
1118	!lon and lat. This is a CF requirement!
1119
1120	call mpp_write_meta(ig_units(il_port) &
1121	,mpp_get_id(field(il_port)) &
1122	, 'coordinates',cval='lon lat')
1123
1124	!Flush the metadata information to the file connected to unit
1125	!ig_units(il_port).
1126
1127	call mpp_write(ig_units(il_port), x_axis(il_port) )
1128	call mpp_write(ig_units(il_port), y_axis(il_port) )
1129	!
1130	!rv Metadata should be flushed from the struct c_axis
1131	!rv into the NetCDF file
1132	!
1133	call mpp_write(ig_units(il_port), c_axis(il_port) )
1134	!
1135	IF (mpp_pe() .EQ. mpp_root_pe()) THEN
1136	WRITE(nulprt,*)'mod_psmile_io - - check grids.nc file:'
1137	cg_clim_cgrdnamnc=cg_clim_cgrdnam//'.nc'
1138	INQUIRE(FILE = cg_clim_cgrdnamnc, EXIST = lg_exist)
1139	IF (.NOT. lg_exist) THEN
1140	WRITE(nulprt,*) ' - no grids.nc file found'
1141	WRITE(nulprt,*) 'lats and lons will not be written to output files'
1142	CALL flush(nulprt)
1143	ELSE
1144	WRITE(nulprt,*)'mod_psmile_io - - grids.nc file found'
1145	! Read the lats and the lons
1146	CALL mpp_open(nc_grdid, 'grids', &
1147	action=MPP_RDONLY, &
1148	form=MPP_NETCDF, &
1149	threading=MPP_SINGLE)
1150	CALL mpp_get_info(nc_grdid, il_ndim, il_nvar, &
1151	il_natt, il_ntime )
1152	ALLOCATE(il_vars(1:il_nvar))
1153	CALL mpp_get_fields (nc_grdid,il_vars(:))
1154	!
1155	! Longitudes
1156	!
1157	clstrg = cga_clim_locator(il_port)//cg_clim_lonsuf
1158	il_i = 0
1159	ll_notfound = .TRUE.
1160	DO WHILE (ll_notfound)
1161	il_i = il_i + 1
1162	IF (il_i .GT. il_nvar) THEN
1163	WRITE(nulprt,*)'psmile_def_metdata: Pb in file grids.nc: No longitudes' &
1164	,'for grid', cga_clim_locator(il_port)
1165	WRITE(nulprt,*)'Calling MPI_ABORT '
1166	CALL flush(nulprt)
1167	CALL MPI_ABORT(mpi_comm,0,mpi_err)
1168	ENDIF
1169	CALL mpp_get_atts(il_vars(il_i),name=cl_varname,siz=lens_lon)
1170	!
1171	!RV,sgi ALLOCATE(lon(1:il_xend-il_xbegin+1,1:il_yend-il_ybegin+1))
1172	!
1173	IF(TRIM(cl_varname).EQ.TRIM(clstrg) ) THEN
1174	ll_notfound = .FALSE.
1175	#ifdef __VERBOSE
1176	write(nulprt,*)'psmile_def_metadata: cl_varname= ' &
1177	,trim(cl_varname)
1178	call flush(nulprt)
1179	#endif
1180	IF(.NOT.ALLOCATED(lon)) &
1181	ALLOCATE(lon(1:il_xend-il_xbegin+1,1:il_yend-il_ybegin+1))
1182	ALLOCATE(rl_lons(lens_lon(1),lens_lon(2)))
1183
1184	!NOTE: I am using tindex = 1 in order to circumvent a
1185	!peculiarity of mpp_io: grids.nc contains the lons/lats
1186	!declared as lon(i,j)/lat(i,j). mpp_io interprets i and j as
1187	!axes. However, the axis i and j have no associated data.
1188	!Therefore, i and j are NOT taken as spatial axes but
1189	!as something related to a time stamp. Until that "feature"
1190	!is corrected use the call below.
1191	!
1192	CALL mpp_read(nc_grdid,il_vars(il_i),rl_lons,tindex=1)
1193	!
1194	!rv,sgi<
1195	!If the user has specified INVERT/REVERSE and
1196	!$CORLAT=NORSUD and/or $CORLON=ESTWST in the
1197	!preprocessing section of the field we INVERT the
1198	!ordering of longitudes such that data and grid are
1199	!are matching. The assertion is that
1200	!grid informations are initialized SUDNOR and WSTEST.
1201
1202	if((ig_def_invert(il_port).gt.1).or. &
1203	(ig_def_reverse(il_port).gt.1)) then
1204	llxlon=ig_def_invert(il_port).eq.3 &
1205	.or.ig_def_invert(il_port).eq.5 &
1206	.or.ig_def_reverse(il_port).eq.3 &
1207	.or.ig_def_reverse(il_port).eq.5
1208	llylat=ig_def_invert(il_port).eq.2 &
1209	.or.ig_def_invert(il_port).eq.5 &
1210	.or.ig_def_reverse(il_port).eq.2 &
1211	.or.ig_def_reverse(il_port).eq.5
1212	call inv2d (rl_lons,lens_lon(1),lens_lon(2) &
1213	,llxlon,llylat)
1214	endif
1215	!rv,sgi>
1216	!Reshape the read array such that it matches the current
1217	!axes of the output file.
1218	!
1219	lon=RESHAPE(rl_lons,shape=SHAPE(lon))
1220	!
1221	!Dump it ot the output file.
1222	!
1223	CALL mpp_write(ig_units(il_port),lonij(il_port),lon)
1224	DEALLOCATE(rl_lons)
1225	ENDIF
1226	ENDDO
1227	!
1228	!Latitudes
1229	!
1230	clstrg = cga_clim_locator(il_port)//cg_clim_latsuf
1231	!
1232	il_i = 0
1233	ll_notfound = .TRUE.
1234	DO WHILE (ll_notfound)
1235	il_i = il_i + 1
1236	IF (il_i .GT. il_nvar) THEN
1237	WRITE(nulprt,*)'psmile_def_metdata: Pb in file grids.nc: No latitudes' &
1238	,'for grid', cga_clim_locator(il_port)
1239	WRITE(nulprt,*)'Calling MPI_ABORT '
1240	CALL flush(nulprt)
1241	CALL MPI_ABORT(mpi_comm,0,mpi_err)
1242	ENDIF
1243
1244	!
1245	CALL mpp_get_atts(il_vars(il_i),name=cl_varname,siz=lens_lat)
1246	IF(TRIM(cl_varname).EQ.TRIM(clstrg) ) THEN
1247	ll_notfound = .FALSE.
1248	#ifdef __VERBOSE
1249	write(nulprt,*)'psmile_def_metadata: cl_varname= ' &
1250	,trim(cl_varname)
1251	call flush(nulprt)
1252	#endif
1253	ALLOCATE(lat(1:il_xend-il_xbegin+1,1:il_yend-il_ybegin+1))
1254	ALLOCATE(rl_lats(lens_lat(1),lens_lat(2)))
1255
1256	CALL mpp_read(nc_grdid,il_vars(il_i),rl_lats,tindex=1)
1257	!rv,sgi<
1258	if((ig_def_invert(il_port).gt.1).or. &
1259	(ig_def_reverse(il_port).gt.1)) then
1260	llxlon=ig_def_invert(il_port).eq.3 &
1261	.or.ig_def_invert(il_port).eq.5 &
1262	.or.ig_def_reverse(il_port).eq.3 &
1263	.or.ig_def_reverse(il_port).eq.5
1264	llylat=ig_def_invert(il_port).eq.2 &
1265	.or.ig_def_invert(il_port).eq.5 &
1266	.or.ig_def_reverse(il_port).eq.2 &
1267	.or.ig_def_reverse(il_port).eq.5
1268	call inv2d (rl_lats,lens_lat(1),lens_lat(2) &
1269	,llxlon,llylat)
1270	endif
1271	!rv,sgi>
1272	lat=reshape(rl_lats,shape=shape(lat))
1273	call mpp_write(ig_units(il_port),latij(il_port),lat)
1274	DEALLOCATE(rl_lats)
1275	ENDIF
1276	ENDDO
1277	!
1278	! Longitude corners
1279	!
1280	clstrg = cga_clim_locator(il_port)//crn_clim_lonsuf
1281
1282	ll_notfound = .TRUE.
1283
1284	DO il_i=1,il_nvar
1285
1286	CALL mpp_get_atts(il_vars(il_i),name=cl_varname,siz=lens_lon)
1287	IF(TRIM(cl_varname).EQ.TRIM(clstrg) ) THEN
1288	ll_notfound = .FALSE.
1289	#ifdef __VERBOSE
1290	write(nulprt,*) &
1291	'psmile_def_metada, lon, clstrg: ',clstrg,lens_lon
1292	#endif
1293	nbr_corners(il_port)=lens_lon(3)
1294	ALLOCATE(rl_lonscrn(lens_lon(1),lens_lon(2) &
1295	,nbr_corners(il_port)))
1296	ALLOCATE(loncrn(nbr_corners(il_port) &
1297	,1:il_xend-il_xbegin+1 &
1298	,1:il_yend-il_ybegin+1))
1299	lens_lon(1)=il_xend-il_xbegin+1
1300	lens_lon(2)=il_yend-il_ybegin+1
1301	lens_lon(3)=nbr_corners(il_port)
1302	!
1303	CALL mpp_read(nc_grdid,il_vars(il_i),rl_lonscrn,tindex=1)
1304	!rv,sgi<
1305	if((ig_def_invert(il_port).gt.1).or. &
1306	(ig_def_reverse(il_port).gt.1)) then
1307	llxlon=ig_def_invert(il_port).eq.3 &
1308	.or.ig_def_invert(il_port).eq.5 &
1309	.or.ig_def_reverse(il_port).eq.3 &
1310	.or.ig_def_reverse(il_port).eq.5
1311	llylat=ig_def_invert(il_port).eq.2 &
1312	.or.ig_def_invert(il_port).eq.5 &
1313	.or.ig_def_reverse(il_port).eq.2 &
1314	.or.ig_def_reverse(il_port).eq.5
1315	call inv2dcrn (rl_lonscrn,lens_lon(1),lens_lon(2) &
1316	,nbr_corners(il_port),llxlon,llylat)
1317	endif
1318	!rv,sgi>
1319
1320	!
1321	!Reshape the read array such that it matches the current
1322	!axes of the output file.
1323	!
1324	call transp_crn2cf(rl_lonscrn,loncrn,lens_lon)
1325	!rv,sgi<
1326	!
1327	! Check for cyclic boundaries !
1328	! For plot routines using a mesh-fill method
1329	! the distances of the cell centers to the corners
1330	! have to be less equal 360. degrees.
1331	!
1332	do il_crn=1,nbr_corners(il_port)
1333	do il_jj=1,il_yend-il_ybegin+1
1334	do il_ii=1,il_xend-il_xbegin+1
1335	difflon=lon(il_ii,il_jj)- &
1336	loncrn(il_crn,il_ii,il_jj)
1337	if(difflon.gt.300.0) then
1338	loncrn(il_crn,il_ii,il_jj)= &
1339	loncrn(il_crn,il_ii,il_jj)+threesixty
1340	elseif(difflon.lt.-300.0) then
1341	loncrn(il_crn,il_ii,il_jj)= &
1342	loncrn(il_crn,il_ii,il_jj)-threesixty
1343	endif
1344	enddo
1345	enddo
1346	enddo
1347	!rv,sgi>
1348	!
1349	!Dump it to the output file.
1350	!
1351	CALL mpp_write(ig_units(il_port),crnlonij(il_port),loncrn)
1352	DEALLOCATE(rl_lonscrn)
1353	DEALLOCATE(loncrn)
1354	ENDIF
1355	ENDDO
1356
1357	IF (ll_notfound) WRITE(nulprt,*)'psmile_def_metdata: No corner longitudes found' &
1358	,'for grid', cga_clim_locator(il_port)
1359
1360	!rv,sgi DEALLOCATE(lon)
1361	!
1362	!Latitude corners
1363	!
1364	!See comments above.
1365	!
1366	clstrg = cga_clim_locator(il_port)//crn_clim_latsuf
1367	ll_notfound = .TRUE.
1368
1369	DO il_i=1,il_nvar
1370
1371
1372	CALL mpp_get_atts(il_vars(il_i),name=cl_varname,siz=lens_lat)
1373	IF(TRIM(cl_varname).EQ.TRIM(clstrg) ) THEN
1374	ll_notfound = .FALSE.
1375	#ifdef __VERBOSE
1376	write(nulprt,*) &
1377	'psmile_def_metada, lat, clstrg: ',clstrg,lens_lat
1378	#endif
1379	nbr_corners(il_port)=lens_lat(3)
1380	ALLOCATE(rl_latscrn(lens_lat(1),lens_lat(2) &
1381	,nbr_corners(il_port)))
1382	ALLOCATE(latcrn(nbr_corners(il_port) &
1383	,1:il_xend-il_xbegin+1 &
1384	,1:il_yend-il_ybegin+1))
1385	lens_lat(1)=il_xend-il_xbegin+1
1386	lens_lat(2)=il_yend-il_ybegin+1
1387	lens_lat(3)=nbr_corners(il_port)
1388
1389	CALL mpp_read(nc_grdid,il_vars(il_i),rl_latscrn,tindex=1)
1390	!rv,sgi<
1391	if((ig_def_invert(il_port).gt.1).or. &
1392	(ig_def_reverse(il_port).gt.1)) then
1393	llxlon=ig_def_invert(il_port).eq.3 &
1394	.or.ig_def_invert(il_port).eq.5 &
1395	.or.ig_def_reverse(il_port).eq.3 &
1396	.or.ig_def_reverse(il_port).eq.5
1397	llylat=ig_def_invert(il_port).eq.2 &
1398	.or.ig_def_invert(il_port).eq.5 &
1399	.or.ig_def_reverse(il_port).eq.2 &
1400	.or.ig_def_reverse(il_port).eq.5
1401	call inv2dcrn(rl_latscrn,lens_lat(1),lens_lat(2) &
1402	,nbr_corners(il_port),llxlon,llylat)
1403	endif
1404	!rv,sgi>
1405
1406	call transp_crn2cf(rl_latscrn,latcrn,lens_lat)
1407	call mpp_write(ig_units(il_port),crnlatij(il_port),latcrn)
1408	DEALLOCATE(rl_latscrn)
1409	DEALLOCATE(latcrn)
1410	ENDIF
1411	ENDDO
1412	IF (ll_notfound) WRITE(nulprt,*)'psmile_def_metdata: No corner latitudes found' &
1413	,'for grid', cga_clim_locator(il_port)
1414	DEALLOCATE(lat)
1415
1416	!rv,sgi<
1417	if(ALLOCATED(lon)) DEALLOCATE(lon)
1418	!rv,sgi>
1419	DEALLOCATE(il_vars)
1420	ENDIF
1421
1422	CALL mpp_flush(ig_units(il_port))
1423
1424	ENDIF
1425
1426	ENDIF
1427	#ifdef __VERBOSE
1428	write(nulprt,*)'End -- psmile_def_metada'
1429	call flush(nulprt)
1430	#endif
1431	ENDDO
1432	end subroutine psmile_def_metadata
1433
1434	subroutine psmile_close_files(id_error)
1435	!------------------------------------------------------------------------------
1436	! This subroutine closes all opened psmile_io files
1437	! Supposed to be called in prism_terminate before psmile_io_cleanup
1438	! and before shutting down MPI!
1439	!-----------------------------------------------------------------------
1440	!
1441	use mod_kinds_model
1442	use mod_prism_proto
1443	use mod_comprism_proto
1444	use mod_psmile_io
1445	implicit none
1446
1447	#include <mpif.h>
1448	integer(kind=ip_intwp_p),intent(out)::id_error
1449	!Local variables
1450	integer(kind=ip_intwp_p)::il_port
1451
1452	if( ig_color.ne.1000*ig_mynummod) return
1453	#ifdef __VERBOSE
1454	write(nulprt,*)'psmile_close_files: Start'
1455	#endif
1456	id_error=CLIM_Ok
1457
1458	do il_port=1,nports
1459	! IF (ig_def_state(il_port) .eq. ip_ignout .or. &
1460	! ig_def_state(il_port) .eq. ip_expout .or. &
1461	! ig_def_state(il_port) .eq. ip_output .or. &
1462	! ig_def_state(il_port) .eq. ip_input ) THEN
1463	!RV:22.01.2003
1464	!Just in case that ig_def_state has been deallocated during termination.
1465	!This way has less interference with rest of prism_terminate.
1466	IF (ig_units(il_port).gt.0) THEN
1467	IF (allocated(dg_times).and. &
1468	associated(dg_times(il_port)%rl_times)) THEN
1469
1470	deallocate(dg_times(il_port)%rl_times)
1471
1472	ENDIF
1473
1474	call mpp_close(ig_units(il_port))
1475
1476	ENDIF
1477	enddo
1478
1479	!RV: 18.12.2002: Moved these lines into psmile_io_cleanup
1480	! call mpp_io_exit()
1481	! call mpp_domains_exit()
1482	! call mpp_exit()
1483
1484	#ifdef __VERBOSE
1485	write(nulprt,*)'psmile_close_files: End'
1486	#endif
1487
1488	end subroutine psmile_close_files
1489
1490	subroutine psmile_io_cleanup(id_error)
1491	!-----------------------------------------------------------------------
1492	! Supposed to be called in prism_terminate after psmile_close_files
1493	!-----------------------------------------------------------------------
1494	use mod_kinds_model
1495	use mod_prism_proto
1496	use mod_comprism_proto
1497	use mod_psmile_io
1498
1499	implicit none
1500
1501	integer(kind=ip_intwp_p),intent(out)::id_error
1502
1503	if( ig_color.ne.1000*ig_mynummod) return
1504	!Global arrays to specify IO domains
1505	if(allocated(domain))deallocate(domain)
1506	if(allocated(xdom))deallocate(xdom)
1507	if(allocated(ydom))deallocate(ydom)
1508	if(allocated(lg_is_odd_distribution))deallocate(lg_is_odd_distribution)
1509
1510	!Global arrarys for opening files)
1511	if(allocated(ig_units))deallocate(ig_units)
1512	if(allocated(ig_action))deallocate(ig_action)
1513	if(allocated(ig_form))deallocate(ig_form)
1514	if(allocated(ig_threading))deallocate(ig_threading)
1515	if(allocated(ig_fileset))deallocate(ig_fileset)
1516
1517	!Global arrays for netcdf metadata informations)
1518	if(allocated(ig_var_index))deallocate(ig_var_index)
1519	if(allocated(ig_ntimes))deallocate(ig_ntimes)
1520	if(allocated(dg_times))deallocate( dg_times)
1521	if(allocated(ig_vars))deallocate(ig_vars)
1522	if(allocated(field))deallocate(field)
1523	if(allocated(latij))deallocate(latij)
1524	if(allocated(lonij))deallocate(lonij)
1525	if(allocated(x_axis))deallocate(x_axis)
1526	if(allocated(y_axis))deallocate(y_axis)
1527	if(allocated(t_axis))deallocate(t_axis)
1528
1529	!Was not allocated. Oasis sends 2D arrays. In case of bundles or 3d arrays
1530	!we have to declare z_axis + plus one for bundles.
1531	if(allocated(z_axis))deallocate(z_axis)
1532
1533	! These arrays are supposed to be part of a different module!)
1534	! However, no information is available in Oasis 3.0. So, let them live here)
1535	! for a while.)
1536
1537	if(allocated(cports_lgname)) deallocate(cports_lgname)
1538	if(allocated(cports_units)) deallocate(cports_units)
1539
1540	call mpp_io_exit()
1541	call mpp_domains_exit()
1542	call mpp_exit()
1543
1544
1545	end subroutine psmile_io_cleanup
1546
1547	subroutine psmile_read_8(id_port_id,rd_field,id_newtime)
1548	!--------------------------------------------------------------------------
1549	! Called in rism_get
1550	use mod_kinds_model
1551	use mod_prism_proto
1552	use mod_comprism_proto
1553	use mod_psmile_io
1554
1555	implicit none
1556
1557	integer(kind=ip_intwp_p),intent(in)::id_newtime,id_port_id
1558	REAL(kind=ip_double_p),DIMENSION(myport(4,id_port_id)), &
1559	intent(out) :: rd_field
1560	!local
1561	real(kind=ip_double_p)::rl_time
1562	integer(kind=ip_intwp_p)::il_record_no
1563	real(kind=ip_double_p)::diff_first,diff_i
1564	integer(kind=ip_intwp_p)::il_i,i
1565
1566	#ifdef __VERBOSE
1567	write(nulprt,*) 'psmile_read_8: Start'
1568	call flush(nulprt)
1569	#endif
1570
1571	!RV:22.01.2003
1572
1573	#ifdef __READ_BY_LOOKUP
1574	!I try to find the record number of timestamp id_newtime by scrolling through
1575	!the time stamps of the file.
1576
1577	rl_time=id_newtime
1578	if(.not.associated(dg_times(id_port_id)%rl_times))then
1579
1580	allocate(dg_times(id_port_id)%rl_times(1:ig_ntimes(id_port_id)))
1581	call mpp_get_times( ig_units(id_port_id) &
1582	,dg_times(id_port_id)%rl_times(:))
1583	endif
1584	! il_record_no=minloc(abs(dg_times(id_port_id)%rl_times-rl_time),dim=1)
1585	il_record_no=1
1586	diff_first=abs(dg_times(id_port_id)%rl_times(1)-rl_time)
1587	do i=1,ig_ntimes(id_port_id)
1588	diff_i=abs(dg_times(id_port_id)%rl_times(i)-rl_time)
1589	if(abs(diff_i).lt.diff_first) then
1590	diff_first=diff_i
1591	il_record_no=i
1592	endif
1593	enddo
1594	if(abs(dg_times(id_port_id)%rl_times(il_record_no)-rl_time).gt.1.d-8) &
1595	then
1596	write(nulprt,*)'psmile_read_8: The time stamp ',id_newtime, &
1597	' is not available in file ',cg_def_inpfile(id_port_id)
1598	call flush(nulprt)
1599	call MPI_ABORT(mpi_comm,0,mpi_err)
1600	endif
1601
1602	#ifdef __VERBOSE
1603	write(nulprt,*)'psmile_read_8: timestamp, record no.',id_newtime, &
1604	il_record_no
1605	call flush(nulprt)
1606	#endif
1607	#else
1608	!Alternative:
1609	!Assuming that time stamps are starting from zero
1610	il_record_no=id_newtime/ig_def_freq(id_port_id) +1
1611	#endif
1612
1613	call mpp_read(ig_units(id_port_id),ig_vars(id_port_id) &
1614	,domain(id_port_id),rd_field,il_record_no)
1615
1616	#ifdef __VERBOSE
1617	write(nulprt,*) 'psmile_read_8: End'
1618	call flush(nulprt)
1619	#endif
1620	end subroutine psmile_read_8
1621
1622	subroutine psmile_read_4(id_port_id,rd_field,id_newtime)
1623	!--------------------------------------------------------------------------
1624	! Called in prism_get
1625	use mod_kinds_model
1626	use mod_prism_proto
1627	use mod_comprism_proto
1628	use mod_psmile_io
1629
1630	implicit none
1631
1632	integer(kind=ip_intwp_p),intent(in)::id_newtime,id_port_id
1633	REAL(kind=ip_single_p), DIMENSION(myport(4,id_port_id)) :: rd_field
1634	!local
1635	REAL(kind=ip_double_p), DIMENSION(myport(4,id_port_id)) :: rl_field
1636	real(kind=ip_double_p)::rl_time
1637	integer(kind=ip_intwp_p)::il_record_no,il_i
1638
1639	#ifdef __VERBOSE
1640	write(nulprt,*) 'psmile_read_4: Start'
1641	#endif
1642
1643	!RV:22.01.2003
1644	#ifdef __READ_BY_LOOKUP
1645	!I try to find the record number of timestamp id_newtime by scrolling through
1646	!the time stamps of the file.
1647	rl_time=id_newtime
1648
1649	if(.not.associated(dg_times(id_port_id)%rl_times))then
1650
1651	allocate(dg_times(id_port_id)%rl_times(1:ig_ntimes(id_port_id)))
1652	call mpp_get_times( ig_units(id_port_id), &
1653	dg_times(id_port_id)%rl_times(:))
1654
1655	endif
1656
1657	il_record_no=minloc(abs(dg_times(id_port_id)%rl_times-rl_time),dim=1)
1658
1659	if(abs(dg_times(id_port_id)%rl_times(il_record_no)-rl_time).gt.1.d-8) &
1660	then
1661	write(nulprt,*)'psmile_read_4: The time stamp ',id_newtime, &
1662	' is not available in file ',cg_def_inpfile(id_port_id)
1663	call flush(nulprt)
1664	call MPI_ABORT(mpi_comm,0,mpi_err)
1665	endif
1666	#ifdef __VERBOSE
1667	write(nulprt,*)'psmile_read_4: timestamp, record no.',id_newtime, &
1668	il_record_no
1669	call flush(nulprt)
1670	#endif
1671	#else
1672	!Alternative:
1673	!Assuming that time stamps are starting from zero
1674	il_record_no=id_newtime/ig_def_freq(id_port_id) +1
1675	#endif
1676
1677	call mpp_read(ig_units(id_port_id),ig_vars(id_port_id) &
1678	,domain(id_port_id),rl_field,il_record_no)
1679	rd_field=rl_field
1680
1681	#ifdef __VERBOSE
1682	write(nulprt,*) 'psmile_read_4: End'
1683	call flush(nulprt)
1684	#endif
1685	end subroutine psmile_read_4
1686
1687	subroutine psmile_write_4(id_port_id,rd_field,id_newtime)
1688	!--------------------------------------------------------------------------
1689	! Called in prism_get or prism_put
1690	use mod_kinds_model
1691	use mod_prism_proto
1692	use mod_comprism_proto
1693	use mod_psmile_io
1694
1695	implicit none
1696	integer(kind=ip_intwp_p),intent(in)::id_newtime,id_port_id
1697	REAL(kind=ip_single_p), DIMENSION(myport(4,id_port_id)) :: rd_field
1698	!Local
1699	real(kind=ip_double_p)::rl_time
1700	real(ip_double_p),DIMENSION(myport(4,id_port_id)) ::rl_field
1701
1702	#ifdef __VERBOSE
1703	write(nulprt,*)' psmile_write_4: Start'
1704	#endif
1705	rl_time=id_newtime
1706	rl_field=rd_field
1707	call mpp_write(ig_units(id_port_id),field(id_port_id) &
1708	,domain(id_port_id),rl_field,rl_time)
1709	#ifdef __VERBOSE
1710	write(nulprt,*)' psmile_write_4: End'
1711	#endif
1712
1713	end subroutine psmile_write_4
1714
1715	subroutine psmile_write_8(id_port_id,rd_field,id_newtime)
1716	!--------------------------------------------------------------------------
1717	! Called in prism_get or prism_put
1718	use mod_kinds_model
1719	use mod_prism_proto
1720	use mod_comprism_proto
1721	use mod_psmile_io
1722
1723	implicit none
1724	integer(kind=ip_intwp_p),intent(in)::id_newtime,id_port_id
1725	integer(kind=ip_intwp_p)::il_stat
1726	REAL(kind=ip_double_p), &
1727	DIMENSION(1:myport(4,id_port_id)),intent(inout) :: rd_field
1728	!Local
1729	real(kind=ip_double_p)::rl_time
1730	integer(kind=ip_intwp_p) :: ij
1731
1732	! write(70+mpi_rank+(ig_color/1000-1)*4,'(4e12.5)')rd_field
1733	#ifdef __VERBOSE
1734	write(nulprt,*)' psmile_write_8: Start'
1735	call flush(nulprt)
1736	#endif
1737	rl_time=id_newtime
1738	#ifdef __VERBOSE
1739	write(nulprt,*)'rl_time ',rl_time
1740	#endif
1741	! do ij =1,4608,4608/10
1742	! write(nulprt,*)'rd_field bef mpp_write ',rd_field(ij)
1743	! enddo
1744	! call flush(nulprt)
1745	call mpp_write(ig_units(id_port_id),field(id_port_id) &
1746	,domain(id_port_id),rd_field,rl_time)
1747
1748	#ifdef __VERBOSE
1749	write(nulprt,*)' psmile_write_8: End'
1750	call flush(nulprt)
1751	#endif
1752	end subroutine psmile_write_8
1753
1754	subroutine indexi(n,arr,indx)
1755	!---------------------------------------------------------------------
1756	! Generates a list indx which sorts arr in ascending order.
1757	! This is called sorting by indexing.
1758	! Code taken from 'Numerical Recipes'
1759	!---------------------------------------------------------------------
1760	use mod_kinds_model
1761	INTEGER(KIND=IP_INTWP_P),intent(in):: n
1762	INTEGER(KIND=IP_INTWP_P),intent(out)::indx(n)
1763	integer(kind=ip_intwp_p),intent(in):: arr(n)
1764	!Local
1765	integer(kind=ip_intwp_p),PARAMETER:: M=7,NSTACK=128
1766	INTEGER(KIND=IP_INTWP_P):: i,indxt,ir,itemp,j,jstack,k,l,istack(NSTACK)
1767	integer(kind=ip_intwp_p):: a
1768
1769	do 11 j=1,n
1770	indx(j)=j
1771	11 continue
1772	jstack=0
1773	l=1
1774	ir=n
1775	1 if(ir-l.lt.M)then
1776	do 13 j=l+1,ir
1777	indxt=indx(j)
1778	a=arr(indxt)
1779	do 12 i=j-1,1,-1
1780	if(arr(indx(i)).le.a)goto 2
1781	indx(i+1)=indx(i)
1782	12 continue
1783	i=0
1784	2 indx(i+1)=indxt
1785	13 continue
1786	if(jstack.eq.0)return
1787	ir=istack(jstack)
1788	l=istack(jstack-1)
1789	jstack=jstack-2
1790	else
1791	k=(l+ir)/2
1792	itemp=indx(k)
1793	indx(k)=indx(l+1)
1794	indx(l+1)=itemp
1795	if(arr(indx(l+1)).gt.arr(indx(ir)))then
1796	itemp=indx(l+1)
1797	indx(l+1)=indx(ir)
1798	indx(ir)=itemp
1799	endif
1800	if(arr(indx(l)).gt.arr(indx(ir)))then
1801	itemp=indx(l)
1802	indx(l)=indx(ir)
1803	indx(ir)=itemp
1804	endif
1805	if(arr(indx(l+1)).gt.arr(indx(l)))then
1806	itemp=indx(l+1)
1807	indx(l+1)=indx(l)
1808	indx(l)=itemp
1809	endif
1810	i=l+1
1811	j=ir
1812	indxt=indx(l)
1813	a=arr(indxt)
1814	3 continue
1815	i=i+1
1816	if(arr(indx(i)).lt.a)goto 3
1817	4 continue
1818	j=j-1
1819	if(arr(indx(j)).gt.a)goto 4
1820	if(j.lt.i)goto 5
1821	itemp=indx(i)
1822	indx(i)=indx(j)
1823	indx(j)=itemp
1824	goto 3
1825	5 indx(l)=indx(j)
1826	indx(j)=indxt
1827	jstack=jstack+2
1828	if(jstack.gt.NSTACK)pause 'NSTACK too small in indexx'
1829	if(ir-i+1.ge.j-l)then
1830	istack(jstack)=ir
1831	istack(jstack-1)=i
1832	ir=j-1
1833	else
1834	istack(jstack)=j-1
1835	istack(jstack-1)=l
1836	l=i
1837	endif
1838	endif
1839	goto 1
1840	END subroutine indexi
1841
1842	subroutine combine_with_date(cd_in,cd_mode,id_initial_date,cd_on)
1843	!-------------------------------------------------------------
1844	use mod_kinds_model
1845	character(len=*),intent(in)::cd_in
1846	character(len=*),intent(in)::cd_mode
1847	character(len=*),intent(out)::cd_on
1848	integer(kind=ip_intwp_p),intent(in)::id_initial_date(:)
1849	!Local
1850	integer(kind=ip_intwp_p),parameter::il_nsuffixes=2
1851	integer(kind=ip_intwp_p)::leni,leno,lens,lenb,il_i
1852	character(len=8)::suffixes(il_nsuffixes),suffix
1853	character(len=len(cd_in))::cl_basename
1854	character(len=4)::cl_year
1855	character(len=2)::cl_date(1:size(id_initial_date)-1)
1856
1857	suffixes(1)='.nc'
1858	suffixes(2)='.grib'
1859
1860	leni=len(trim(cd_in))
1861	leno=len(cd_on)
1862
1863	do il_i=1,il_nsuffixes
1864
1865	suffix=suffixes(il_i)
1866	lens=len(trim(suffix))
1867
1868	if(cd_in(leni-lens+1:leni).eq.suffix(1:lens)) then
1869	cl_basename(1:leni-lens)=cd_in(1:leni-lens)
1870	lenb=leni-lens
1871	EXIT
1872	else
1873	cl_basename(1:leni)=cd_in(1:leni)
1874	lenb=leni
1875	suffix='.nc'
1876	endif
1877
1878	enddo
1879
1880
1881	write(cl_year,'(i4.4)')id_initial_date(1)
1882
1883	do il_i=2,size(id_initial_date)
1884
1885	write(cl_date(il_i-1),'(i2.2)')id_initial_date(il_i)
1886
1887	enddo
1888
1889	if(trim(cd_mode) .eq. 'since' ) then
1890	cd_on=trim(cl_basename(1:lenb))//' '//trim(cd_mode)//' ' &
1891	//cl_year//'-'//cl_date(1)//'-' &
1892	//cl_date(2)//' ' &
1893	//cl_date(3)//':' &
1894	//cl_date(4)//':' &
1895	//cl_date(5)
1896	else
1897	cd_on=trim(cl_basename(1:lenb))//'_'//trim(cd_mode)//'.' &
1898	//cl_year//'-'//cl_date(1)//'-' &
1899	//cl_date(2)//'T' &
1900	//cl_date(3)//':' &
1901	//cl_date(4)//':' &
1902	//cl_date(5)//suffix
1903	endif
1904
1905	return
1906	end subroutine combine_with_date
1907
1908	subroutine transp_crn2cf(x,y,vshape)
1909	use mod_kinds_model
1910	implicit none
1911	! Input parameters
1912	Integer, Intent(in):: vshape(3)
1913	REAL(kind=ip_realwp_p),Intent(in)::x(vshape(1),vshape(2),vshape(3))
1914	! Ouput parameters
1915	REAL(kind=ip_realwp_p),Intent(out)::y(vshape(3),vshape(1),vshape(2))
1916	!Local
1917	integer::il_i,il_j,il_k
1918
1919	do il_k=1,vshape(3)
1920	do il_j=1,vshape(2)
1921	do il_i=1,vshape(1)
1922	y(il_k, il_i, il_j)=x(il_i ,il_j, il_k)
1923	enddo
1924	enddo
1925	enddo
1926
1927	end subroutine transp_crn2cf
1928
1929	!rv,sgi<
1930	subroutine inv2d(x,ilon,ilat,lxlon,lylat)
1931	use mod_kinds_model
1932	use mod_comprism_proto
1933	implicit none
1934	!
1935	! Input Parameters
1936	!
1937	Integer,Intent(In)::ilon,ilat
1938	REAL(kind=ip_realwp_p)::x(ilon,ilat)
1939	Logical,Intent(In)::lxlon,lylat
1940	!
1941	! Local Variables
1942	!
1943	REAL(kind=ip_realwp_p)::xtmp
1944	Integer(kind=ip_realwp_p):: ijmed,iimed,ji,jj
1945	#ifdef __VERBOSE
1946	write(nulprt,*)'inv2d start:',ilon,ilat,lxlon,lylat
1947	call flush(nulprt)
1948	#endif
1949	!
1950	! Reorder lattitude related values stored as the y-direction.
1951	!
1952	if(lylat) then
1953	ijmed = ilat/2
1954	DO jj = 1, ijmed
1955	DO ji = 1, ilon
1956	xtmp = x(ji,ilat + 1 - jj)
1957	x(ji,ilat + 1 - jj) = x(ji,jj)
1958	x(ji,jj) = xtmp
1959	ENDDO
1960	ENDDO
1961	endif
1962
1963	!
1964	! Reorder longitude related values stored as the x-direction.
1965	!
1966	if(lxlon) then
1967	iimed = ilon/2
1968	DO jj = 1, ilat
1969	DO ji = 1, iimed
1970	xtmp = x(ilon + 1 - ji,jj)
1971	x(ilon + 1 - ji,jj) = x(ji,jj)
1972	x(ji,jj) = xtmp
1973	ENDDO
1974	ENDDO
1975	endif
1976
1977	#ifdef __VERBOSE
1978	write(nulprt,*)'inv2d end'
1979	call flush(nulprt)
1980	#endif
1981
1982	end subroutine inv2d
1983
1984	subroutine inv2dcrn(x,ilon,ilat,icrn,lxlon,lylat)
1985	use mod_kinds_model
1986	use mod_comprism_proto
1987	implicit none
1988	!
1989	! Input Parameters
1990	!
1991	Integer,Intent(In)::ilon,ilat,icrn
1992	REAL(kind=ip_realwp_p)::x(ilon,ilat,icrn)
1993	Logical,Intent(In)::lxlon,lylat
1994	!
1995	! Local Variables
1996	!
1997	REAL(kind=ip_realwp_p)::xtmp
1998	Integer(kind=ip_realwp_p):: ijmed,iimed,ji,jj,jc
1999	#ifdef __VERBOSE
2000	write(nulprt,*)'inv2dcrn start:',ilon,ilat,icrn,lxlon,lylat
2001	call flush(nulprt)
2002	#endif
2003	!
2004	! Reorder corners of lattitude
2005	!
2006	if(lylat) then
2007	ijmed = ilat/2
2008	DO jc=1,icrn
2009	DO jj = 1, ijmed
2010	DO ji = 1, ilon
2011	xtmp = x(ji,ilat + 1 - jj,jc)
2012	x(ji,ilat + 1 - jj,jc) = x(ji,jj,jc)
2013	x(ji,jj,jc) = xtmp
2014	ENDDO
2015	ENDDO
2016	ENDDO
2017	endif
2018
2019	!
2020	! Reorder corners of longitude
2021	!
2022	if(lxlon) then
2023	iimed = ilon/2
2024	DO jc=1,icrn
2025	DO jj = 1, ilat
2026	DO ji = 1, iimed
2027	xtmp = x(ilon + 1 - ji,jj,jc)
2028	x(ilon + 1 - ji,jj,jc) = x(ji,jj,jc)
2029	x(ji,jj,jc) = xtmp
2030	ENDDO
2031	ENDDO
2032	ENDDO
2033	endif
2034
2035
2036	#ifdef __VERBOSE
2037	write(nulprt,*)'inv2dcrn end'
2038	call flush(nulprt)
2039	#endif
2040	end subroutine inv2dcrn
2041	!rv,sgi>

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

Context Navigation

source: CPL/oasis3/trunk/src/lib/psmile/src/mod_psmile_io.F90 @ 1677

Download in other formats: