/[lmdze]/trunk/bibio/Writefield/writefield_gen.f
ViewVC logotype

Annotation of /trunk/bibio/Writefield/writefield_gen.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 109 - (hide annotations)
Wed Sep 17 10:08:00 2014 UTC (9 years, 8 months ago) by guez
File size: 847 byte(s)
Moved a call to writefield from guide to tau2alpha. (dxdys does not
change with itau.) So dxdys does not need to be a module variable any
longer. Other variables of modules tau2alpha_m downgraded to local
variables of tau2alpha, since they were not used elsewhere.

Procedures write_field[13]d and formcoord were never called. Could
then remove int2str.

Inline writefield_gen into writefield.

CreateNewField takes an integer array argument instead of 3 scalar
integers. CreateNewField now creates a number of dimensions adapted to
the rank of the output field, instead of always 4 dimensions.

Changed names of variables of module write_field: fieldid to
ncid, fieldindex to record, fieldvarid to varid.

In writefield_gen, if index == -1, no need to call GetFieldIndex
again, we know that the result is nbfield.

In guide, moved calls to writefield for some variables inside if
first_call: those variables do not change with time. Removed ztau:
computed only to be output, does not seem meaningful. Removed
writefield for aire: does not change with time and is already in
"grilles_gcm.nc".

1 guez 108 module writefield_gen_m
2    
3     implicit none
4    
5     contains
6    
7     subroutine WriteField_gen(name, Field, dimx, dimy, dimz)
8    
9     use CreateNewField_m, only: CreateNewField
10     use GetFieldIndex_m, only: GetFieldIndex
11     USE netcdf95, ONLY: nf95_put_var
12 guez 109 USE write_field, ONLY: ncid, record, varid, nbfield
13 guez 108
14     character(len=*) :: name
15     integer :: dimx, dimy, dimz
16     real, dimension(dimx, dimy, dimz) :: Field
17     integer :: index
18    
19     !------------------------------------------------------------
20    
21     Index=GetFieldIndex(name)
22     if (Index==-1) then
23 guez 109 call CreateNewField(name, (/dimx, dimy, dimz/))
24     Index=nbfield
25 guez 108 else
26 guez 109 Record(Index)=Record(Index)+1
27 guez 108 endif
28    
29 guez 109 call NF95_PUT_VAR(Ncid(Index), Varid(Index), Field, &
30     start = (/1, 1, 1, Record(Index)/))
31 guez 108
32     end subroutine WriteField_gen
33    
34     end module writefield_gen_m

  ViewVC Help
Powered by ViewVC 1.1.21