Context Navigation

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

Last change on this file since 156 was 151, checked in by ymipsl, 11 years ago

Implementation of mixte parallelism MPI/OpenMP into src directory

YM

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

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

Download in other formats: