source: src/htune_plot.R @ 176

Last change on this file since 176 was 101, checked in by htune, 7 years ago

New syntax for csv metrics
Fredho

##############################################################################
# Computing metrics from LES and SCM
# Result in csv format
# Auteur: F Couvreux, R Honnert, F Hourdin, N Villefranque, C Rio, n'co
##############################################################################
#
# metrics are specified in list_case
# metrics names follow the syntax
#
# CASE_SUBCASE_METRICS_T1_T2
# T1 and T2 are initial and final time step for time averages
#
# available METRICS :
# ===================
#
# 1/ zav-400-600-var -> variable "var" averaged between 400 and 600 m
#
# 2/ Ay-var -> integral ( min ( var -var(t=0) ) dz ) / integral ( dz )
# integral taken from 0 to zmax
#
# 3/ nebzave, neb2zave, neb4zave : Effective cloud height
# = int ( nebp z dz ) / int ( nebp dz ) with p=1, 2 or 4
#
# 4/ nebmax : maximum cloud fraction on the column
#
# 5/ nebzmin, nebzmax : minimum/maximum cloud height
#
# 6/ lwp : liquid water path
#
# TBD :
# =====
# 1/ integrals are computed assuming rho=1 because rho is not systematically avalble
# 2/ the time average is coded for zav metrics only
#
##############################################################################

File size: 4.0 KB
Line 
1###############################################################################
2# Auteurs : F. Hourdin
3# Modif 2017/12/15 : F. Hourdin
4# Reads z-t files from LES and SCM stored in DATA/
5#    SMC001.nc, SCM002.nc, ... and same for LES
6# tout_tracer plots vertical profiles of a variable N at time 
7###############################################################################
8
9
10###############################################################################
11# Plots
12###############################################################################
13
14tout_tracer <- function(varname,case_name,subcase_name,time_scm,time_les,pltsu) {
15  # if les = 1, compute zf accordingly
16  # else, it is 1D
17  print('debut tout tracer')
18  print(c(time_scm,time_les))
19  print(c("pltsu=",pltsu))
20  first=1
21  #pdf(file='myplot.pdf')
22  for (k in 1:NSCMS) {
23    file= paste("WAVE1/",case_name,"/",subcase_name,"/SCM_1-",sprintf("%3.3i",k),".nc",sep="")
24    une_courbe(file,varname,1,time_scm,first,"grey")
25    first=0
26  }
27  print("ON PASSE AU TRACER DES LES")
28  for (k in 0:(NLES-1)) {
29     file= paste("LES/",case_name,"/",subcase_name,"/LES",print(k),".nc",sep="")
30     une_courbe(file,varname,0,time_les,first,"blue")
31  }
32  #dev.off()
33  print("ET ON EN SORT")
34}
35
36############################################################################
37# Plots vertical profiles at a given time
38############################################################################
39une_courbe <- function(file,varname,timez,timeval,first,color) {
40       ncid = nc_open(file)
41       zval = ncvar_get(ncid,"zf")
42       var=ncvar_get(ncid,varname)
43       if (timez==0) { zzz<-zval } else { zzz<-zval[,timez]}
44       if (first==1) {
45        xmin=pltsu[1]
46        xmax=pltsu[2]
47        zmax=pltsu[3]
48        plot(var[,timeval],zzz,col=color,xlim=c(xmin,xmax),ylim=c(0,zmax),type="l",xlab=varname,ylab="z")
49       }
50       else {
51         lines(var[,timeval],zzz,col=color)
52       }
53       nc_close(ncid)   
54     }
55#tracer d'une serie temporelle pour une variable donnee (soit un niveau z soit une integrale
56trace_serie_s <- function(case_name,subcase_name,varname,targetvar) {
57    first=1
58    #pdf(file='myplot2.pdf')
59    for (k in 1:NSCMS) {
60      file= paste("WAVE1/",case_name,"/",subcase_name,"/SCM_1-",sprintf("%3.3i",k),".nc",sep="")
61      print(c("fichier SCM",k,file))
62      trace_serie(file,varname,first,"grey",targetvar)
63      first=0
64    }
65    print("ON PASSE AU TRACER DES LES")
66    for (k in 0:(NLES-1)) {
67      first=0
68      file= paste("LES/",case_name,"/",subcase_name,"/LES",print(k),".nc",sep="")
69      trace_serie(file,varname,first,"blue",targetvar)
70    }
71    #dev.off()
72    print("ET ON EN SORT")
73}
74
75trace_serie <- function(file,varname,first,color,targetvar) {
76  nc=nc_open(file)
77  neb = ncvar_get(nc,varname)
78  lm=dim(neb)[2]
79  km=dim(neb)[1]
80  zval=rep(0,km)
81  # rep(x,n) replicates the value x ntimes = une maniere d'initialiser un tableau de 0
82  zf = ncvar_get(nc,"zf")
83  if ( length(c(dim(zf))) == 1 ) { zval[] = zf[] } else { zval[] = zf[,1] }
84  nc_close(nc)
85  # la on veut peut etre modifier calcul z pour avoir soit quelque chose d'integrer soit une variable extraite a un niveau donn????
86  if ( varname=="rneb"){
87     if ( (targetvar == "nebzmin" )| (targetvar == "nebzmax") | (targetvar == "nebzave") ) {
88      zmaxint=4000.
89      xmin=0.
90      xmax=2500.
91      varnew<-compute_zhneb(zval,neb,zmaxint,targetvar)
92     } else {
93      zmaxint=4000.
94      xmin=0.
95      xmax=1.
96      varnew<-compute_zhneb(zval,neb,zmaxint,targetvar)
97     }
98  } else {
99    varnew=rep(0,lm)
100    zsel=500.
101    if (varname=="theta") {
102      xmin=pltsu[1]
103      xmax=pltsu[2]
104    }
105    if (varname=="hur") {
106      zsel=6000.
107      xmin=0.
108      xmax=1.
109    }
110    if (varname=="qv") {
111      xmin=pltsu[1]
112      xmax=pltsu[2]
113    }
114    for (k in 2:(km-1)) {
115      if ((zval[k]-zsel)*(zval[k+1]-zsel) <=0.) { ksel = k }
116    }     
117    for (l in 0:(lm-1)) {
118      varnew[l]=neb[ksel,l]
119    }
120  }
121  if ( first == 1 ) {
122    plot(varnew,ylim=c(xmin,xmax),type="l",col=color,ylab=varname)
123  } else {
124    lines(varnew,col=color)
125  }
126  }
Note: See TracBrowser for help on using the repository browser.