Context Navigation

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

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

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

Download in other formats: