source: LMDZ5/branches/LMDZ5V2.0-dev/libf/phylmd/inistats.F @ 5440

Last change on this file since 5440 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.1 KB
Line 
1      subroutine inistats(ierr)
2
3      implicit none
4
5#include "dimensions.h"
6#include "paramet.h"
7#include "comgeom.h"
8#include "comvert.h"
9#include "comconst.h"
10#include "statto.h"
11#include "netcdf.inc"
12
13      integer,intent(out) :: ierr
14      integer :: nid
15      integer :: l,nsteppd
16      real, dimension(llm) ::  sig_s
17      integer :: idim_lat,idim_lon,idim_llm,idim_llmp1,idim_time
18      real, dimension(istime) :: lt
19      integer :: nvarid
20      real, dimension(llm) :: pseudoalt
21
22      write (*,*)
23      write (*,*) '                        || STATS ||'
24      write (*,*)
25      write (*,*) 'daysec',daysec
26      write (*,*) 'dtphys',dtphys
27      nsteppd=nint(daysec/dtphys)
28      write (*,*) 'nsteppd=',nsteppd
29      if (abs(float(nsteppd)-daysec/dtphys).gt.1.e-8*daysec)
30     &   stop'Dans Instat:  1jour .ne. n pas physiques'
31
32      if(mod(nsteppd,istime).ne.0)
33     &   stop'Dans Instat:  1jour .ne. n*istime pas physiques'
34
35      istats=nsteppd/istime
36      write (*,*) 'istats=',istats
37      write (*,*) 'Storing ',istime,'times per day'
38      write (*,*) 'thus every ',istats,'physical timestep '
39      write (*,*)
40
41      do l= 1, llm
42         sig_s(l)=((ap(l)+ap(l+1))/preff+bp(l)+bp(l+1))/2.
43         pseudoalt(l)=log(preff/presnivs(l))*8. 
44      enddo
45
46      ierr = NF_CREATE("stats.nc",NF_CLOBBER,nid)
47      if (ierr.ne.NF_NOERR) then
48         write (*,*) NF_STRERROR(ierr)
49         stop ""
50      endif
51
52      ierr = NF_DEF_DIM (nid, "latitude", jjp1, idim_lat)
53      ierr = NF_DEF_DIM (nid, "longitude", iip1, idim_lon)
54      ierr = NF_DEF_DIM (nid, "altitude", llm, idim_llm)
55      ierr = NF_DEF_DIM (nid, "llmp1", llm+1, idim_llmp1)
56      ierr = NF_DEF_DIM (nid, "Time", NF_UNLIMITED, idim_time)
57
58      ierr = NF_ENDDEF(nid)
59      call def_var_stats(nid,"Time","Time",
60     &            "hours since 0000-00-0 00:00:00",1,
61     &            idim_time,nvarid,ierr)
62! Time is initialised later by mkstats subroutine
63
64      call def_var_stats(nid,"latitude","latitude",
65     &            "degrees_north",1,idim_lat,nvarid,ierr)
66#ifdef NC_DOUBLE
67      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,rlatu/pi*180)
68#else
69      ierr = NF_PUT_VAR_REAL (nid,nvarid,rlatu/pi*180)
70#endif
71      call def_var_stats(nid,"longitude","East longitude",
72     &            "degrees_east",1,idim_lon,nvarid,ierr)
73#ifdef NC_DOUBLE
74      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,rlonv/pi*180)
75#else
76      ierr = NF_PUT_VAR_REAL (nid,nvarid,rlonv/pi*180)
77#endif
78
79! Niveaux verticaux, aps et bps
80      ierr = NF_REDEF (nid)
81! presnivs
82#ifdef NC_DOUBLE
83      ierr = NF_DEF_VAR (nid,"presnivs", NF_DOUBLE, 1,idim_llm,nvarid)
84#else
85      ierr = NF_DEF_VAR (nid,"presnivs", NF_FLOAT, 1,idim_llm,nvarid)
86#endif
87      ierr = NF_PUT_ATT_TEXT (nid,nvarid,"long_name",15,
88     &                        "Vertical levels")
89      ierr = NF_PUT_ATT_TEXT (nid,nvarid,'units',2,"Pa")
90      ierr = NF_PUT_ATT_TEXT (nid,nvarid,'positive',4,"down")
91      ierr = NF_ENDDEF(nid)
92#ifdef NC_DOUBLE
93      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,
94     &                          presnivs(1:llm))
95#else
96      ierr = NF_PUT_VAR_REAL (nid,nvarid,
97     &                        presnivs(1:llm))
98#endif
99! Pseudo alts
100#ifdef NC_DOUBLE
101      ierr = NF_DEF_VAR (nid,"altitude", NF_DOUBLE, 1,idim_llm,nvarid)
102#else
103      ierr = NF_DEF_VAR (nid,"altitude", NF_FLOAT, 1,idim_llm,nvarid)
104#endif
105      ierr = NF_PUT_ATT_TEXT (nid,nvarid,"long_name",8,"altitude")
106      ierr = NF_PUT_ATT_TEXT (nid,nvarid,'units',2,"km")
107      ierr = NF_PUT_ATT_TEXT (nid,nvarid,'positive',2,"up")
108      ierr = NF_ENDDEF(nid)
109#ifdef NC_DOUBLE
110      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,pseudoalt)
111#else
112      ierr = NF_PUT_VAR_REAL (nid,nvarid,pseudoalt)
113#endif
114!      call def_var_stats(nid,"aps","hybrid pressure at midlayers"," ",
115!     &            1,idim_llm,nvarid,ierr)
116!#ifdef NC_DOUBLE
117!      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,aps)
118!#else
119!      ierr = NF_PUT_VAR_REAL (nid,nvarid,aps)
120!#endif
121
122!      call def_var_stats(nid,"bps","hybrid sigma at midlayers"," ",
123!     &            1,idim_llm,nvarid,ierr)
124!#ifdef NC_DOUBLE
125!      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,bps)
126!#else
127!      ierr = NF_PUT_VAR_REAL (nid,nvarid,bps)
128!#endif
129
130      ierr=NF_CLOSE(nid)
131
132      end subroutine inistats
133
Note: See TracBrowser for help on using the repository browser.