source: LMDZ5/branches/LMDZ5V2.0-dev/libf/phylmd/mkstat.F90 @ 3770

Last change on this file since 3770 was 1486, checked in by Laurent Fairhead, 14 years ago

Inclusion de la routine "stats" du modele martien qui permet de sortir
le cycle diurne moyen de différentes variables et l'écart-type.
Pour activer cette routine mettre
callstats=y
dans config.def. Les résultats seront placés dans le fichier stats.nc
Pour rajouter des variables: allez à la fin de physiq.F et prendre modèle sur
les lignes call wstats(...) et ne pas oublier de modifier le nombre de variables
sorties dans statto.h (n2dvar et n2dvar).
Attention: la routine n'est pas optimale (elle stocke ses résultats intermédiaires
dans le fichier stats.nc et les lit et écrit à chaque appel de la physique)
et n'a été testée dans LMDZ pour l'instant que sur PC/gfortran
et sur SX8/monoprocesseur. Encore du travail donc. Mais la routine est utilisée
régulièrement avec le modèle martien :-)


Inclusion of the martian model "stats" routine that outputs the mean diurnal
cycle of variables as well as their standard deviation.
To activate the routine, one needs to put
callstats=y
in the config.def file. Results will be output in a stats.nc file
To add more output variables: go to the end of physiq.F and use the call wstats(...)
lines as a template and do not forget to modify the number of output variables in
statto.h (n2dvar and n3dvar)
Caution: the routine is not optimal (as it stocks preliminary results in the stats.nc
file reading and writing at each physics step) and has only been tested in LMDZ
with a PC/gfortran setup and on a NEC-SX8/monoprocessor setup. Some work then still to
be done. But it is used regularly with the martial model :-)

