Context Navigation

source: codes/icosagcm/trunk/src/write_field.f90 @ 107

Last change on this file since 107 was 107, checked in by ymipsl, 12 years ago

Memory leak fix in writefield

YM

File size: 65.2 KB

Line
1	module write_field
2	USE genmod
3	implicit none
4
5	INTEGER,SAVE :: ncprec
6
7	TYPE ncvar
8	INTEGER :: size
9	INTEGER,POINTER :: nc_id(:)
10	INTEGER :: displ
11	END TYPE ncvar
12
13	integer, parameter :: MaxWriteField = 1000
14	integer, dimension(MaxWriteField),save :: FieldId
15	TYPE(ncvar), dimension(MaxWriteField),save :: FieldVarId
16	integer, dimension(MaxWriteField),save :: FieldIndex
17	character(len=255), dimension(MaxWriteField) :: FieldName
18
19	integer,save :: NbField = 0
20
21	contains
22
23	SUBROUTINE init_writeField
24	USE ioipsl
25	use netcdf_mod
26	IMPLICIT NONE
27	CHARACTER(LEN=255) :: netcdf_prec
28
29	netcdf_prec='float'
30	CALL getin("netcdf_prec",netcdf_prec)
31
32	SELECT CASE(TRIM(netcdf_prec))
33	CASE('float')
34	ncprec=NF90_FLOAT
35	CASE('double')
36	ncprec=NF90_DOUBLE
37	CASE default
38	PRINT*,'Bad selector for variable netcdf_prec : <', TRIM(netcdf_prec),"> options are <float>, <double>"
39	STOP
40	END SELECT
41	END SUBROUTINE init_writeField
42
43	function GetFieldIndex(name)
44	implicit none
45	integer :: GetFieldindex
46	character(len=*) :: name
47
48	character(len=255) :: TrueName
49	integer :: i
50
51
52	TrueName=TRIM(ADJUSTL(name))
53
54	GetFieldIndex=-1
55	do i=1,NbField
56	if (TrueName==FieldName(i)) then
57	GetFieldIndex=i
58	exit
59	endif
60	enddo
61	end function GetFieldIndex
62
63	SUBROUTINE Writefield(name_in,field,nind,once)
64	USE domain_mod
65	USE field_mod
66	USE transfert_mpi_mod
67	USE dimensions
68	USE mpipara
69	IMPLICIT NONE
70	CHARACTER(LEN=*),INTENT(IN) :: name_in
71	TYPE(t_field),POINTER :: field(:)
72	INTEGER,OPTIONAL,INTENT(IN) :: nind
73	LOGICAL,OPTIONAL,INTENT(IN) :: once
74
75	TYPE(t_field),POINTER :: field_glo(:)
76
77	IF (field(1)%ndim==2) THEN
78	CALL allocate_field_glo(field_glo,field(1)%field_type,field(1)%data_type)
79	ELSE IF (field(1)%ndim==3) THEN
80	CALL allocate_field_glo(field_glo,field(1)%field_type,field(1)%data_type,field(1)%dim3)
81	ELSE IF (field(1)%ndim==4) THEN
82	CALL allocate_field_glo(field_glo,field(1)%field_type,field(1)%data_type,field(1)%dim3,field(1)%dim4)
83	ENDIF
84
85	CALL gather_field(field,field_glo)
86
87	IF (mpi_rank==0) THEN
88	IF (PRESENT(nind)) THEN
89	CALL writefield_gen(name_in,field_glo,domain_glo,nind,once=once)
90	ELSE
91	CALL writefield_gen(name_in,field_glo,domain_glo,1,ndomain_glo,once=once)
92	ENDIF
93	ENDIF
94
95	CALL deallocate_field_glo(field_glo)
96
97	END SUBROUTINE writefield
98
99	! SUBROUTINE Writefield(name_in,field,nind)
100	! USE netcdf
101	! USE domain_mod
102	! use field_mod
103	! USE dimensions
104	! USE geometry
105	! IMPLICIT NONE
106	! CHARACTER(LEN=*),INTENT(IN) :: name_in
107	! TYPE(t_field),POINTER :: field(:)
108	! INTEGER,OPTIONAL,INTENT(IN) :: nind
109	! REAL(r8),ALLOCATABLE :: field_val2d(:)
110	! REAL(r8),ALLOCATABLE :: field_val3d(:,:)
111	! REAL(r8),ALLOCATABLE :: field_val4d(:,:,:)
112	! TYPE(t_domain),POINTER :: d
113	! INTEGER :: Index
114	! INTEGER :: ind,i,j,k,n,ncell,q
115	! INTEGER :: iie,jje,iin,jjn
116	! INTEGER :: status
117	! CHARACTER(len=255) :: name
118	! CHARACTER(len=255) :: str_ind
119	! INTEGER :: ind_b,ind_e
120	! INTEGER :: halo_size
121	! LOGICAL :: single
122	!
123	!
124	! name=TRIM(ADJUSTL(name_in))
125
126	! IF (PRESENT(nind)) THEN
127	! name=TRIM(name)//"_"//TRIM(int2str(nind))
128	! PRINT *,"NAME",nind,int2str(nind),name
129	! ind_b=nind
130	! ind_e=nind
131	! halo_size=1
132	! single=.TRUE.
133	! ELSE
134	! ind_b=1
135	! ind_e=ndomain
136	! halo_size=0
137	! single=.FALSE.
138	! ENDIF
139
140	! Index=GetFieldIndex(name)
141	! if (Index==-1) then
142	! call create_header(name,field,nind)
143	! Index=GetFieldIndex(name)
144	! else
145	! FieldIndex(Index)=FieldIndex(Index)+1.
146	! endif
147	!
148	! IF (Field(ind_b)%field_type==field_T) THEN
149	! ncell=1
150	! DO ind=ind_b,ind_e
151	! d=>domain(ind)
152	! IF (Field(ind)%field_type/=field_T) THEN
153	! PRINT *,"Writefield, grille non geree"
154	! RETURN
155	! ENDIF
156
157	! n=0
158	! DO j=d%jj_begin-halo_size,d%jj_end+halo_size
159	! DO i=d%ii_begin-halo_size,d%ii_end+halo_size
160	! IF (d%own(i,j) .OR. single) n=n+1
161	! ENDDO
162	! ENDDO
163
164	! IF (field(ind)%ndim==2) THEN
165	! ALLOCATE(Field_val2d(n))
166	! n=0
167	! DO j=d%jj_begin-halo_size,d%jj_end+halo_size
168	! DO i=d%ii_begin-halo_size,d%ii_end+halo_size
169	! k=d%iim*(j-1)+i
170	! IF (d%own(i,j) .OR. single) THEN
171	! n=n+1
172	! Field_val2d(n)=field(ind)%rval2d(k)
173	! ENDIF
174	! ENDDO
175	! ENDDO
176	! status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(1),Field_val2d, &
177	! start=(/ ncell,FieldIndex(Index) /),count=(/n,1 /))
178	! DEALLOCATE(field_val2d)
179	! ELSE IF (field(ind)%ndim==3) THEN
180	! ALLOCATE(Field_val3d(n,size(field(ind)%rval3d,2)))
181	! n=0
182	! DO j=d%jj_begin-halo_size,d%jj_end+halo_size
183	! DO i=d%ii_begin-halo_size,d%ii_end+halo_size
184	! k=d%iim*(j-1)+i
185	! IF (d%own(i,j) .OR. single) THEN
186	! n=n+1
187	! Field_val3d(n,:)=field(ind)%rval3d(k,:)
188	! ENDIF
189	! ENDDO
190	! ENDDO
191	! status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(1),Field_val3d,start=(/ ncell,1,FieldIndex(Index) /), &
192	! count=(/n,size(field(1)%rval3d,2),1 /))
193	! DEALLOCATE(field_val3d)
194	! ELSE IF (field(1)%ndim==4) THEN
195
196	! DO q=1,FieldVarId(index)%size
197	!
198	! ALLOCATE(Field_val3d(n,size(field(ind)%rval4d,2)))
199	! n=0
200	! DO j=d%jj_begin-halo_size,d%jj_end+halo_size
201	! DO i=d%ii_begin-halo_size,d%ii_end+halo_size
202	! k=d%iim*(j-1)+i
203	! IF (d%own(i,j) .OR. single) THEN
204	! n=n+1
205	! Field_val3d(n,:)=field(ind)%rval4d(k,:,q)
206	! ENDIF
207	! ENDDO
208	! ENDDO
209	! status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(q),Field_val3d,start=(/ ncell,1,FieldIndex(Index) /), &
210	! count=(/n,size(field(1)%rval4d,2),1 /))
211	! DEALLOCATE(field_val3d)
212	! ENDDO
213	! ENDIF
214	!
215	! PRINT *,NF90_STRERROR(status)
216	! ncell=ncell+n
217
218	! ENDDO
219	!
220	! ELSE IF (Field(ind_b)%field_type==field_Z) THEN
221	! ncell=1
222	! n=0
223	! DO ind=ind_b,ind_e
224	! d=>domain(ind)
225	! CALL swap_geometry(ind)
226	! CALL swap_dimensions(ind)
227	!
228	! n=0
229	! DO j=jj_begin+1,jj_end
230	! DO i=ii_begin,ii_end-1
231	! n=n+1
232	! ENDDO
233	! ENDDO
234
235	! DO j=jj_begin,jj_end-1
236	! DO i=ii_begin+1,ii_end
237	! n=n+1
238	! ENDDO
239	! ENDDO
240
241	! IF (field(ind)%ndim==2) THEN
242	! ALLOCATE(Field_val2d(n))
243
244	! n=0
245	! DO j=jj_begin+1,jj_end
246	! DO i=ii_begin,ii_end-1
247	! n=n+1
248	! k=iim*(j-1)+i
249	! Field_val2d(n)=field(ind)%rval2d(k+z_down)
250	! ENDDO
251	! ENDDO
252
253	! DO j=jj_begin,jj_end-1
254	! DO i=ii_begin+1,ii_end
255	! n=n+1
256	! k=iim*(j-1)+i
257	! Field_val2d(n)=field(ind)%rval2d(k+z_up)
258	! ENDDO
259	! ENDDO
260
261	! status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(1), &
262	! Field_val2d,start=(/ ncell,FieldIndex(Index) /),count=(/n,1 /))
263	! DEALLOCATE(field_val2d)
264
265	! ELSE IF (field(ind)%ndim==3) THEN
266	! ALLOCATE(Field_val3d(n,size(field(ind)%rval3d,2)))
267	! n=0
268	! DO j=jj_begin+1,jj_end
269	! DO i=ii_begin,ii_end-1
270	! n=n+1
271	! k=iim*(j-1)+i
272	! Field_val3d(n,:)=field(ind)%rval3d(k+z_down,:)
273	! ENDDO
274	! ENDDO
275
276	! DO j=jj_begin,jj_end-1
277	! DO i=ii_begin+1,ii_end
278	! n=n+1
279	! k=iim*(j-1)+i
280	! Field_val3d(n,:)=field(ind)%rval3d(k+z_up,:)
281	! ENDDO
282	! ENDDO
283	! status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(1),Field_val3d,start=(/ ncell,1,FieldIndex(Index) /), &
284	! count=(/n,size(field(1)%rval3d,2),1 /))
285	! DEALLOCATE(field_val3d)
286	! ELSE IF (field(1)%ndim==4) THEN
287
288	! DO q=1,FieldVarId(index)%size
289	! ALLOCATE(Field_val3d(n,size(field(ind)%rval4d,2)))
290	! n=0
291	! DO j=jj_begin+1,jj_end
292	! DO i=ii_begin,ii_end-1
293	! n=n+1
294	! k=iim*(j-1)+i
295	! Field_val3d(n,:)=field(ind)%rval4d(k+z_down,:,q)
296	! ENDDO
297	! ENDDO
298
299	! DO j=jj_begin,jj_end-1
300	! DO i=ii_begin+1,ii_end
301	! n=n+1
302	! k=iim*(j-1)+i
303	! Field_val3d(n,:)=field(ind)%rval4d(k+z_up,:,q)
304	! ENDDO
305	! ENDDO
306
307	! status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(q),Field_val3d,start=(/ ncell,1,1,FieldIndex(Index) /), &
308	! count=(/n,size(field(1)%rval4d,2),1 /))
309	! DEALLOCATE(field_val3d)
310	! ENDDO
311	! ENDIF
312	!
313	! PRINT *,NF90_STRERROR(status)
314	! ncell=ncell+n
315
316	! ENDDO
317	!
318	! ENDIF
319	! status=NF90_SYNC(FieldId(Index))
320	!
321	! END SUBROUTINE Writefield
322
323
324	SUBROUTINE Writefield_gen(name_in, field, domain_type, ind_b_in, ind_e_in,once )
325	USE netcdf_mod
326	USE domain_mod
327	USE field_mod
328	! USE dimensions
329	! USE geometry
330	IMPLICIT NONE
331	CHARACTER(LEN=*),INTENT(IN) :: name_in
332	TYPE(t_field), POINTER :: field(:)
333	TYPE(t_domain),INTENT(IN),TARGET :: domain_type(:)
334	INTEGER,OPTIONAL,INTENT(IN) :: ind_b_in
335	INTEGER,OPTIONAL,INTENT(IN) :: ind_e_in
336	REAL(r8),ALLOCATABLE :: field_val2d(:)
337	REAL(r8),ALLOCATABLE :: field_val3d(:,:)
338	REAL(r8),ALLOCATABLE :: field_val4d(:,:,:)
339	LOGICAL,OPTIONAL, INTENT(IN) :: once
340	TYPE(t_domain),POINTER :: d
341	INTEGER :: Index
342	INTEGER :: ind,i,j,k,n,ncell,q
343	INTEGER :: iie,jje,iin,jjn
344	INTEGER :: status
345	CHARACTER(len=255) :: name
346	CHARACTER(len=255) :: str_ind
347	INTEGER :: ind_b,ind_e
348	INTEGER :: halo_size
349	LOGICAL :: single
350
351
352	name=TRIM(ADJUSTL(name_in))
353
354	IF (PRESENT(ind_b_in) .AND. .NOT. PRESENT(ind_e_in)) THEN
355	name=TRIM(name)//"_"//TRIM(int2str(ind_b))
356	ind_b=ind_b_in
357	ind_e=ind_b_in
358	halo_size=1
359	single=.TRUE.
360	ELSE IF (.NOT. PRESENT(ind_b_in) .AND. PRESENT(ind_e_in)) THEN
361	name=TRIM(name)//"_"//TRIM(int2str(ind_e))
362	ind_b=ind_e_in
363	ind_e=ind_e_in
364	halo_size=1
365	single=.TRUE.
366	ELSE IF ( PRESENT(ind_b_in) .AND. PRESENT(ind_e_in)) THEN
367	ind_b=ind_b_in
368	ind_e=ind_e_in
369	halo_size=0
370	single=.FALSE.
371	ELSE
372	halo_size=0
373	ind_b=1
374	ind_e=ndomain
375	single=.FALSE.
376	ENDIF
377
378	Index=GetFieldIndex(name)
379	if (Index==-1) then
380	call create_header_gen(name_in,field,domain_type,ind_b_in,ind_e_in,once)
381	Index=GetFieldIndex(name)
382	else
383	FieldIndex(Index)=FieldIndex(Index)+1.
384	endif
385
386	IF (Field(ind_b)%field_type==field_T) THEN
387
388	ncell=0
389	DO ind=ind_b,ind_e
390	d=>domain_type(ind)
391	DO j=d%jj_begin-halo_size,d%jj_end+halo_size
392	DO i=d%ii_begin-halo_size,d%ii_end+halo_size
393	IF (d%assign_domain(i,j)==ind .OR. single) ncell=ncell+1
394	ENDDO
395	ENDDO
396	ENDDO
397
398	IF (field(1)%ndim==2) THEN
399	ALLOCATE(Field_val2d(ncell))
400	n=0
401	DO ind=ind_b,ind_e
402	d=>domain_type(ind)
403	DO j=d%jj_begin-halo_size,d%jj_end+halo_size
404	DO i=d%ii_begin-halo_size,d%ii_end+halo_size
405	k=d%iim*(j-1)+i
406	IF (d%assign_domain(i,j)==ind .OR. single) THEN
407	n=n+1
408	Field_val2d(n)=field(ind)%rval2d(k)
409	ENDIF
410	ENDDO
411	ENDDO
412	ENDDO
413
414	IF (PRESENT(once) .AND. once) THEN
415	status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(1),Field_val2d, &
416	start=(/ 1 /),count=(/ncell /))
417	ELSE
418	status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(1),Field_val2d, &
419	start=(/ 1,FieldIndex(Index) /),count=(/ncell,1 /))
420	ENDIF
421
422	DEALLOCATE(field_val2d)
423	ELSE IF (field(1)%ndim==3) THEN
424	ALLOCATE(Field_val3d(ncell,size(field(1)%rval3d,2)))
425	n=0
426	DO ind=ind_b,ind_e
427	d=>domain_type(ind)
428	DO j=d%jj_begin-halo_size,d%jj_end+halo_size
429	DO i=d%ii_begin-halo_size,d%ii_end+halo_size
430	k=d%iim*(j-1)+i
431	IF (d%assign_domain(i,j)==ind .OR. single) THEN
432	n=n+1
433	Field_val3d(n,:)=field(ind)%rval3d(k,:)
434	ENDIF
435	ENDDO
436	ENDDO
437	ENDDO
438
439	IF (PRESENT(once) .AND. once) THEN
440	status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(1),Field_val3d,start=(/ 1,1 /), &
441	count=(/ncell,size(field(1)%rval3d,2) /))
442	ELSE
443	status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(1),Field_val3d,start=(/ 1,1,FieldIndex(Index) /), &
444	count=(/ncell,size(field(1)%rval3d,2),1 /))
445	ENDIF
446
447	DEALLOCATE(field_val3d)
448	ELSE IF (field(1)%ndim==4) THEN
449
450	DO q=1,FieldVarId(index)%size
451
452	ALLOCATE(Field_val3d(ncell,size(field(1)%rval4d,2)))
453	n=0
454	DO ind=ind_b,ind_e
455	d=>domain_type(ind)
456	DO j=d%jj_begin-halo_size,d%jj_end+halo_size
457	DO i=d%ii_begin-halo_size,d%ii_end+halo_size
458	k=d%iim*(j-1)+i
459	IF (d%assign_domain(i,j)==ind .OR. single) THEN
460	n=n+1
461	Field_val3d(n,:)=field(ind)%rval4d(k,:,q)
462	ENDIF
463	ENDDO
464	ENDDO
465	ENDDO
466
467	IF (PRESENT(once) .AND. once) THEN
468	status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(q),Field_val3d,start=(/ 1,1 /), &
469	count=(/ncell,size(field(1)%rval4d,2),1 /))
470	ELSE
471	status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(q),Field_val3d,start=(/ 1,1,FieldIndex(Index) /), &
472	count=(/ncell,size(field(1)%rval4d,2),1 /))
473	ENDIF
474	DEALLOCATE(field_val3d)
475	ENDDO
476	ENDIF
477
478	! PRINT *,NF90_STRERROR(status)
479	! ncell=ncell+n
480
481	! ENDDO
482
483	ELSE IF (Field(ind_b)%field_type==field_Z) THEN
484
485	ncell=0
486	DO ind=ind_b,ind_e
487	d=>domain_type(ind)
488
489	DO j=d%jj_begin+1,d%jj_end
490	DO i=d%ii_begin,d%ii_end-1
491	ncell=ncell+1
492	ENDDO
493	ENDDO
494
495	DO j=d%jj_begin,d%jj_end-1
496	DO i=d%ii_begin+1,d%ii_end
497	ncell=ncell+1
498	ENDDO
499	ENDDO
500	ENDDO
501
502	IF (field(1)%ndim==2) THEN
503	ALLOCATE(Field_val2d(ncell))
504
505	n=0
506
507	DO ind=ind_b,ind_e
508	d=>domain_type(ind)
509	DO j=d%jj_begin+1,d%jj_end
510	DO i=d%ii_begin,d%ii_end-1
511	n=n+1
512	k=d%iim*(j-1)+i
513	Field_val2d(n)=field(ind)%rval2d(k+d%z_down)
514	ENDDO
515	ENDDO
516
517	DO j=d%jj_begin,d%jj_end-1
518	DO i=d%ii_begin+1,d%ii_end
519	n=n+1
520	k=d%iim*(j-1)+i
521	Field_val2d(n)=field(ind)%rval2d(k+d%z_up)
522	ENDDO
523	ENDDO
524	ENDDO
525
526	IF (PRESENT(once) .AND. once) THEN
527	status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(1), &
528	Field_val2d,start=(/ 1 /),count=(/ncell /))
529	ELSE
530	status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(1), &
531	Field_val2d,start=(/ 1,FieldIndex(Index) /),count=(/ncell,1 /))
532	ENDIF
533	DEALLOCATE(field_val2d)
534
535	ELSE IF (field(1)%ndim==3) THEN
536	ALLOCATE(Field_val3d(ncell,size(field(1)%rval3d,2)))
537	n=0
538	DO ind=ind_b,ind_e
539	d=>domain_type(ind)
540	DO j=d%jj_begin+1,d%jj_end
541	DO i=d%ii_begin,d%ii_end-1
542	n=n+1
543	k=d%iim*(j-1)+i
544	Field_val3d(n,:)=field(ind)%rval3d(k+d%z_down,:)
545	ENDDO
546	ENDDO
547
548	DO j=d%jj_begin,d%jj_end-1
549	DO i=d%ii_begin+1,d%ii_end
550	n=n+1
551	k=d%iim*(j-1)+i
552	Field_val3d(n,:)=field(ind)%rval3d(k+d%z_up,:)
553	ENDDO
554	ENDDO
555	ENDDO
556
557	IF (PRESENT(once) .AND. once) THEN
558	status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(1),Field_val3d,start=(/ 1,1 /), &
559	count=(/ncell,size(field(1)%rval3d,2) /))
560	ELSE
561	status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(1),Field_val3d,start=(/ 1,1,FieldIndex(Index) /), &
562	count=(/ncell,size(field(1)%rval3d,2),1 /))
563	ENDIF
564	DEALLOCATE(field_val3d)
565
566	ELSE IF (field(1)%ndim==4) THEN
567
568	DO q=1,FieldVarId(index)%size
569	ALLOCATE(Field_val3d(ncell,size(field(1)%rval4d,2)))
570	n=0
571	DO ind=ind_b,ind_e
572	d=>domain_type(ind)
573	DO j=d%jj_begin+1,d%jj_end
574	DO i=d%ii_begin,d%ii_end-1
575	n=n+1
576	k=d%iim*(j-1)+i
577	Field_val3d(n,:)=field(ind)%rval4d(k+d%z_down,:,q)
578	ENDDO
579	ENDDO
580
581	DO j=d%jj_begin,d%jj_end-1
582	DO i=d%ii_begin+1,d%ii_end
583	n=n+1
584	k=d%iim*(j-1)+i
585	Field_val3d(n,:)=field(ind)%rval4d(k+d%z_up,:,q)
586	ENDDO
587	ENDDO
588	ENDDO
589
590	IF (PRESENT(once) .AND. once) THEN
591	status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(q),Field_val3d,start=(/ 1,1,1 /), &
592	count=(/ncell,size(field(1)%rval4d,2) /))
593	ELSE
594	status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(q),Field_val3d,start=(/ 1,1,1,FieldIndex(Index) /), &
595	count=(/ncell,size(field(1)%rval4d,2),1 /))
596	ENDIF
597	DEALLOCATE(field_val3d)
598	ENDDO
599	ENDIF
600
601	! PRINT *,NF90_STRERROR(status)
602	! ncell=ncell+n
603	!
604	! ENDDO
605
606	ENDIF
607	status=NF90_SYNC(FieldId(Index))
608
609	END SUBROUTINE Writefield_gen
610
611
612
613	SUBROUTINE Writefield_mpi(name_in,field,nind)
614	USE netcdf_mod
615	USE domain_mod
616	use field_mod
617	USE dimensions
618	USE geometry
619	IMPLICIT NONE
620	CHARACTER(LEN=*),INTENT(IN) :: name_in
621	TYPE(t_field),POINTER :: field(:)
622	INTEGER,OPTIONAL,INTENT(IN) :: nind
623	REAL(r8),ALLOCATABLE :: field_val2d(:)
624	REAL(r8),ALLOCATABLE :: field_val3d(:,:)
625	REAL(r8),ALLOCATABLE :: field_val4d(:,:,:)
626	TYPE(t_domain),POINTER :: d
627	INTEGER :: Index
628	INTEGER :: ind,i,j,l,k,n,ncell,q
629	INTEGER :: iie,jje,iin,jjn
630	INTEGER :: status
631	CHARACTER(len=255) :: name
632	CHARACTER(len=255) :: str_ind
633	INTEGER :: ind_b,ind_e
634	INTEGER :: halo_size
635	LOGICAL :: single
636	INTEGER :: displ
637
638
639	name=TRIM(ADJUSTL(name_in))
640
641	IF (PRESENT(nind)) THEN
642	name=TRIM(name)//"_"//TRIM(int2str(nind))
643	PRINT *,"NAME",nind,int2str(nind),name
644	ind_b=nind
645	ind_e=nind
646	halo_size=1
647	single=.TRUE.
648	ELSE
649	ind_b=1
650	ind_e=ndomain
651	halo_size=0
652	single=.FALSE.
653	ENDIF
654
655	Index=GetFieldIndex(name)
656	if (Index==-1) then
657	call create_header_mpi(name,field,nind)
658	Index=GetFieldIndex(name)
659	else
660	FieldIndex(Index)=FieldIndex(Index)+1.
661	endif
662
663	IF (Field(ind_b)%field_type==field_T) THEN
664	ncell=1
665	DO ind=ind_b,ind_e
666	d=>domain(ind)
667	IF (Field(ind)%field_type/=field_T) THEN
668	PRINT *,"Writefield, grille non geree"
669	RETURN
670	ENDIF
671
672	n=0
673	DO j=d%jj_begin-halo_size,d%jj_end+halo_size
674	DO i=d%ii_begin-halo_size,d%ii_end+halo_size
675	IF (d%own(i,j) .OR. single) n=n+1
676	ENDDO
677	ENDDO
678
679	displ=FieldVarId(index)%displ
680
681	IF (field(ind)%ndim==2) THEN
682	ALLOCATE(Field_val2d(n))
683	n=0
684	DO j=d%jj_begin-halo_size,d%jj_end+halo_size
685	DO i=d%ii_begin-halo_size,d%ii_end+halo_size
686	k=d%iim*(j-1)+i
687	IF (d%own(i,j) .OR. single) THEN
688	n=n+1
689	Field_val2d(n)=field(ind)%rval2d(k)
690	ENDIF
691	ENDDO
692	ENDDO
693	status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(1),Field_val2d, &
694	start=(/ displ+ncell,FieldIndex(Index) /),count=(/n,1 /))
695	DEALLOCATE(field_val2d)
696	ELSE IF (field(ind)%ndim==3) THEN
697	ALLOCATE(Field_val3d(n,size(field(ind)%rval3d,2)))
698	n=0
699	DO j=d%jj_begin-halo_size,d%jj_end+halo_size
700	DO i=d%ii_begin-halo_size,d%ii_end+halo_size
701	k=d%iim*(j-1)+i
702	IF (d%own(i,j) .OR. single) THEN
703	n=n+1
704	Field_val3d(n,:)=field(ind)%rval3d(k,:)
705	ENDIF
706	ENDDO
707	ENDDO
708	! DO l=1,size(field(ind)%rval3d,2)
709	status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(1),Field_val3d, &
710	start=(/ displ+ncell,1,FieldIndex(Index) /), count=(/n,size(field(ind)%rval3d,2),1 /))
711	! ENDDO
712	DEALLOCATE(field_val3d)
713	ELSE IF (field(1)%ndim==4) THEN
714
715	DO q=1,FieldVarId(index)%size
716
717	ALLOCATE(Field_val3d(n,size(field(ind)%rval4d,2)))
718	n=0
719	DO j=d%jj_begin-halo_size,d%jj_end+halo_size
720	DO i=d%ii_begin-halo_size,d%ii_end+halo_size
721	k=d%iim*(j-1)+i
722	IF (d%own(i,j) .OR. single) THEN
723	n=n+1
724	Field_val3d(n,:)=field(ind)%rval4d(k,:,q)
725	ENDIF
726	ENDDO
727	ENDDO
728	! DO l=1,size(field(ind)%rval4d,2)
729	status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(q),Field_val3d(:,l), &
730	start=(/ displ+ncell,l,FieldIndex(Index) /), count=(/n,size(field(ind)%rval4d,2),1 /))
731	! ENDDO
732	DEALLOCATE(field_val3d)
733	ENDDO
734	ENDIF
735
736	! PRINT *,NF90_STRERROR(status)
737	ncell=ncell+n
738
739	ENDDO
740
741	ELSE IF (Field(ind_b)%field_type==field_Z) THEN
742	ncell=1
743	n=0
744	DO ind=ind_b,ind_e
745	d=>domain(ind)
746	CALL swap_geometry(ind)
747	CALL swap_dimensions(ind)
748
749	n=0
750	DO j=jj_begin+1,jj_end
751	DO i=ii_begin,ii_end-1
752	n=n+1
753	ENDDO
754	ENDDO
755
756	DO j=jj_begin,jj_end-1
757	DO i=ii_begin+1,ii_end
758	n=n+1
759	ENDDO
760	ENDDO
761
762	displ=FieldVarId(index)%displ
763
764	IF (field(ind)%ndim==2) THEN
765	ALLOCATE(Field_val2d(n))
766
767	n=0
768	DO j=jj_begin+1,jj_end
769	DO i=ii_begin,ii_end-1
770	n=n+1
771	k=iim*(j-1)+i
772	Field_val2d(n)=field(ind)%rval2d(k+z_down)
773	ENDDO
774	ENDDO
775
776	DO j=jj_begin,jj_end-1
777	DO i=ii_begin+1,ii_end
778	n=n+1
779	k=iim*(j-1)+i
780	Field_val2d(n)=field(ind)%rval2d(k+z_up)
781	ENDDO
782	ENDDO
783
784	status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(1), &
785	Field_val2d,start=(/ displ+ncell,FieldIndex(Index) /),count=(/n,1 /))
786	DEALLOCATE(field_val2d)
787
788	ELSE IF (field(ind)%ndim==3) THEN
789	ALLOCATE(Field_val3d(n,size(field(ind)%rval3d,2)))
790	n=0
791	DO j=jj_begin+1,jj_end
792	DO i=ii_begin,ii_end-1
793	n=n+1
794	k=iim*(j-1)+i
795	Field_val3d(n,:)=field(ind)%rval3d(k+z_down,:)
796	ENDDO
797	ENDDO
798
799	DO j=jj_begin,jj_end-1
800	DO i=ii_begin+1,ii_end
801	n=n+1
802	k=iim*(j-1)+i
803	Field_val3d(n,:)=field(ind)%rval3d(k+z_up,:)
804	ENDDO
805	ENDDO
806	! DO l=1,size(field(ind)%rval3d,2)
807	status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(1),Field_val3d, &
808	start=(/ displ+ncell,1,FieldIndex(Index) /), count=(/n,size(field(ind)%rval3d,2),1 /))
809	! ENDDO
810	DEALLOCATE(field_val3d)
811	ELSE IF (field(1)%ndim==4) THEN
812
813	DO q=1,FieldVarId(index)%size
814	ALLOCATE(Field_val3d(n,size(field(ind)%rval4d,2)))
815	n=0
816	DO j=jj_begin+1,jj_end
817	DO i=ii_begin,ii_end-1
818	n=n+1
819	k=iim*(j-1)+i
820	Field_val3d(n,:)=field(ind)%rval4d(k+z_down,:,q)
821	ENDDO
822	ENDDO
823
824	DO j=jj_begin,jj_end-1
825	DO i=ii_begin+1,ii_end
826	n=n+1
827	k=iim*(j-1)+i
828	Field_val3d(n,:)=field(ind)%rval4d(k+z_up,:,q)
829	ENDDO
830	ENDDO
831
832	! DO l=1,size(field(ind)%rval4d,2)
833
834	status=NF90_PUT_VAR(FieldId(Index),FieldVarId(index)%nc_id(q),Field_val3d, &
835	start=(/ displ+ncell,1,FieldIndex(Index) /), count=(/n,size(field(ind)%rval4d,2),1 /))
836	! ENDDO
837	DEALLOCATE(field_val3d)
838	ENDDO
839	ENDIF
840
841	! PRINT *,NF90_STRERROR(status)
842	ncell=ncell+n
843
844	ENDDO
845
846	ENDIF
847	status=NF90_SYNC(FieldId(Index))
848
849	END SUBROUTINE Writefield_mpi
850
851
852
853	! SUBROUTINE Create_header(name,field,nind)
854	! USE netcdf
855	! USE field_mod
856	! USE domain_mod
857	! USE spherical_geom_mod
858	! USE dimensions
859	! USE geometry
860	! IMPLICIT NONE
861	! CHARACTER(LEN=*) :: name
862	! TYPE(t_field),POINTER :: field(:)
863	! INTEGER,OPTIONAL,INTENT(IN) :: nind
864	! INTEGER :: ncell
865	! INTEGER :: nvert
866	! REAL(rstd),ALLOCATABLE :: lon(:),lat(:),bounds_lon(:,:),bounds_lat(:,:)
867	! TYPE(t_domain),POINTER :: d
868	! INTEGER :: nvertId,ncid,lonId,latId,bounds_lonId,bounds_latId,timeId,ncellId
869	! INTEGER :: dim3id,dim4id
870	! INTEGER :: status
871	! INTEGER :: ind,i,j,k,n,q
872	! INTEGER :: iie,jje,iin,jjn
873	! INTEGER :: ind_b,ind_e
874	! INTEGER :: halo_size
875	! LOGICAL :: single
876	! INTEGER :: nij
877	!
878	! NbField=NbField+1
879	! FieldName(NbField)=TRIM(ADJUSTL(name))
880	! FieldIndex(NbField)=1
881	!
882	! IF (PRESENT(nind)) THEN
883	! ind_b=nind
884	! ind_e=nind
885	! halo_size=1
886	! single=.TRUE.
887	! ELSE
888	! ind_b=1
889	! ind_e=ndomain
890	! halo_size=0
891	! single=.FALSE.
892	! ENDIF
893	!
894	! ncell=0
895	!
896	! IF (Field(ind_b)%field_type==field_T) THEN
897	! nvert=6
898	!
899	! DO ind=ind_b,ind_e
900	! d=>domain(ind)
901	!
902	! DO j=d%jj_begin-halo_size,d%jj_end+halo_size
903	! DO i=d%ii_begin-halo_size,d%ii_end+halo_size
904	! IF (single .OR. domain(ind)%own(i,j)) ncell=ncell+1
905	! ENDDO
906	! ENDDO
907
908	! END DO
909	!
910	! status = NF90_CREATE(TRIM(ADJUSTL(name))//'.nc', NF90_CLOBBER, ncid)
911	! FieldId(NbField)=ncid
912	! status = NF90_DEF_DIM(ncid,'cell',ncell,ncellId)
913	! status = NF90_DEF_DIM(ncid,'nvert',nvert,nvertid)
914
915	! IF (Field(ind_b)%ndim==2) THEN
916	! FieldVarId(NbField)%size=1
917	! ALLOCATE(FieldVarId(NbField)%nc_id(1))
918	! ELSE IF (Field(ind_b)%ndim==3) THEN
919	! FieldVarId(NbField)%size=1
920	! ALLOCATE(FieldVarId(NbField)%nc_id(1))
921	! status = NF90_DEF_DIM(ncid,'Z',size(field(ind_b)%rval3d,2),dim3id)
922	! ELSE IF (Field(1)%ndim==4) THEN
923	! FieldVarId(NbField)%size=size(field(ind_b)%rval4d,3)
924	! ALLOCATE(FieldVarId(NbField)%nc_id(FieldVarId(NbField)%size))
925	! status = NF90_DEF_DIM(ncid,'Z',size(field(ind_b)%rval4d,2),dim3id)
926	! status = NF90_DEF_DIM(ncid,'Q',size(field(ind_b)%rval4d,3),dim4id)
927	! ENDIF
928	!
929	! status = NF90_DEF_DIM(ncid,'time_counter',NF90_UNLIMITED,timeId)
930	!
931	! status = NF90_DEF_VAR(ncid,'lon',NF90_DOUBLE,(/ ncellId /),lonId)
932	! status = NF90_PUT_ATT(ncid,lonId,"long_name","longitude")
933	! status = NF90_PUT_ATT(ncid,lonId,"units","degrees_east")
934	! status = NF90_PUT_ATT(ncid,lonId,"bounds","bounds_lon")
935	! status = NF90_DEF_VAR(ncid,'lat',NF90_DOUBLE,(/ ncellId /),latId)
936	! status = NF90_PUT_ATT(ncid,latId,"long_name","latitude")
937	! status = NF90_PUT_ATT(ncid,latId,"units","degrees_north")
938	! status = NF90_PUT_ATT(ncid,latId,"bounds","bounds_lat")
939	! status = NF90_DEF_VAR(ncid,'bounds_lon',NF90_DOUBLE,(/ nvertid,ncellId /),bounds_lonId)
940	! status = NF90_DEF_VAR(ncid,'bounds_lat',NF90_DOUBLE,(/ nvertid,ncellId /),bounds_latId)
941
942	! IF (Field(ind_b)%ndim==2) THEN
943	! status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name)),NF90_DOUBLE,(/ ncellId,timeId /),FieldVarId(NbField)%nc_id(1))
944	! status = NF90_PUT_ATT(ncid,FieldVarId(NbField)%nc_id(1),"coordinates","lon lat")
945	! ELSE IF (Field(ind_b)%ndim==3) THEN
946	! status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name)),NF90_DOUBLE,(/ ncellId,dim3id,timeId /),FieldVarId(NbField)%nc_id(1))
947	! status = NF90_PUT_ATT(ncid,FieldVarId(NbField)%nc_id(1),"coordinates","lon lat")
948	! ELSE IF (Field(ind_b)%ndim==4) THEN
949	! DO i=1,FieldVarId(NbField)%size
950	! status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name))//int2str(i),NF90_DOUBLE,(/ ncellId,dim3id,timeId /), &
951	! FieldVarId(NbField)%nc_id(i))
952	! status = NF90_PUT_ATT(ncid,FieldVarId(NbField)%nc_id(i),"coordinates","lon lat")
953	! ENDDO
954	! ENDIF
955	!
956	!
957	! status = NF90_ENDDEF(ncid)
958
959	! ncell=1
960	! DO ind=ind_b,ind_e
961	! d=>domain(ind)
962	!
963	! n=0
964	! DO j=d%jj_begin-halo_size,d%jj_end+halo_size
965	! DO i=d%ii_begin-halo_size,d%ii_end+halo_size
966	! IF (single .OR. d%own(i,j)) n=n+1
967	! ENDDO
968	! ENDDO
969	!
970	! ALLOCATE(lon(n),lat(n),bounds_lon(0:nvert-1,n),bounds_lat(0:nvert-1,n))
971	!
972	! n=0
973	! DO j=d%jj_begin-halo_size,d%jj_end+halo_size
974	! DO i=d%ii_begin-halo_size,d%ii_end+halo_size
975	! IF (d%own(i,j) .OR. single) THEN
976	! n=n+1
977	! CALL xyz2lonlat(d%xyz(:,i,j),lon(n),lat(n))
978	! lon(n)=lon(n)*180/Pi
979	! lat(n)=lat(n)*180/Pi
980	! DO k=0,5
981	! CALL xyz2lonlat(d%vertex(:,k,i,j),bounds_lon(k,n), bounds_lat(k,n))
982	! bounds_lat(k,n)=bounds_lat(k,n)*180/Pi
983	! bounds_lon(k,n)=bounds_lon(k,n)*180/Pi
984	! ENDDO
985	! ENDIF
986	! ENDDO
987	! ENDDO
988	! status=NF90_PUT_VAR(ncid,lonid,REAL(lon,r8),start=(/ ncell /),count=(/ n /))
989	! status=NF90_PUT_VAR(ncid,latid,REAL(lat,r8),start=(/ ncell /),count=(/ n /))
990	! status=NF90_PUT_VAR(ncid,bounds_lonId,REAL(bounds_lon,r8),start=(/ 1,ncell /),count=(/ nvert,n /))
991	! status=NF90_PUT_VAR(ncid,bounds_latId,REAL(bounds_lat,r8),start=(/ 1,ncell /),count=(/ nvert,n /))
992	!
993	! ncell=ncell+n
994	! DEALLOCATE(lon,lat,bounds_lon,bounds_lat)
995	! END DO
996
997	! ELSE IF (Field(ind_b)%field_type==field_Z) THEN
998	! nvert=3
999	! DO ind=ind_b,ind_e
1000	! d=>domain(ind)
1001	!
1002	! DO j=d%jj_begin+1,d%jj_end
1003	! DO i=d%ii_begin,d%ii_end-1
1004	! ncell=ncell+1
1005	! ENDDO
1006	! ENDDO
1007
1008	! DO j=d%jj_begin,d%jj_end-1
1009	! DO i=d%ii_begin+1,d%ii_end
1010	! ncell=ncell+1
1011	! ENDDO
1012	! ENDDO
1013
1014	! END DO
1015	!
1016	! status = NF90_CREATE(TRIM(ADJUSTL(name))//'.nc', NF90_CLOBBER, ncid)
1017	! FieldId(NbField)=ncid
1018	! status = NF90_DEF_DIM(ncid,'cell',ncell,ncellId)
1019	! status = NF90_DEF_DIM(ncid,'nvert',nvert,nvertid)
1020
1021	! IF (Field(ind_b)%ndim==2) THEN
1022	! FieldVarId(NbField)%size=1
1023	! ALLOCATE(FieldVarId(NbField)%nc_id(1))
1024	! ELSE IF (Field(ind_b)%ndim==3) THEN
1025	! FieldVarId(NbField)%size=1
1026	! ALLOCATE(FieldVarId(NbField)%nc_id(1))
1027	! status = NF90_DEF_DIM(ncid,'Z',size(field(ind_b)%rval3d,2),dim3id)
1028	! ELSE IF (Field(1)%ndim==4) THEN
1029	! FieldVarId(NbField)%size=size(field(ind_b)%rval4d,3)
1030	! ALLOCATE(FieldVarId(NbField)%nc_id(FieldVarId(NbField)%size))
1031	! status = NF90_DEF_DIM(ncid,'Z',size(field(ind_b)%rval4d,2),dim3id)
1032	! ENDIF
1033
1034
1035	!
1036	! status = NF90_DEF_DIM(ncid,'time_counter',NF90_UNLIMITED,timeId)
1037	!
1038	! status = NF90_DEF_VAR(ncid,'lon',NF90_DOUBLE,(/ ncellId /),lonId)
1039	! status = NF90_PUT_ATT(ncid,lonId,"long_name","longitude")
1040	! status = NF90_PUT_ATT(ncid,lonId,"units","degrees_east")
1041	! status = NF90_PUT_ATT(ncid,lonId,"bounds","bounds_lon")
1042	! status = NF90_DEF_VAR(ncid,'lat',NF90_DOUBLE,(/ ncellId /),latId)
1043	! status = NF90_PUT_ATT(ncid,latId,"long_name","latitude")
1044	! status = NF90_PUT_ATT(ncid,latId,"units","degrees_north")
1045	! status = NF90_PUT_ATT(ncid,latId,"bounds","bounds_lat")
1046	! status = NF90_DEF_VAR(ncid,'bounds_lon',NF90_DOUBLE,(/ nvertid,ncellId /),bounds_lonId)
1047	! status = NF90_DEF_VAR(ncid,'bounds_lat',NF90_DOUBLE,(/ nvertid,ncellId /),bounds_latId)
1048
1049
1050	! IF (Field(ind_b)%ndim==2) THEN
1051	! status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name)),NF90_DOUBLE,(/ ncellId,timeId /),FieldVarId(NbField)%nc_id(1))
1052	! status = NF90_PUT_ATT(ncid,FieldVarId(NbField)%nc_id(1),"coordinates","lon lat")
1053	! ELSE IF (Field(ind_b)%ndim==3) THEN
1054	! status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name)),NF90_DOUBLE,(/ ncellId,dim3id,timeId /),FieldVarId(NbField)%nc_id(1))
1055	! status = NF90_PUT_ATT(ncid,FieldVarId(NbField)%nc_id(1),"coordinates","lon lat")
1056	! ELSE IF (Field(ind_b)%ndim==4) THEN
1057	! DO q=1,FieldVarId(NbField)%size
1058	! status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name)//int2str(q)),NF90_DOUBLE, &
1059	! (/ ncellId,dim3id,timeId /),FieldVarId(NbField)%nc_id(q))
1060	! status = NF90_PUT_ATT(ncid,FieldVarId(NbField)%nc_id(q),"coordinates","lon lat")
1061	! ENDDO
1062	! ENDIF
1063	!
1064	! status = NF90_ENDDEF(ncid)
1065
1066	! ncell=1
1067	! DO ind=ind_b,ind_e
1068	! d=>domain(ind)
1069	! CALL swap_geometry(ind)
1070	! CALL swap_dimensions(ind)
1071	!
1072	! n=0
1073	! DO j=jj_begin+1,jj_end
1074	! DO i=ii_begin,ii_end-1
1075	! n=n+1
1076	! ENDDO
1077	! ENDDO
1078
1079	! DO j=jj_begin,jj_end-1
1080	! DO i=ii_begin+1,ii_end
1081	! n=n+1
1082	! ENDDO
1083	! ENDDO
1084
1085	! ALLOCATE(lon(n),lat(n),bounds_lon(0:nvert-1,n),bounds_lat(0:nvert-1,n))
1086	!
1087	! n=0
1088	!
1089	! DO j=jj_begin+1,jj_end
1090	! DO i=ii_begin,ii_end-1
1091	! nij=(j-1)*iim+i
1092	! n=n+1
1093	! CALL xyz2lonlat(xyz_v(nij+z_down,:)/radius,lon(n),lat(n))
1094	! lon(n)=lon(n)*180/Pi
1095	!! IF (lon(n)<0) lon(n)=lon(n)+360
1096	! lat(n)=lat(n)*180/Pi
1097	! CALL xyz2lonlat(xyz_i(nij,:)/radius,bounds_lon(0,n), bounds_lat(0,n))
1098	! CALL xyz2lonlat(xyz_i(nij+t_ldown,:)/radius,bounds_lon(1,n), bounds_lat(1,n))
1099	! CALL xyz2lonlat(xyz_i(nij+t_rdown,:)/radius,bounds_lon(2,n), bounds_lat(2,n))
1100
1101	! DO k=0,2
1102	! bounds_lat(k,n)=bounds_lat(k,n)*180/Pi
1103	! bounds_lon(k,n)=bounds_lon(k,n)*180/Pi
1104	! IF (bounds_lon(k,n)<0) bounds_lon(k,n)=bounds_lon(k,n)+360
1105	! ENDDO
1106	! ENDDO
1107	! ENDDO
1108
1109	! DO j=jj_begin,jj_end-1
1110	! DO i=ii_begin+1,ii_end
1111	! nij=(j-1)*iim+i
1112	! n=n+1
1113	! CALL xyz2lonlat(xyz_v(nij+z_up,:)/radius,lon(n),lat(n))
1114	! lon(n)=lon(n)*180/Pi
1115	! IF (lon(n)<0) lon(n)=lon(n)+360
1116	! lat(n)=lat(n)*180/Pi
1117	! CALL xyz2lonlat(xyz_i(nij,:)/radius,bounds_lon(0,n), bounds_lat(0,n))
1118	! CALL xyz2lonlat(xyz_i(nij+t_rup,:)/radius,bounds_lon(1,n), bounds_lat(1,n))
1119	! CALL xyz2lonlat(xyz_i(nij+t_lup,:)/radius,bounds_lon(2,n), bounds_lat(2,n))
1120
1121	! DO k=0,2
1122	! bounds_lat(k,n)=bounds_lat(k,n)*180/Pi
1123	! bounds_lon(k,n)=bounds_lon(k,n)*180/Pi
1124	! IF (bounds_lon(k,n)<0) bounds_lon(k,n)=bounds_lon(k,n)+360
1125	! ENDDO
1126	! ENDDO
1127	! ENDDO
1128	!
1129	!
1130	! status=NF90_PUT_VAR(ncid,lonid,REAL(lon,r8),start=(/ ncell /),count=(/ n /))
1131	! status=NF90_PUT_VAR(ncid,latid,REAL(lat,r8),start=(/ ncell /),count=(/ n /))
1132	! status=NF90_PUT_VAR(ncid,bounds_lonId,REAL(bounds_lon,r8),start=(/ 1,ncell /),count=(/ nvert,n /))
1133	! status=NF90_PUT_VAR(ncid,bounds_latId,REAL(bounds_lat,r8),start=(/ 1,ncell /),count=(/ nvert,n /))
1134	!
1135	! ncell=ncell+n
1136	! DEALLOCATE(lon,lat,bounds_lon,bounds_lat)
1137	! END DO
1138	! ENDIF
1139
1140
1141	!
1142	! status = NF90_CLOSE(ncid)
1143
1144	! END SUBROUTINE Create_Header
1145
1146
1147
1148
1149	SUBROUTINE Create_header_gen(name_in,field,domain_type,ind_b_in,ind_e_in,once)
1150	USE netcdf_mod
1151	USE field_mod
1152	USE domain_mod
1153	USE metric
1154	USE spherical_geom_mod
1155	! USE dimensions
1156	! USE geometry
1157	IMPLICIT NONE
1158	CHARACTER(LEN=*),INTENT(IN) :: name_in
1159	TYPE(t_field),POINTER :: field(:)
1160	TYPE(t_domain),INTENT(IN),TARGET :: domain_type(:)
1161	INTEGER,OPTIONAL,INTENT(IN) :: ind_b_in
1162	INTEGER,OPTIONAL,INTENT(IN) :: ind_e_in
1163	LOGICAL,OPTIONAL,INTENT(IN) :: once
1164
1165	INTEGER :: ncell
1166	INTEGER :: nvert
1167	REAL(rstd),ALLOCATABLE :: lon(:),lat(:),bounds_lon(:,:),bounds_lat(:,:)
1168	TYPE(t_domain),POINTER :: d
1169	INTEGER :: nvertId,ncid,lonId,latId,bounds_lonId,bounds_latId,timeId,ncellId
1170	INTEGER :: dim3id,dim4id
1171	INTEGER :: status
1172	INTEGER :: ind,i,j,k,n,q
1173	INTEGER :: iie,jje,iin,jjn
1174	INTEGER :: ind_b,ind_e
1175	INTEGER :: halo_size
1176	LOGICAL :: single
1177	INTEGER :: nij
1178	CHARACTER(LEN=255) :: name
1179
1180
1181	name=TRIM(ADJUSTL(name_in))
1182
1183	IF (PRESENT(ind_b_in) .AND. .NOT. PRESENT(ind_e_in)) THEN
1184	name=TRIM(name)//"_"//TRIM(int2str(ind_b))
1185	ind_b=ind_b_in
1186	ind_e=ind_b_in
1187	halo_size=1
1188	single=.TRUE.
1189	ELSE IF (.NOT. PRESENT(ind_b_in) .AND. PRESENT(ind_e_in)) THEN
1190	name=TRIM(name)//"_"//TRIM(int2str(ind_e))
1191	ind_b=ind_e_in
1192	ind_e=ind_e_in
1193	halo_size=1
1194	single=.TRUE.
1195	ELSE IF ( PRESENT(ind_b_in) .AND. PRESENT(ind_e_in)) THEN
1196	ind_b=ind_b_in
1197	ind_e=ind_e_in
1198	halo_size=0
1199	single=.FALSE.
1200	ELSE
1201	halo_size=0
1202	ind_b=1
1203	ind_e=ndomain
1204	single=.FALSE.
1205	ENDIF
1206
1207	NbField=NbField+1
1208	FieldName(NbField)=TRIM(ADJUSTL(name))
1209	FieldIndex(NbField)=1
1210
1211	ncell=0
1212
1213	IF (Field(ind_b)%field_type==field_T) THEN
1214	nvert=6
1215
1216	DO ind=ind_b,ind_e
1217	d=>domain_type(ind)
1218
1219	DO j=d%jj_begin-halo_size,d%jj_end+halo_size
1220	DO i=d%ii_begin-halo_size,d%ii_end+halo_size
1221	IF (single .OR. d%assign_domain(i,j)==ind) ncell=ncell+1
1222	ENDDO
1223	ENDDO
1224
1225	END DO
1226
1227	status = NF90_CREATE(TRIM(ADJUSTL(name))//'.nc', NF90_CLOBBER, ncid)
1228	FieldId(NbField)=ncid
1229	status = NF90_DEF_DIM(ncid,'cell_i',ncell,ncellId)
1230	status = NF90_DEF_DIM(ncid,'nvert_i',nvert,nvertid)
1231
1232	IF (Field(ind_b)%ndim==2) THEN
1233	FieldVarId(NbField)%size=1
1234	ALLOCATE(FieldVarId(NbField)%nc_id(1))
1235	ELSE IF (Field(ind_b)%ndim==3) THEN
1236	FieldVarId(NbField)%size=1
1237	ALLOCATE(FieldVarId(NbField)%nc_id(1))
1238	status = NF90_DEF_DIM(ncid,'lev',size(field(ind_b)%rval3d,2),dim3id)
1239	ELSE IF (Field(1)%ndim==4) THEN
1240	FieldVarId(NbField)%size=size(field(ind_b)%rval4d,3)
1241	ALLOCATE(FieldVarId(NbField)%nc_id(FieldVarId(NbField)%size))
1242	status = NF90_DEF_DIM(ncid,'lev',size(field(ind_b)%rval4d,2),dim3id)
1243	! status = NF90_DEF_DIM(ncid,'Q',size(field(ind_b)%rval4d,3),dim4id)
1244	ENDIF
1245
1246	status = NF90_DEF_DIM(ncid,'time_counter',NF90_UNLIMITED,timeId)
1247
1248	status = NF90_DEF_VAR(ncid,'lon_i',NF90_DOUBLE,(/ ncellId /),lonId)
1249	status = NF90_PUT_ATT(ncid,lonId,"long_name","longitude")
1250	status = NF90_PUT_ATT(ncid,lonId,"units","degrees_east")
1251	status = NF90_PUT_ATT(ncid,lonId,"bounds","bounds_lon_i")
1252	status = NF90_DEF_VAR(ncid,'lat_i',NF90_DOUBLE,(/ ncellId /),latId)
1253	status = NF90_PUT_ATT(ncid,latId,"long_name","latitude")
1254	status = NF90_PUT_ATT(ncid,latId,"units","degrees_north")
1255	status = NF90_PUT_ATT(ncid,latId,"bounds","bounds_lat_i")
1256	status = NF90_DEF_VAR(ncid,'bounds_lon_i',NF90_DOUBLE,(/ nvertid,ncellId /),bounds_lonId)
1257	status = NF90_DEF_VAR(ncid,'bounds_lat_i',NF90_DOUBLE,(/ nvertid,ncellId /),bounds_latId)
1258
1259	IF (Field(ind_b)%ndim==2) THEN
1260	IF (PRESENT(once) .AND. once) THEN
1261	status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name)),ncprec,(/ ncellId /),FieldVarId(NbField)%nc_id(1))
1262	ELSE
1263	status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name)),ncprec,(/ ncellId,timeId /),FieldVarId(NbField)%nc_id(1))
1264	ENDIF
1265	status = NF90_PUT_ATT(ncid,FieldVarId(NbField)%nc_id(1),"coordinates","lon_i lat_i")
1266	ELSE IF (Field(ind_b)%ndim==3) THEN
1267	IF (PRESENT(once) .AND. once) THEN
1268	status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name)),ncprec,(/ ncellId,dim3id /),FieldVarId(NbField)%nc_id(1))
1269	ELSE
1270	status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name)),ncprec,(/ ncellId,dim3id,timeId /),FieldVarId(NbField)%nc_id(1))
1271	ENDIF
1272	status = NF90_PUT_ATT(ncid,FieldVarId(NbField)%nc_id(1),"coordinates","lon_i lat_i")
1273	ELSE IF (Field(ind_b)%ndim==4) THEN
1274	DO i=1,FieldVarId(NbField)%size
1275	IF (PRESENT(once) .AND. once) THEN
1276	status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name))//int2str(i),ncprec,(/ ncellId,dim3id /), &
1277	FieldVarId(NbField)%nc_id(i))
1278	ELSE
1279	status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name))//int2str(i),ncprec,(/ ncellId,dim3id,timeId /), &
1280	FieldVarId(NbField)%nc_id(i))
1281	ENDIF
1282	status = NF90_PUT_ATT(ncid,FieldVarId(NbField)%nc_id(i),"coordinates","lon_i lat_i")
1283	ENDDO
1284	ENDIF
1285
1286
1287	status = NF90_ENDDEF(ncid)
1288
1289	! ncell=1
1290	! DO ind=ind_b,ind_e
1291	! d=>domain_type(ind)
1292
1293	! n=0
1294	! DO j=d%jj_begin-halo_size,d%jj_end+halo_size
1295	! DO i=d%ii_begin-halo_size,d%ii_end+halo_size
1296	! IF (single .OR. d%assign_domain(i,j)==ind) n=n+1
1297	! ENDDO
1298	! ENDDO
1299
1300	ALLOCATE(lon(ncell),lat(ncell),bounds_lon(0:nvert-1,ncell),bounds_lat(0:nvert-1,ncell))
1301
1302	n=0
1303	DO ind=ind_b,ind_e
1304	d=>domain_type(ind)
1305	DO j=d%jj_begin-halo_size,d%jj_end+halo_size
1306	DO i=d%ii_begin-halo_size,d%ii_end+halo_size
1307	IF (d%assign_domain(i,j)==ind .OR. single) THEN
1308	n=n+1
1309	CALL xyz2lonlat(d%xyz(:,i,j),lon(n),lat(n))
1310	lon(n)=lon(n)*180/Pi
1311	lat(n)=lat(n)*180/Pi
1312	DO k=0,5
1313	CALL xyz2lonlat(d%vertex(:,k,i,j),bounds_lon(k,n), bounds_lat(k,n))
1314	bounds_lat(k,n)=bounds_lat(k,n)*180/Pi
1315	bounds_lon(k,n)=bounds_lon(k,n)*180/Pi
1316	ENDDO
1317	ENDIF
1318	ENDDO
1319	ENDDO
1320	ENDDO
1321	status=NF90_PUT_VAR(ncid,lonid,REAL(lon,r8),start=(/ 1 /),count=(/ ncell /))
1322	status=NF90_PUT_VAR(ncid,latid,REAL(lat,r8),start=(/ 1 /),count=(/ ncell /))
1323	status=NF90_PUT_VAR(ncid,bounds_lonId,REAL(bounds_lon,r8),start=(/ 1,1 /),count=(/ nvert,ncell /))
1324	status=NF90_PUT_VAR(ncid,bounds_latId,REAL(bounds_lat,r8),start=(/ 1,1 /),count=(/ nvert,ncell /))
1325
1326	! ncell=ncell+n
1327	DEALLOCATE(lon,lat,bounds_lon,bounds_lat)
1328	! END DO
1329
1330	ELSE IF (Field(ind_b)%field_type==field_Z) THEN
1331	nvert=3
1332	DO ind=ind_b,ind_e
1333	d=>domain_type(ind)
1334
1335	DO j=d%jj_begin+1,d%jj_end
1336	DO i=d%ii_begin,d%ii_end-1
1337	ncell=ncell+1
1338	ENDDO
1339	ENDDO
1340
1341	DO j=d%jj_begin,d%jj_end-1
1342	DO i=d%ii_begin+1,d%ii_end
1343	ncell=ncell+1
1344	ENDDO
1345	ENDDO
1346
1347	END DO
1348
1349	status = NF90_CREATE(TRIM(ADJUSTL(name))//'.nc', NF90_CLOBBER, ncid)
1350	FieldId(NbField)=ncid
1351	status = NF90_DEF_DIM(ncid,'cell_v',ncell,ncellId)
1352	status = NF90_DEF_DIM(ncid,'nvert_v',nvert,nvertid)
1353
1354	IF (Field(ind_b)%ndim==2) THEN
1355	FieldVarId(NbField)%size=1
1356	ALLOCATE(FieldVarId(NbField)%nc_id(1))
1357	ELSE IF (Field(ind_b)%ndim==3) THEN
1358	FieldVarId(NbField)%size=1
1359	ALLOCATE(FieldVarId(NbField)%nc_id(1))
1360	status = NF90_DEF_DIM(ncid,'lev',size(field(ind_b)%rval3d,2),dim3id)
1361	ELSE IF (Field(1)%ndim==4) THEN
1362	FieldVarId(NbField)%size=size(field(ind_b)%rval4d,3)
1363	ALLOCATE(FieldVarId(NbField)%nc_id(FieldVarId(NbField)%size))
1364	status = NF90_DEF_DIM(ncid,'lev',size(field(ind_b)%rval4d,2),dim3id)
1365	ENDIF
1366
1367
1368
1369	status = NF90_DEF_DIM(ncid,'time_counter',NF90_UNLIMITED,timeId)
1370
1371	status = NF90_DEF_VAR(ncid,'lon_v',NF90_DOUBLE,(/ ncellId /),lonId)
1372	status = NF90_PUT_ATT(ncid,lonId,"long_name","longitude")
1373	status = NF90_PUT_ATT(ncid,lonId,"units","degrees_east")
1374	status = NF90_PUT_ATT(ncid,lonId,"bounds","bounds_lon_v")
1375	status = NF90_DEF_VAR(ncid,'lat_v',NF90_DOUBLE,(/ ncellId /),latId)
1376	status = NF90_PUT_ATT(ncid,latId,"long_name","latitude")
1377	status = NF90_PUT_ATT(ncid,latId,"units","degrees_north")
1378	status = NF90_PUT_ATT(ncid,latId,"bounds","bounds_lat_v")
1379	status = NF90_DEF_VAR(ncid,'bounds_lon_v',NF90_DOUBLE,(/ nvertid,ncellId /),bounds_lonId)
1380	status = NF90_DEF_VAR(ncid,'bounds_lat_v',NF90_DOUBLE,(/ nvertid,ncellId /),bounds_latId)
1381
1382
1383	IF (Field(ind_b)%ndim==2) THEN
1384	IF (PRESENT(once) .AND. once) THEN
1385	status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name)),ncprec,(/ ncellId /),FieldVarId(NbField)%nc_id(1))
1386	ELSE
1387	status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name)),ncprec,(/ ncellId,timeId /),FieldVarId(NbField)%nc_id(1))
1388	ENDIF
1389	status = NF90_PUT_ATT(ncid,FieldVarId(NbField)%nc_id(1),"coordinates","lon_v lat_v")
1390	ELSE IF (Field(ind_b)%ndim==3) THEN
1391	IF (PRESENT(once) .AND. once) THEN
1392	status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name)),ncprec,(/ ncellId,dim3id /),FieldVarId(NbField)%nc_id(1))
1393	ELSE
1394	status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name)),ncprec,(/ ncellId,dim3id,timeId /),FieldVarId(NbField)%nc_id(1))
1395	ENDIF
1396	status = NF90_PUT_ATT(ncid,FieldVarId(NbField)%nc_id(1),"coordinates","lon_v lat_v")
1397	ELSE IF (Field(ind_b)%ndim==4) THEN
1398	DO q=1,FieldVarId(NbField)%size
1399	IF (PRESENT(once) .AND. once) THEN
1400	status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name)//int2str(q)),ncprec, &
1401	(/ ncellId,dim3id /),FieldVarId(NbField)%nc_id(q))
1402	ELSE
1403	status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name)//int2str(q)),ncprec, &
1404	(/ ncellId,dim3id,timeId /),FieldVarId(NbField)%nc_id(q))
1405	ENDIF
1406	status = NF90_PUT_ATT(ncid,FieldVarId(NbField)%nc_id(q),"coordinates","lon_v lat_v")
1407	ENDDO
1408	ENDIF
1409
1410	status = NF90_ENDDEF(ncid)
1411
1412	! ncell=1
1413	! DO ind=ind_b,ind_e
1414	! d=>domain_type(ind)
1415	!
1416	! n=0
1417	! DO j=d%jj_begin+1,d%jj_end
1418	! DO i=d%ii_begin,d%ii_end-1
1419	! n=n+1
1420	! ENDDO
1421	! ENDDO
1422	!
1423	! DO j=d%jj_begin,d%jj_end-1
1424	! DO i=d%ii_begin+1,d%ii_end
1425	! n=n+1
1426	! ENDDO
1427	! ENDDO
1428
1429	! ALLOCATE(lon(n),lat(n),bounds_lon(0:nvert-1,n),bounds_lat(0:nvert-1,n))
1430	ALLOCATE(lon(ncell),lat(ncell),bounds_lon(0:nvert-1,ncell),bounds_lat(0:nvert-1,ncell))
1431
1432	n=0
1433
1434	DO ind=ind_b,ind_e
1435	d=>domain_type(ind)
1436	DO j=d%jj_begin+1,d%jj_end
1437	DO i=d%ii_begin,d%ii_end-1
1438	nij=(j-1)*d%iim+i
1439	n=n+1
1440	CALL xyz2lonlat(d%vertex(:,vdown,i,j),lon(n),lat(n))
1441	lon(n)=lon(n)*180/Pi
1442	! IF (lon(n)<0) lon(n)=lon(n)+360
1443	lat(n)=lat(n)*180/Pi
1444	! CALL xyz2lonlat(xyz_i(nij,:)/radius,bounds_lon(0,n), bounds_lat(0,n))
1445	! CALL xyz2lonlat(xyz_i(nij+t_ldown,:)/radius,bounds_lon(1,n), bounds_lat(1,n))
1446	! CALL xyz2lonlat(xyz_i(nij+t_rdown,:)/radius,bounds_lon(2,n), bounds_lat(2,n))
1447
1448	CALL xyz2lonlat(d%xyz(:,i,j),bounds_lon(0,n), bounds_lat(0,n))
1449	CALL xyz2lonlat(d%xyz(:,i,j-1),bounds_lon(1,n), bounds_lat(1,n))
1450	CALL xyz2lonlat(d%xyz(:,i+1,j-1),bounds_lon(2,n), bounds_lat(2,n))
1451
1452	DO k=0,2
1453	bounds_lat(k,n)=bounds_lat(k,n)*180/Pi
1454	bounds_lon(k,n)=bounds_lon(k,n)*180/Pi
1455	! IF (bounds_lon(k,n)<0) bounds_lon(k,n)=bounds_lon(k,n)+360
1456	ENDDO
1457	ENDDO
1458	ENDDO
1459
1460	DO j=d%jj_begin,d%jj_end-1
1461	DO i=d%ii_begin+1,d%ii_end
1462	nij=(j-1)*d%iim+i
1463	n=n+1
1464	! CALL xyz2lonlat(xyz_v(nij+z_up,:)/radius,lon(n),lat(n))
1465	CALL xyz2lonlat(d%vertex(:,vup,i,j),lon(n),lat(n))
1466	lon(n)=lon(n)*180/Pi
1467	! IF (lon(n)<0) lon(n)=lon(n)+360
1468	lat(n)=lat(n)*180/Pi
1469	! CALL xyz2lonlat(xyz_i(nij,:)/radius,bounds_lon(0,n), bounds_lat(0,n))
1470	! CALL xyz2lonlat(xyz_i(nij+t_rup,:)/radius,bounds_lon(1,n), bounds_lat(1,n))
1471	! CALL xyz2lonlat(xyz_i(nij+t_lup,:)/radius,bounds_lon(2,n), bounds_lat(2,n))
1472	CALL xyz2lonlat(d%xyz(:,i,j),bounds_lon(0,n), bounds_lat(0,n))
1473	CALL xyz2lonlat(d%xyz(:,i,j+1),bounds_lon(1,n), bounds_lat(1,n))
1474	CALL xyz2lonlat(d%xyz(:,i-1,j+1),bounds_lon(2,n), bounds_lat(2,n))
1475
1476	DO k=0,2
1477	bounds_lat(k,n)=bounds_lat(k,n)*180/Pi
1478	bounds_lon(k,n)=bounds_lon(k,n)*180/Pi
1479	! IF (bounds_lon(k,n)<0) bounds_lon(k,n)=bounds_lon(k,n)+360
1480	ENDDO
1481	ENDDO
1482	ENDDO
1483	ENDDO
1484
1485	status=NF90_PUT_VAR(ncid,lonid,REAL(lon,r8),start=(/ 1 /),count=(/ ncell /))
1486	status=NF90_PUT_VAR(ncid,latid,REAL(lat,r8),start=(/ 1 /),count=(/ ncell /))
1487	status=NF90_PUT_VAR(ncid,bounds_lonId,REAL(bounds_lon,r8),start=(/ 1,1 /),count=(/ nvert,ncell /))
1488	status=NF90_PUT_VAR(ncid,bounds_latId,REAL(bounds_lat,r8),start=(/ 1,1 /),count=(/ nvert,ncell /))
1489
1490	! ncell=ncell+n
1491	DEALLOCATE(lon,lat,bounds_lon,bounds_lat)
1492	! END DO
1493	ENDIF
1494
1495
1496
1497	! status = NF90_CLOSE(ncid)
1498
1499	END SUBROUTINE Create_Header_gen
1500
1501	SUBROUTINE Create_header_mpi(name,field,nind)
1502	USE netcdf_mod
1503	USE field_mod
1504	USE domain_mod
1505	USE spherical_geom_mod
1506	USE dimensions
1507	USE geometry
1508	USE mpi_mod
1509	USE mpipara
1510	IMPLICIT NONE
1511	CHARACTER(LEN=*) :: name
1512	TYPE(t_field),POINTER :: field(:)
1513	INTEGER,OPTIONAL,INTENT(IN) :: nind
1514	INTEGER :: ncell
1515	INTEGER :: nvert
1516	REAL(rstd),ALLOCATABLE :: lon(:),lat(:),bounds_lon(:,:),bounds_lat(:,:)
1517	TYPE(t_domain),POINTER :: d
1518	INTEGER :: nvertId,ncid,lonId,latId,bounds_lonId,bounds_latId,timeId,ncellId
1519	INTEGER :: dim3id,dim4id
1520	INTEGER :: status
1521	INTEGER :: ind,i,j,k,n,q
1522	INTEGER :: iie,jje,iin,jjn
1523	INTEGER :: ind_b,ind_e
1524	INTEGER :: halo_size
1525	LOGICAL :: single
1526	INTEGER :: nij
1527	INTEGER :: ncell_glo(0:mpi_size-1)
1528	INTEGER :: displ, ncell_tot
1529
1530
1531	NbField=NbField+1
1532	FieldName(NbField)=TRIM(ADJUSTL(name))
1533	FieldIndex(NbField)=1
1534
1535	IF (PRESENT(nind)) THEN
1536	ind_b=nind
1537	ind_e=nind
1538	halo_size=1
1539	single=.TRUE.
1540	ELSE
1541	ind_b=1
1542	ind_e=ndomain
1543	halo_size=0
1544	single=.FALSE.
1545	ENDIF
1546
1547	ncell=0
1548
1549	IF (Field(ind_b)%field_type==field_T) THEN
1550	nvert=6
1551
1552	DO ind=ind_b,ind_e
1553	d=>domain(ind)
1554
1555	DO j=d%jj_begin-halo_size,d%jj_end+halo_size
1556	DO i=d%ii_begin-halo_size,d%ii_end+halo_size
1557	IF (single .OR. domain(ind)%own(i,j)) ncell=ncell+1
1558	ENDDO
1559	ENDDO
1560
1561	END DO
1562
1563	CALL MPI_ALLGATHER(ncell,1,MPI_INTEGER,ncell_glo,1,MPI_INTEGER,comm_icosa,ierr)
1564
1565	displ=0
1566	DO i=1,mpi_rank
1567	displ=displ+ncell_glo(i-1)
1568	ENDDO
1569	FieldVarId(NbField)%displ=displ
1570	ncell_tot=sum(ncell_glo(:))
1571
1572	! status = NF90_CREATE_PAR(TRIM(ADJUSTL(name))//'.nc', IOR(NF90_NETCDF4, NF90_MPIIO), comm_icosa, MPI_INFO_NULL, ncid)
1573	FieldId(NbField)=ncid
1574	status = NF90_DEF_DIM(ncid,'cell',ncell_tot,ncellId)
1575	status = NF90_DEF_DIM(ncid,'nvert',nvert,nvertid)
1576
1577	IF (Field(ind_b)%ndim==2) THEN
1578	FieldVarId(NbField)%size=1
1579	ALLOCATE(FieldVarId(NbField)%nc_id(1))
1580	ELSE IF (Field(ind_b)%ndim==3) THEN
1581	FieldVarId(NbField)%size=1
1582	ALLOCATE(FieldVarId(NbField)%nc_id(1))
1583	status = NF90_DEF_DIM(ncid,'lev',size(field(ind_b)%rval3d,2),dim3id)
1584	ELSE IF (Field(1)%ndim==4) THEN
1585	FieldVarId(NbField)%size=size(field(ind_b)%rval4d,3)
1586	ALLOCATE(FieldVarId(NbField)%nc_id(FieldVarId(NbField)%size))
1587	status = NF90_DEF_DIM(ncid,'lev',size(field(ind_b)%rval4d,2),dim3id)
1588	! status = NF90_DEF_DIM(ncid,'Q',size(field(ind_b)%rval4d,3),dim4id)
1589	ENDIF
1590
1591	status = NF90_DEF_DIM(ncid,'time_counter',NF90_UNLIMITED,timeId)
1592
1593	status = NF90_DEF_VAR(ncid,'lon',NF90_DOUBLE,(/ ncellId /),lonId)
1594	status = NF90_PUT_ATT(ncid,lonId,"long_name","longitude")
1595	status = NF90_PUT_ATT(ncid,lonId,"units","degrees_east")
1596	status = NF90_PUT_ATT(ncid,lonId,"bounds","bounds_lon")
1597	status = NF90_DEF_VAR(ncid,'lat',NF90_DOUBLE,(/ ncellId /),latId)
1598	status = NF90_PUT_ATT(ncid,latId,"long_name","latitude")
1599	status = NF90_PUT_ATT(ncid,latId,"units","degrees_north")
1600	status = NF90_PUT_ATT(ncid,latId,"bounds","bounds_lat")
1601	status = NF90_DEF_VAR(ncid,'bounds_lon',NF90_DOUBLE,(/ nvertid,ncellId /),bounds_lonId)
1602	status = NF90_DEF_VAR(ncid,'bounds_lat',NF90_DOUBLE,(/ nvertid,ncellId /),bounds_latId)
1603
1604	IF (Field(ind_b)%ndim==2) THEN
1605	status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name)),ncprec,(/ ncellId,timeId /),FieldVarId(NbField)%nc_id(1))
1606	status = NF90_PUT_ATT(ncid,FieldVarId(NbField)%nc_id(1),"coordinates","lon lat")
1607	status = NF90_DEF_VAR_CHUNKING(ncid, FieldVarId(NbField)%nc_id(1), NF90_CHUNKED, (/ncell_tot,1/))
1608	ELSE IF (Field(ind_b)%ndim==3) THEN
1609	status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name)),ncprec,(/ ncellId,dim3id,timeId /),FieldVarId(NbField)%nc_id(1))
1610	status = NF90_PUT_ATT(ncid,FieldVarId(NbField)%nc_id(1),"coordinates","lon lat")
1611	status = NF90_DEF_VAR_CHUNKING(ncid, FieldVarId(NbField)%nc_id(1), NF90_CHUNKED, &
1612	(/ncell_tot,size(field(ind_b)%rval3d,2),1/))
1613	ELSE IF (Field(ind_b)%ndim==4) THEN
1614	DO i=1,FieldVarId(NbField)%size
1615	status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name))//int2str(i),ncprec,(/ ncellId,dim3id,timeId /), &
1616	FieldVarId(NbField)%nc_id(i))
1617	status = NF90_PUT_ATT(ncid,FieldVarId(NbField)%nc_id(i),"coordinates","lon lat")
1618	status = NF90_DEF_VAR_CHUNKING(ncid, FieldVarId(NbField)%nc_id(q), NF90_CHUNKED, &
1619	(/ncell_tot,size(field(ind_b)%rval4d,2),1/))
1620	ENDDO
1621	ENDIF
1622
1623
1624	status = NF90_ENDDEF(ncid)
1625
1626	ncell=1
1627	DO ind=ind_b,ind_e
1628	d=>domain(ind)
1629
1630	n=0
1631	DO j=d%jj_begin-halo_size,d%jj_end+halo_size
1632	DO i=d%ii_begin-halo_size,d%ii_end+halo_size
1633	IF (single .OR. d%own(i,j)) n=n+1
1634	ENDDO
1635	ENDDO
1636
1637	ALLOCATE(lon(n),lat(n),bounds_lon(0:nvert-1,n),bounds_lat(0:nvert-1,n))
1638
1639	n=0
1640	DO j=d%jj_begin-halo_size,d%jj_end+halo_size
1641	DO i=d%ii_begin-halo_size,d%ii_end+halo_size
1642	IF (d%own(i,j) .OR. single) THEN
1643	n=n+1
1644	CALL xyz2lonlat(d%xyz(:,i,j),lon(n),lat(n))
1645	lon(n)=lon(n)*180/Pi
1646	lat(n)=lat(n)*180/Pi
1647	DO k=0,5
1648	CALL xyz2lonlat(d%vertex(:,k,i,j),bounds_lon(k,n), bounds_lat(k,n))
1649	bounds_lat(k,n)=bounds_lat(k,n)*180/Pi
1650	bounds_lon(k,n)=bounds_lon(k,n)*180/Pi
1651	ENDDO
1652	ENDIF
1653	ENDDO
1654	ENDDO
1655	status=NF90_PUT_VAR(ncid,lonid,REAL(lon,r8),start=(/ displ+ncell /),count=(/ n /))
1656	status=NF90_PUT_VAR(ncid,latid,REAL(lat,r8),start=(/ displ+ncell /),count=(/ n /))
1657	status=NF90_PUT_VAR(ncid,bounds_lonId,REAL(bounds_lon,r8),start=(/ 1,displ+ncell /),count=(/ nvert,n /))
1658	status=NF90_PUT_VAR(ncid,bounds_latId,REAL(bounds_lat,r8),start=(/ 1,displ+ncell /),count=(/ nvert,n /))
1659
1660	ncell=ncell+n
1661	DEALLOCATE(lon,lat,bounds_lon,bounds_lat)
1662	END DO
1663
1664	ELSE IF (Field(ind_b)%field_type==field_Z) THEN
1665	nvert=3
1666	DO ind=ind_b,ind_e
1667	d=>domain(ind)
1668
1669	DO j=d%jj_begin+1,d%jj_end
1670	DO i=d%ii_begin,d%ii_end-1
1671	ncell=ncell+1
1672	ENDDO
1673	ENDDO
1674
1675	DO j=d%jj_begin,d%jj_end-1
1676	DO i=d%ii_begin+1,d%ii_end
1677	ncell=ncell+1
1678	ENDDO
1679	ENDDO
1680
1681	END DO
1682
1683	CALL MPI_ALLGATHER(ncell,1,MPI_INTEGER,ncell_glo,1,MPI_INTEGER,comm_icosa,ierr)
1684
1685	displ=0
1686	DO i=1,mpi_rank
1687	displ=displ+ncell_glo(i-1)
1688	ENDDO
1689	FieldVarId(NbField)%displ=displ
1690	ncell_tot=sum(ncell_glo(:))
1691
1692	! status = NF90_CREATE_PAR(TRIM(ADJUSTL(name))//'.nc',IOR(NF90_NETCDF4, NF90_MPIIO), comm_icosa, MPI_INFO_NULL, ncid)
1693	! status = NF90_CREATE(TRIM(ADJUSTL(name))//'.nc', NF90_CLOBBER, ncid)
1694	FieldId(NbField)=ncid
1695	status = NF90_DEF_DIM(ncid,'cell',ncell_tot,ncellId)
1696	status = NF90_DEF_DIM(ncid,'nvert',nvert,nvertid)
1697
1698	IF (Field(ind_b)%ndim==2) THEN
1699	FieldVarId(NbField)%size=1
1700	ALLOCATE(FieldVarId(NbField)%nc_id(1))
1701	ELSE IF (Field(ind_b)%ndim==3) THEN
1702	FieldVarId(NbField)%size=1
1703	ALLOCATE(FieldVarId(NbField)%nc_id(1))
1704	status = NF90_DEF_DIM(ncid,'lev',size(field(ind_b)%rval3d,2),dim3id)
1705	ELSE IF (Field(1)%ndim==4) THEN
1706	FieldVarId(NbField)%size=size(field(ind_b)%rval4d,3)
1707	ALLOCATE(FieldVarId(NbField)%nc_id(FieldVarId(NbField)%size))
1708	status = NF90_DEF_DIM(ncid,'lev',size(field(ind_b)%rval4d,2),dim3id)
1709	ENDIF
1710
1711
1712
1713	status = NF90_DEF_DIM(ncid,'time_counter',NF90_UNLIMITED,timeId)
1714
1715	status = NF90_DEF_VAR(ncid,'lon',NF90_DOUBLE,(/ ncellId /),lonId)
1716	status = NF90_PUT_ATT(ncid,lonId,"long_name","longitude")
1717	status = NF90_PUT_ATT(ncid,lonId,"units","degrees_east")
1718	status = NF90_PUT_ATT(ncid,lonId,"bounds","bounds_lon")
1719	status = NF90_DEF_VAR(ncid,'lat',NF90_DOUBLE,(/ ncellId /),latId)
1720	status = NF90_PUT_ATT(ncid,latId,"long_name","latitude")
1721	status = NF90_PUT_ATT(ncid,latId,"units","degrees_north")
1722	status = NF90_PUT_ATT(ncid,latId,"bounds","bounds_lat")
1723	status = NF90_DEF_VAR(ncid,'bounds_lon',NF90_DOUBLE,(/ nvertid,ncellId /),bounds_lonId)
1724	status = NF90_DEF_VAR(ncid,'bounds_lat',NF90_DOUBLE,(/ nvertid,ncellId /),bounds_latId)
1725
1726
1727	IF (Field(ind_b)%ndim==2) THEN
1728	status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name)),ncprec,(/ ncellId,timeId /),FieldVarId(NbField)%nc_id(1))
1729	status = NF90_PUT_ATT(ncid,FieldVarId(NbField)%nc_id(1),"coordinates","lon lat")
1730	status = NF90_DEF_VAR_CHUNKING(ncid, FieldVarId(NbField)%nc_id(1), NF90_CHUNKED, (/ncell_tot,1/))
1731	ELSE IF (Field(ind_b)%ndim==3) THEN
1732	status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name)),ncprec,(/ ncellId,dim3id,timeId /),FieldVarId(NbField)%nc_id(1))
1733	status = NF90_PUT_ATT(ncid,FieldVarId(NbField)%nc_id(1),"coordinates","lon lat")
1734	status = NF90_DEF_VAR_CHUNKING(ncid, FieldVarId(NbField)%nc_id(1), NF90_CHUNKED, &
1735	(/ncell_tot,size(field(ind_b)%rval3d,2),1/))
1736	ELSE IF (Field(ind_b)%ndim==4) THEN
1737	DO q=1,FieldVarId(NbField)%size
1738	status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(name)//int2str(q)),ncprec, &
1739	(/ ncellId,dim3id,timeId /),FieldVarId(NbField)%nc_id(q))
1740	status = NF90_PUT_ATT(ncid,FieldVarId(NbField)%nc_id(q),"coordinates","lon lat")
1741	status = NF90_DEF_VAR_CHUNKING(ncid, FieldVarId(NbField)%nc_id(q), NF90_CHUNKED, &
1742	(/ncell_tot,size(field(ind_b)%rval4d,2),1/))
1743	ENDDO
1744	ENDIF
1745
1746	status = NF90_ENDDEF(ncid)
1747
1748	ncell=1
1749	DO ind=ind_b,ind_e
1750	d=>domain(ind)
1751	CALL swap_geometry(ind)
1752	CALL swap_dimensions(ind)
1753
1754	n=0
1755	DO j=jj_begin+1,jj_end
1756	DO i=ii_begin,ii_end-1
1757	n=n+1
1758	ENDDO
1759	ENDDO
1760
1761	DO j=jj_begin,jj_end-1
1762	DO i=ii_begin+1,ii_end
1763	n=n+1
1764	ENDDO
1765	ENDDO
1766
1767	ALLOCATE(lon(n),lat(n),bounds_lon(0:nvert-1,n),bounds_lat(0:nvert-1,n))
1768
1769	n=0
1770
1771	DO j=jj_begin+1,jj_end
1772	DO i=ii_begin,ii_end-1
1773	nij=(j-1)*iim+i
1774	n=n+1
1775	CALL xyz2lonlat(xyz_v(nij+z_down,:)/radius,lon(n),lat(n))
1776	lon(n)=lon(n)*180/Pi
1777	! IF (lon(n)<0) lon(n)=lon(n)+360
1778	lat(n)=lat(n)*180/Pi
1779	CALL xyz2lonlat(xyz_i(nij,:)/radius,bounds_lon(0,n), bounds_lat(0,n))
1780	CALL xyz2lonlat(xyz_i(nij+t_ldown,:)/radius,bounds_lon(1,n), bounds_lat(1,n))
1781	CALL xyz2lonlat(xyz_i(nij+t_rdown,:)/radius,bounds_lon(2,n), bounds_lat(2,n))
1782
1783	DO k=0,2
1784	bounds_lat(k,n)=bounds_lat(k,n)*180/Pi
1785	bounds_lon(k,n)=bounds_lon(k,n)*180/Pi
1786	! IF (bounds_lon(k,n)<0) bounds_lon(k,n)=bounds_lon(k,n)+360
1787	ENDDO
1788	ENDDO
1789	ENDDO
1790
1791	DO j=jj_begin,jj_end-1
1792	DO i=ii_begin+1,ii_end
1793	nij=(j-1)*iim+i
1794	n=n+1
1795	CALL xyz2lonlat(xyz_v(nij+z_up,:)/radius,lon(n),lat(n))
1796	lon(n)=lon(n)*180/Pi
1797	! IF (lon(n)<0) lon(n)=lon(n)+360
1798	lat(n)=lat(n)*180/Pi
1799	CALL xyz2lonlat(xyz_i(nij,:)/radius,bounds_lon(0,n), bounds_lat(0,n))
1800	CALL xyz2lonlat(xyz_i(nij+t_rup,:)/radius,bounds_lon(1,n), bounds_lat(1,n))
1801	CALL xyz2lonlat(xyz_i(nij+t_lup,:)/radius,bounds_lon(2,n), bounds_lat(2,n))
1802
1803	DO k=0,2
1804	bounds_lat(k,n)=bounds_lat(k,n)*180/Pi
1805	bounds_lon(k,n)=bounds_lon(k,n)*180/Pi
1806	! IF (bounds_lon(k,n)<0) bounds_lon(k,n)=bounds_lon(k,n)+360
1807	ENDDO
1808	ENDDO
1809	ENDDO
1810
1811
1812	status=NF90_PUT_VAR(ncid,lonid,REAL(lon,r8),start=(/ displ+ncell /),count=(/ n /))
1813	status=NF90_PUT_VAR(ncid,latid,REAL(lat,r8),start=(/ displ+ncell /),count=(/ n /))
1814	status=NF90_PUT_VAR(ncid,bounds_lonId,REAL(bounds_lon,r8),start=(/ 1,displ+ncell /),count=(/ nvert,n /))
1815	status=NF90_PUT_VAR(ncid,bounds_latId,REAL(bounds_lat,r8),start=(/ 1,displ+ncell /),count=(/ nvert,n /))
1816
1817	ncell=ncell+n
1818	DEALLOCATE(lon,lat,bounds_lon,bounds_lat)
1819	END DO
1820	ENDIF
1821
1822
1823
1824	! status = NF90_CLOSE(ncid)
1825
1826	end subroutine Create_Header_mpi
1827
1828	SUBROUTINE Close_files
1829	USE netcdf
1830	IMPLICIT NONE
1831	INTEGER :: i,k,status
1832
1833	DO i=1,NbField
1834	status=NF90_CLOSE(FieldId(i))
1835	ENDDO
1836	END SUBROUTINE Close_files
1837
1838
1839	function int2str(int)
1840	implicit none
1841	integer, parameter :: MaxLen=10
1842	integer,intent(in) :: int
1843	character(len=MaxLen) :: int2str
1844	logical :: flag
1845	integer :: i
1846	flag=.true.
1847
1848	i=int
1849
1850	int2str=''
1851	do while (flag)
1852	int2str=CHAR(MOD(i,10)+48)//int2str
1853	i=i/10
1854	if (i==0) flag=.false.
1855	enddo
1856	end function int2str
1857
1858	end module write_field
1859

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

Download in other formats: