/[lmdze]/trunk/Sources/IOIPSL/Histcom/histvert.f
ViewVC logotype

Contents of /trunk/Sources/IOIPSL/Histcom/histvert.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 62 - (show annotations)
Thu Jul 26 14:37:37 2012 UTC (11 years, 11 months ago) by guez
Original Path: trunk/libf/IOIPSL/Histcom/histvert.f90
File size: 4661 byte(s)
Changed handling of compiler in compilation system.

Removed the prefix letters "y", "p", "t" or "z" in some names of variables.

Replaced calls to NetCDF by calls to NetCDF95.

Extracted "ioget_calendar" procedures from "calendar.f90" into a
separate file.

Extracted to a separate file, "mathop2.f90", procedures that were not
part of the generic interface "mathop" in "mathop.f90".

Removed computation of "dq" in "bilan_dyn", which was not used.

In "iniadvtrac", removed schemes 20 Slopes and 30 Prather. Was not
compatible with declarations of array sizes.

In "clcdrag", "ustarhb", "vdif_kcay", "yamada4" and "coefkz", changed
the size of some arrays from "klon" to "knon".

Removed possible call to "conema3" in "physiq".

Removed unused argument "cd" in "yamada".

1 module histvert_m
2
3 implicit none
4
5 contains
6
7 SUBROUTINE histvert(pfileid, pzaxname, pzaxtitle, pzaxunit, pzsize, &
8 pzvalues, pzaxid, pdirect)
9
10 ! This subroutine defines a vertical axis and returns its id. It
11 ! gives the user the possibility to define many different vertical
12 ! axes. For each variable defined with histdef a vertical axis can
13 ! be specified by its ID.
14
15 USE find_str_m, ONLY: find_str
16 USE strlowercase_m, ONLY: strlowercase
17 USE errioipsl, ONLY: histerr
18 USE histcom_var, ONLY: nb_zax, nb_zax_max, ncdf_ids, zax_ids, &
19 zax_name, zax_name_length, zax_size
20 USE netcdf, ONLY: nf90_def_dim, nf90_def_var, nf90_enddef, &
21 nf90_float, nf90_put_att, nf90_put_var, nf90_redef
22
23 INTEGER, INTENT(IN):: pfileid
24 ! ID of the file the variable should be archived in
25
26 INTEGER, INTENT(IN):: pzsize ! size of the vertical axis
27 CHARACTER(len=*), INTENT(IN):: pzaxname ! name of the vertical axis
28 CHARACTER(len=*), INTENT(IN):: pzaxunit ! units of the vertical axis
29 CHARACTER(len=*), INTENT(IN):: pzaxtitle ! title of the vertical axis
30 REAL, INTENT(IN):: pzvalues(pzsize) ! coordinate values of the vetical axis
31 INTEGER, INTENT(OUT):: pzaxid ! ID of the axis (not the netCDF ID)
32
33 CHARACTER(len=*), INTENT(IN), OPTIONAL:: pdirect
34 ! positive direction of the axis: up or down
35
36 ! Local:
37 INTEGER:: pos, iv, nb, zdimid, zaxid_tmp
38 CHARACTER(len=20):: str20, tab_str20(nb_zax_max)
39 INTEGER:: tab_str20_length(nb_zax_max)
40 CHARACTER(len=70):: str70, str71, str72
41 CHARACTER(len=80):: str80
42 CHARACTER(len=20):: direction
43 INTEGER:: iret, leng, ncid
44
45 !---------------------------------------------------------------------
46
47 ! 1.0 Verifications:
48 ! Do we have enough space for an extra axis ?
49 ! Is the name already in use ?
50
51 ! - Direction of axis. Can we get if from the user.
52 ! If not we put unknown.
53
54 IF (present(pdirect)) THEN
55 direction = trim(pdirect)
56 CALL strlowercase(direction)
57 ELSE
58 direction = 'unknown'
59 END IF
60
61 ! Check the consistency of the attribute
62
63 IF ((direction/='unknown') .AND. (direction/='up') .AND. &
64 (direction/='down')) THEN
65 direction = 'unknown'
66 str80 = 'The specified axis was: ' // trim(direction)
67 CALL histerr(2, 'histvert', &
68 'The specified direction for the vertical axis is not possible.', &
69 'it is replaced by: unknown', str80)
70 END IF
71
72 IF (nb_zax(pfileid)+1>nb_zax_max) THEN
73 CALL histerr(3, 'histvert', &
74 'Table of vertical axes too small. You should increase ', &
75 'nb_zax_max in M_HISTCOM.f90 in order to accomodate all ', &
76 'these variables ')
77 END IF
78
79 iv = nb_zax(pfileid)
80 IF (iv>1) THEN
81 str20 = pzaxname
82 nb = iv - 1
83 tab_str20(1:nb) = zax_name(pfileid, 1:nb)
84 tab_str20_length(1:nb) = zax_name_length(pfileid, 1:nb)
85 CALL find_str(nb, tab_str20, tab_str20_length, str20, pos)
86 ELSE
87 pos = 0
88 END IF
89
90 IF (pos>0) THEN
91 str70 = 'Vertical axis already exists'
92 WRITE(str71, '("Check variable ", a, " in file", I3)') str20, &
93 pfileid
94 str72 = 'Can also be a wrong file ID in another declaration'
95 CALL histerr(3, 'histvert', str70, str71, str72)
96 END IF
97
98 iv = nb_zax(pfileid) + 1
99
100 ! 2.0 Add the information to the file
101
102 ncid = ncdf_ids(pfileid)
103
104 leng = min(len_trim(pzaxname), 20)
105 iret = nf90_def_dim(ncid, pzaxname(1:leng), pzsize, zaxid_tmp)
106 iret = nf90_def_var(ncid, pzaxname(1:leng), nf90_float, zaxid_tmp, zdimid)
107
108 leng = min(len_trim(pzaxunit), 20)
109 iret = nf90_put_att(ncid, zdimid, 'units', pzaxunit(1:leng))
110 iret = nf90_put_att(ncid, zdimid, 'positive', trim(direction))
111
112 iret = nf90_put_att(ncid, zdimid, 'valid_min', real(minval( &
113 pzvalues(1:pzsize))))
114 iret = nf90_put_att(ncid, zdimid, 'valid_max', real(maxval( &
115 pzvalues(1:pzsize))))
116
117 leng = min(len_trim(pzaxname), 20)
118 iret = nf90_put_att(ncid, zdimid, 'title', pzaxname(1:leng))
119 leng = min(len_trim(pzaxtitle), 80)
120 iret = nf90_put_att(ncid, zdimid, 'long_name', pzaxtitle(1:leng))
121
122 iret = nf90_enddef(ncid)
123
124 iret = nf90_put_var(ncid, zdimid, pzvalues(1:pzsize))
125
126 iret = nf90_redef(ncid)
127
128 ! 3.0 add the information to the common
129
130 nb_zax(pfileid) = iv
131 zax_size(pfileid, iv) = pzsize
132 zax_name(pfileid, iv) = pzaxname
133 zax_name_length(pfileid, iv) = len_trim(pzaxname)
134 zax_ids(pfileid, iv) = zaxid_tmp
135 pzaxid = iv
136
137 END SUBROUTINE histvert
138
139 end module histvert_m

  ViewVC Help
Powered by ViewVC 1.1.21