File size: 4.6 KB
Line 
1subroutine mkstats(ierr)
2
3 
4!
5!  This program writes a stats.nc file from sums and sums of squares
6!  to means and standard deviations and also writes netcdf style
7!  file so that the data can be viewed easily.  The data file is
8!  overwritten in place. 
9!  SRL  21 May 1996
10!  Yann W. july 2003
11
12
13implicit none
14
15#include "dimensions.h"
16#include "statto.h"
17#include "netcdf.inc"
18
19integer,parameter :: iip1=iim+1
20integer,parameter :: jjp1=jjm+1
21integer :: ierr,nid,nbvar,i,ndims,lt,nvarid
22integer, dimension(4) :: id,varid,start,size
23integer, dimension(5) :: dimids
24character (len=50) :: name,nameout,units,title
25real, dimension(iip1,jjp1,llm) :: sum3d,square3d,mean3d,sd3d
26real, dimension(iip1,jjp1) :: sum2d,square2d,mean2d,sd2d
27real, dimension(istime) :: time
28real, dimension(jjp1) :: lat
29real, dimension(iip1) :: lon
30real, dimension(llm) :: alt
31logical :: lcopy=.true.
32!integer :: latid,lonid,altid,timeid
33integer :: meanid,sdid
34!integer, dimension(4) :: dimout
35
36! Incrementation of count for the last step, which is not done in wstats
37count(istime)=count(istime)+1
38
39ierr = NF_OPEN("stats.nc",NF_WRITE,nid)
40
41! We catch the id of dimensions of the stats file
42
43ierr= NF_INQ_DIMID(nid,"latitude",id(1))
44ierr= NF_INQ_DIMID(nid,"longitude",id(2))
45ierr= NF_INQ_DIMID(nid,"altitude",id(3))
46ierr= NF_INQ_DIMID(nid,"Time",id(4))
47
48ierr= NF_INQ_VARID(nid,"latitude",varid(1))
49ierr= NF_INQ_VARID(nid,"longitude",varid(2))
50ierr= NF_INQ_VARID(nid,"altitude",varid(3))
51ierr= NF_INQ_VARID(nid,"Time",varid(4))
52
53! Time initialisation
54
55do i=1,istime
56   time(i)=i*24./istime
57#ifdef NC_DOUBLE
58   ierr= NF_PUT_VARA_DOUBLE(nid,varid(4),i,1,time(i))
59#else
60   ierr= NF_PUT_VARA_REAL(nid,varid(4),i,1,time(i))
61#endif
62enddo
63
64! We catche the values of the variables
65
66#ifdef NC_DOUBLE
67         ierr = NF_GET_VAR_DOUBLE(nid,varid(1),lat)
68         ierr = NF_GET_VAR_DOUBLE(nid,varid(2),lon)
69         ierr = NF_GET_VAR_DOUBLE(nid,varid(3),alt)
70#else
71         ierr = NF_GET_VAR_REAL(nid,varid(1),lat)
72         ierr = NF_GET_VAR_REAL(nid,varid(2),lon)
73         ierr = NF_GET_VAR_REAL(nid,varid(3),alt)
74#endif
75
76! We catch the number of variables in the stats file
77ierr = NF_INQ_NVARS(nid,nbvar)
78
79! to catche the "real" number of variables (without the "additionnal variables")
80nbvar=(nbvar-4)/2
81
82do i=1,nbvar
83   varid=(i-1)*2+5
84
85   ! What's the variable's name?
86   ierr=NF_INQ_VARNAME(nid,varid,name)
87   write(*,*) "OK variable ",name
88   ! Its units?
89   units=" "
90   ierr=NF_GET_ATT_TEXT(nid,varid,"units",units)
91   ! Its title?
92   title=" "
93   ierr=NF_GET_ATT_TEXT(nid,varid,"title",title)
94   ! Its number of dimensions?   
95   ierr=NF_INQ_VARNDIMS(nid,varid,ndims)
96   ! Its values?
97
98   if(ndims==4) then ! lat, lon, alt & time
99
100!      dimout(1)=lonid
101!      dimout(2)=latid
102!      dimout(3)=altid
103!      dimout(4)=timeid
104
105      size=(/iip1,jjp1,llm,1/)
106      do lt=1,istime
107         start=(/1,1,1,lt/)
108         ! Extraction of the "source" variables
109#ifdef NC_DOUBLE
110         ierr = NF_GET_VARA_DOUBLE(nid,varid,start,size,sum3d)
111         ierr = NF_GET_VARA_DOUBLE(nid,varid+1,start,size,square3d)
112#else
113         ierr = NF_GET_VARA_REAL(nid,varid,start,size,sum3d)
114         ierr = NF_GET_VARA_REAL(nid,varid+1,start,size,square3d)
115#endif
116         ! Calculation of these variables
117         mean3d=sum3d/count(lt)
118         sd3d=sqrt(max(0.,square3d/count(lt)-mean3d**2))
119         ! Writing of the variables
120#ifdef NC_DOUBLE
121         ierr = NF_PUT_VARA_DOUBLE(nid,varid,start,size,mean3d)
122         ierr = NF_PUT_VARA_DOUBLE(nid,varid+1,start,size,sd3d)
123#else
124         ierr = NF_PUT_VARA_REAL(nid,varid,start,size,mean3d)
125         ierr = NF_PUT_VARA_REAL(nid,varid+1,start,size,sd3d)
126#endif
127      enddo
128
129    else if (ndims.eq.3) then
130
131!      dimout(1)=lonid
132!      dimout(2)=latid
133!      dimout(3)=timeid
134
135      size=(/iip1,jjp1,1,0/)
136      do lt=1,istime
137         start=(/1,1,lt,0/)
138         ! Extraction of the "source" variables
139#ifdef NC_DOUBLE
140         ierr = NF_GET_VARA_DOUBLE(nid,varid,start,size,sum2d)
141         ierr = NF_GET_VARA_DOUBLE(nid,varid+1,start,size,square2d)
142#else
143         ierr = NF_GET_VARA_REAL(nid,varid,start,size,sum2d)
144         ierr = NF_GET_VARA_REAL(nid,varid+1,start,size,square2d)
145#endif
146         ! Calculation of these variables
147         mean2d=sum2d/count(lt)
148         sd2d=sqrt(max(0.,square2d/count(lt)-mean2d**2))
149         ! Writing of the variables
150#ifdef NC_DOUBLE
151         ierr = NF_PUT_VARA_DOUBLE(nid,varid,start,size,mean2d)
152         ierr = NF_PUT_VARA_DOUBLE(nid,varid+1,start,size,sd2d)
153#else
154         ierr = NF_PUT_VARA_REAL(nid,varid,start,size,mean2d)
155         ierr = NF_PUT_VARA_REAL(nid,varid+1,start,size,sd2d)
156#endif
157      enddo
158
159    endif
160enddo
161
162ierr= NF_CLOSE(nid)
163
164end
Note: See TracBrowser for help on using the repository browser.