From 23303edc189e2d2f6fb86acae746bde3c1eb405a Mon Sep 17 00:00:00 2001 From: Gaelle DELAUTIER <gaelle.delautier@meteo.fr> Date: Fri, 9 Dec 2016 14:27:01 +0100 Subject: [PATCH] M.Leriche F.Brosse 09/12/2016 : add production/destruction terms computation --- MY_RUN/KTEST/009_ICARTT/004_diag/DIAG1.nam1 | 12 + MY_RUN/KTEST/009_ICARTT/004_diag/DIAG1.nam2 | 12 + .../KTEST/009_ICARTT/004_diag/ReLACS_poet.nam | 301 ++++++++++++++++++ .../KTEST/009_ICARTT/004_diag/clean_diag_xyz | 5 + MY_RUN/KTEST/009_ICARTT/004_diag/run_diag_xyz | 24 ++ .../009_ICARTT/006_ncl/plot_ICARTT_budget.ncl | 149 +++++++++ 6 files changed, 503 insertions(+) create mode 100644 MY_RUN/KTEST/009_ICARTT/004_diag/DIAG1.nam1 create mode 100644 MY_RUN/KTEST/009_ICARTT/004_diag/DIAG1.nam2 create mode 100644 MY_RUN/KTEST/009_ICARTT/004_diag/ReLACS_poet.nam create mode 100755 MY_RUN/KTEST/009_ICARTT/004_diag/clean_diag_xyz create mode 100755 MY_RUN/KTEST/009_ICARTT/004_diag/run_diag_xyz create mode 100644 MY_RUN/KTEST/009_ICARTT/006_ncl/plot_ICARTT_budget.ncl diff --git a/MY_RUN/KTEST/009_ICARTT/004_diag/DIAG1.nam1 b/MY_RUN/KTEST/009_ICARTT/004_diag/DIAG1.nam1 new file mode 100644 index 000000000..a38ba8633 --- /dev/null +++ b/MY_RUN/KTEST/009_ICARTT/004_diag/DIAG1.nam1 @@ -0,0 +1,12 @@ +&NAM_CONFIO LCDF4=T, LLFIOUT=T, LLFIREAD=F / +&NAM_CONFZ + ! NZ_VERB=5 , NZ_PROC=0 , NB_PROCIO_R=8 , NB_PROCIO_W=1 +/ +&NAM_DIAG + LVAR_MRW=T + LCHEMDIAG=T, + CSPEC_DIAG='O3,CO', + CSPEC_BU_DIAG='O3,CO' / +&NAM_DIAG_FILE YINIFILE(1) = "ICART.1.SEG01.001" , + YINIFILEPGD(1) = "ICARTT1008_PGD_15km", + YSUFFIX='dg' / diff --git a/MY_RUN/KTEST/009_ICARTT/004_diag/DIAG1.nam2 b/MY_RUN/KTEST/009_ICARTT/004_diag/DIAG1.nam2 new file mode 100644 index 000000000..8d3c6363b --- /dev/null +++ b/MY_RUN/KTEST/009_ICARTT/004_diag/DIAG1.nam2 @@ -0,0 +1,12 @@ +&NAM_CONFIO LCDF4=T, LLFIOUT=T, LLFIREAD=F / +&NAM_CONFZ + ! NZ_VERB=5 , NZ_PROC=0 , NB_PROCIO_R=8 , NB_PROCIO_W=1 +/ +&NAM_DIAG + LVAR_MRW=T, + LCHEMDIAG=T, + CSPEC_DIAG='O3,CO', + CSPEC_BU_DIAG='O3,CO' / +&NAM_DIAG_FILE YINIFILE(1) = "ICART.1.SEG01.002" , + YINIFILEPGD(1) = "ICARTT1008_PGD_15km", + YSUFFIX='dg' / diff --git a/MY_RUN/KTEST/009_ICARTT/004_diag/ReLACS_poet.nam b/MY_RUN/KTEST/009_ICARTT/004_diag/ReLACS_poet.nam new file mode 100644 index 000000000..95e3e203b --- /dev/null +++ b/MY_RUN/KTEST/009_ICARTT/004_diag/ReLACS_poet.nam @@ -0,0 +1,301 @@ +=================================================================== +*** the following section will be read by ch_field_valuen.f90 *** +=================================================================== + +FORMPROF +different norm-profile +3 13 +(F6.0,F7.2,F7.2,F7.2) + 0. 0.30 1.00 1.00 + 300. 0.50 1.00 1.00 + 350. 0.80 1.00 0.70 + 400. 1.00 1.00 0.40 + 500. 1.05 1.00 0.30 + 600. 1.06 1.00 0.20 + 1000. 1.08 1.00 0.10 + 1500. 1.10 1.00 0.05 + 3000. 1.25 1.00 0.01 + 6000. 1.50 1.00 0.01 + 9000. 2.00 0.50 0.01 +10000. 2.25 0.20 0.01 +11000. 6.25 0.01 0.01 + +=================================================================== +*** the following section will be read by ch_field_valuen.f90 *** +=================================================================== + +NORMINIT +initial values (units are par per part = MIX) +MIX +37 +(1X,A12,1X,E17.2) +'O3 ' 40.00E-09 +'H2O2 ' 0.00E-12 +'NO ' 0.00E-09 +'NO2 ' 0.00E+00 +'NO3 ' 0.00E+00 +'N2O5 ' 0.00E+00 +'HONO ' 0.00E+00 +'HNO3 ' 0.00E+00 +'HNO4 ' 0.00E+00 +'SO2 ' 0.00E+00 +'CO ' 100.00E-09 +'HO2 ' 0.00E+00 +'CH4 ' 1700.00E-09 +'ETH ' 0.00E+00 +'ALKA ' 0.00E+00 +'ALKE ' 0.00E+00 +'BIO ' 0.00E+00 +'ARO ' 0.00E+00 +'HCHO ' 0.00E+00 +'ALD ' 0.00E+00 +'KET ' 0.00E+00 +'CARBO ' 0.00E+00 +'ONIT ' 0.00E+00 +'PAN ' 0.00E+00 +'OP1 ' 0.00E+00 +'OP2 ' 0.00E+00 +'ORA2 ' 0.00E+00 +'MO2 ' 0.00E+00 +'ALKAP ' 0.00E+00 +'ALKEP ' 0.00E+00 +'BIOP ' 0.00E+00 +'PHO ' 0.00E+00 +'ADD ' 0.00E+00 +'AROP ' 0.00E+00 +'CARBOP ' 0.00E+00 +'OLN ' 0.00E+00 +'XO2 ' 0.00E+00 + +=================================================================== +*** the following section will be read by ch_field_valuen.f90 *** +=================================================================== + +PROFASSO +norm-profiles to be associated +37 +(1X,A12,1X,I3) +'O3 ' 1 +'H2O2 ' 3 +'NO ' 3 +'NO2 ' 3 +'NO3 ' 3 +'N2O5 ' 3 +'HONO ' 3 +'HNO3 ' 3 +'HNO4 ' 3 +'SO2 ' 3 +'CO ' 2 +'HO2 ' 3 +'CH4 ' 2 +'ETH ' 3 +'ALKA ' 3 +'ALKE ' 3 +'BIO ' 3 +'ARO ' 3 +'HCHO ' 3 +'ALD ' 3 +'KET ' 3 +'CARBO ' 3 +'ONIT ' 3 +'PAN ' 3 +'OP1 ' 3 +'OP2 ' 3 +'ORA2 ' 3 +'MO2 ' 3 +'ALKAP ' 3 +'ALKEP ' 3 +'BIOP ' 3 +'PHO ' 3 +'ADD ' 3 +'AROP ' 3 +'CARBOP ' 3 +'OLN ' 3 +'XO2 ' 3 + + +===================================================================== +*** the following section will be read by ch_init_dep_isban.F90 *** +===================================================================== + +SURF_RES +surface resistances (s/m), refer to Seinfeld and Pandis, 1998, p. 975, Tab.19.2 + 1 +(A32,2E15.5) +NONE 0.0 + +===================================================================== +*** the following section will be read by ch_init_depconst.F90 *** +===================================================================== + +MASS_MOL +molecular mass (in g/mol) for molecular diffusion, from Stockwell et al., 1997 + 41 +(A32,2E15.5) +O3 0.48000E+02 +H2O2 0.34000E+02 +NO 0.30000E+02 +NO2 0.46000E+02 +NO3 0.62000E+02 +N2O5 0.10800E+03 +HONO 0.47000E+02 +HNO3 0.63000E+02 +HNO4 0.79000E+02 +NH3 0.17000E+02 +SO2 0.64000E+02 +SULF 0.98000E+02 +CO 0.28000E+02 +OH 0.17000E+02 +HO2 0.33000E+02 +CH4 0.16000E+02 +ETH 0.30000E+02 +ALKA 0.50861E+02 +ALKE 0.29160E+02 +BIO 0.10200E+03 +ARO 0.10100E+03 +HCHO 0.30000E+02 +ALD 0.44000E+02 +KET 0.72000E+02 +CARBO 0.66900E+02 +ONIT 0.11900E+03 +PAN 0.12100E+03 +OP1 0.48000E+02 +OP2 0.63100E+02 +ORA1 0.46000E+02 +ORA2 0.60000E+02 +MO2 0.47000E+02 +ALKAP 0.81380E+02 +ALKEP 0.83611E+02 +BIOP 0.11700E+03 +PHO 0.10700E+03 +ADD 0.11680E+03 +AROP 0.14867E+03 +CARBOP 0.85434E+02 +OLN 0.13600E+03 +XO2 0.10000E+03 + +===================================================================== +*** the following section will be read by ch_init_depconst.F90 *** +===================================================================== + +REA_FACT +reactivity factor with biology, Seinfeld and Pandis, 1998, p. 975, Tab. 19.3 + 41 +(A32,2E15.5) +O3 0.10000E+01 +H2O2 0.10000E+01 +NO 0.00000E+00 +NO2 0.10000E+00 +NO3 0.10000E+00 +N2O5 0.10000E+00 +HONO 0.10000E+00 +HNO3 0.00000E+00 +HNO4 0.00000E+00 +NH3 0.00000E+00 +SO2 0.00000E+00 +SULF 0.00000E+00 +CO 0.00000E+00 +OH 0.00000E+00 +HO2 0.00000E+00 +CH4 0.00000E+00 +ETH 0.00000E+00 +ALKA 0.00000E+00 +ALKE 0.00000E+00 +BIO 0.00000E+00 +ARO 0.00000E+00 +HCHO 0.00000E+00 +ALD 0.00000E+00 +KET 0.00000E+00 +CARBO 0.00000E+00 +ONIT 0.00000E+00 +PAN 0.10000E+00 +OP1 0.30000E+00 +OP2 0.10000E+00 +ORA1 0.00000E+00 +ORA2 0.00000E+00 +MO2 0.00000E+00 +ALKAP 0.00000E+00 +ALKEP 0.00000E+00 +BIOP 0.00000E+00 +PHO 0.00000E+00 +ADD 0.00000E+00 +AROP 0.00000E+00 +CARBOP 0.00000E+00 +OLN 0.00000E+00 +XO2 0.00000E+00 + +===================================================================== +*** the following section will be read by ch_init_depconst.F90 *** +===================================================================== + +HENRY_SP +Effective Henrys law factor / exponent, See Leriche et al.2013 + 35 +(A32,2E15.5) +O3 1.03000E-02 -0.28300E+04 +H2O2 8.44000E+04 -0.76000E+04 +NO 1.92000E-03 -0.17900E+04 +NO2 1.20000E-02 -0.25160E+04 +NO3 3.80000E-02 -0.87070E+04 +N2O5 2.10000E+00 -0.34000E+04 +HONO 8.38000E+04 -0.31200E+04 +HNO3 1.46000E+13 -0.10500E+05 +HNO4 4.78000E+04 -0.69000E+04 +NH3 3.24000E+03 0.19000E+03 +SO2 5.59000E+04 -0.48950E+04 +SULF 6.64000E+14 -0.87000E+04 +CO 9.81000E-04 -0.17200E+04 +OH 3.90000E+01 -0.00000E+00 +HO2 3.49000E+04 -0.00000E+00 +CH4 1.41000E-03 -0.20400E+04 +ETH 1.88000E-03 -0.28750E+05 +ALKA 0.15000E-02 -0.32750E+04 +ALKE 0.59600E-02 -0.21700E+04 +BIO 0.38500E-01 -0.00000E+00 +ARO 0.18000E+00 -0.41000E+04 +HCHO 3.23000E+03 -0.71960E+04 +ALD 0.12900E+02 -0.58900E+04 +KET 0.27800E+02 -0.55300E+04 +CARBO 0.36000E+06 -0.75450E+04 +ONIT 0.10000E+02 -0.59100E+04 +PAN 0.28000E+01 -0.57300E+04 +OP1 0.30000E+01 -0.52800E+04 +OP2 0.33600E+01 -0.59950E+04 +ORA1 5.07000E+06 -0.59500E+04 +ORA2 2.66000E+05 -0.62000E+04 +MO2 0.24500E+01 -0.23200E+04 +ALKAP 0.26600E+03 -0.60000E+04 +ALKEP 0.83000E+03 -0.60000E+04 +BIOP 0.23110E+04 -0.60000E+04 + +================================================================ +*** the following section will be read by build_emisstabn.F90 *** +*** <called by ch_init_emissionn.F90> *** +*** CON: flux is given in molecules/cm2/s *** +================================================================ + +EMISUNIT +Emission GEIA_POET +CON + +================================================================ +*** the following section will be read by build_pronoslistn.F90 *** +*** <called by ch_init_emissionn.F90> *** +================================================================ + +AGREGATION +Schema reactionnel ReLACS/ inventaire GEIA-poet +NH3 1. NH3DA 1. NH3SF 1. NH3BC 1. NH3S 1. NH3C 1. NH3H 1. NH3FF 1. NH3I +SO2 1. SO2BFI 1. SO2BFT 1. SO2BFR 1. SO2FFI 1. SO2FFP 1. SO2FFR 1. SO2FFRE 1. SO2FFT 1. SO2FFTL 1. SO2FFTI 1. SO2II 1. SO2INF 1. SO2ICH 1. SO2LUW 0.004 COBB +CO 1. COBB 1. COBIO 1. COAN +NO 0.7 NOXBB 0.7 NOXAN +NO2 0.3 NOXBB 0.3 NOXAN +HCHO 1. CH2OBB 1. CH2OAN +BIO 1. ISOPBIO 1. TERPBIO +ETH 1. C2H6BB 1. C2H6BIO 1. C2H6AN +ARO 1. TOLUBB 1. TOLUAN +ALD 1. CH3CHOBB 1. CH3CHOAN +KET 0.33 ACETBB 0.33 ACETBIO 0.33 ACETAN 1.61 MEKBB 1.61 MEKAN +ALKE 0.96 C2H4BB 0.96 C2H4BIO 0.96 C2H4AN 1.04 C3H6BB 1.04 C3H6BIO 1.04 C3H6AN 1.04 BUTEBB 1.04 BUTEAN +ALKA 0.44 C3H8BB 0.44 C3H8BIO 0.44 C3H8AN 0.38 CH3OHBB 0.38 CH3OHBIO 0.38 CH3OHAN 1. BUTABB 1. BUTAAN 1. C2H5OHBB 1. C2H5OHAN +END_AGREGATION diff --git a/MY_RUN/KTEST/009_ICARTT/004_diag/clean_diag_xyz b/MY_RUN/KTEST/009_ICARTT/004_diag/clean_diag_xyz new file mode 100755 index 000000000..f92472050 --- /dev/null +++ b/MY_RUN/KTEST/009_ICARTT/004_diag/clean_diag_xyz @@ -0,0 +1,5 @@ +set -x +rm -f ICART* OUTPUT_LISTING* OUTPUT_TRANSFER pipe* *.tex +rm -f DIAG1.nam file_for_xtransfer +rm -f DATAE1 DATAJ1 DATAS1 + diff --git a/MY_RUN/KTEST/009_ICARTT/004_diag/run_diag_xyz b/MY_RUN/KTEST/009_ICARTT/004_diag/run_diag_xyz new file mode 100755 index 000000000..5c8a553cf --- /dev/null +++ b/MY_RUN/KTEST/009_ICARTT/004_diag/run_diag_xyz @@ -0,0 +1,24 @@ +#!/bin/sh +#MNH_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier +#MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence +#MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt +#MNH_LIC for details. version 1. +set -x +set -e +rm -f ICART* OUT* +# +export CHIMIE_FILES=${CHIMIE_FILES:-"$HOME/CHIMIE_FILES"} + +ln -sf ${CHIMIE_FILES}/tuv50/DATAE1 . +ln -sf ${CHIMIE_FILES}/tuv50/DATAJ1 . +ln -sf ${CHIMIE_FILES}/tuv50/DATAS1 . + +ln -sf ../003_mesonh/ICART.1.SEG01.002.*??? . +ln -sf ../003_mesonh/ICART.1.SEG01.001.*??? . +ln -sf ../001_pgd1/ICARTT1008_PGD_15km.* . + +cp DIAG1.nam1 DIAG1.nam +time ${MPIRUN} DIAG${XYZ} +cp DIAG1.nam2 DIAG1.nam +time ${MPIRUN} DIAG${XYZ} +#exit diff --git a/MY_RUN/KTEST/009_ICARTT/006_ncl/plot_ICARTT_budget.ncl b/MY_RUN/KTEST/009_ICARTT/006_ncl/plot_ICARTT_budget.ncl new file mode 100644 index 000000000..5c5f5e594 --- /dev/null +++ b/MY_RUN/KTEST/009_ICARTT/006_ncl/plot_ICARTT_budget.ncl @@ -0,0 +1,149 @@ +;================================================; +load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl" +load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl" +load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl" +load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl" +; ================================================; +begin +;=================================================; +; open file and read in data +;=================================================; + a = addfile("ICART.1.SEG01.001dg.nc4", "r") + a2 = addfile("ICART.1.SEG01.002dg.nc4", "r") + +;=================================================; +; Get informations on variable sizes +; dims are dims-2 to remove non-physical values +;=================================================; + jphext = a->JPHEXT + mdims = getfilevardimsizes(a,"THT") ; get dimension sizes + nd = dimsizes(mdims) + imax=mdims(nd-1)-2*jphext + jmax=mdims(nd-2)-2*jphext + kmax=mdims(nd-3)-2 + +;-------------------------------------------------; +; Read data. +;-------------------------------------------------; + +; Liste de toutes les réactions impliquant O3 + o3_list_t1=a->O3_CHREACLIST + +; Tableau 4D (reac,Z,Y,X) regroupant les termes de prod. et destr. de O3 +; Niveau k=13 approx. 1250m + o3_budget_t1=a->O3_BUDGET(:,13,jphext:jmax+jphext-1,jphext:imax+jphext-1) + o3_bud_t1=dim_avg_n(o3_budget_t1,(/1,2/)) + o3_bud_t1=o3_bud_t1*1e9*3600 + o3_bud_t1@units="ppbv/h" + +; Liste de toutes les réactions impliquant CO + co_list_t1=a->CO_CHREACLIST + +; Tableau 4D (reac,Z,Y,X) regroupant les termes de prod. et destr. de CO +; Niveau k=13 approx. 1250m + co_budget_t1=a->CO_BUDGET(:,13,jphext:jmax+jphext-1,jphext:imax+jphext-1) + co_bud_t1=dim_avg_n(co_budget_t1,(/1,2/)) + co_bud_t1=co_bud_t1*1e9*3600 + co_bud_t1@units="ppbv/h" + +; Liste de toutes les réactions impliquant O3 + o3_list_t2=a2->O3_CHREACLIST + +; Tableau 4D (reac,Z,Y,X) regroupant les termes de prod. et destr. de O3 +; Niveau k=13 approx. 1250m + o3_budget_t2=a2->O3_BUDGET(:,13,jphext:jmax+jphext-1,jphext:imax+jphext-1) + o3_bud_t2=dim_avg_n(o3_budget_t2,(/1,2/)) + o3_bud_t2=o3_bud_t2*1e9*3600 + o3_bud_t2@units="ppbv/h" + +; Liste de toutes les réactions impliquant CO + co_list_t2=a2->CO_CHREACLIST + +; Tableau 4D (reac,Z,Y,X) regroupant les termes de prod. et destr. de CO +; Niveau k=13 approx. 1250m + co_budget_t2=a2->CO_BUDGET(:,13,jphext:jmax+jphext-1,jphext:imax+jphext-1) + co_bud_t2=dim_avg_n(co_budget_t2,(/1,2/)) + co_bud_t2=co_bud_t2*1e9*3600 + co_bud_t2@units="ppbv/h" + +;=================================================; +; PLOT +;=================================================; +; interpolation des champs a 1250 m + + figname ="zsection_1250_bud" + wks = gsn_open_wks("png",figname) ; open a ncgm file + gsn_define_colormap(wks,"WhBlGrYeRe") ; Choose colormap + + res = True + res@gsnDraw = False ; don't draw yet + res@gsnFrame = False ; don't advance frame yet + +; X-axis title (tiY) + res@tiXAxisFontHeightF = 0.018 ; font height + res@tiXAxisFont = 21 ; font index + res@tiXAxisString = "Chemical reactions" ; string to use as the X-Axis title + +; Y-axis title (tiY) + res@tiYAxisFontHeightF = 0.018 ; font height + res@tiYAxisFont = 21 ; font index + +; Bar plot + res@gsnXYBarChart = True ; turn on bar chat + res@gsnYRefLine = 0. ; reference line + res@gsnAboveYRefLineColor = "red" ; above ref line fill red + res@gsnBelowYRefLineColor = "blue" ; below ref line fill blue + res@xyCurveDrawOrder = "PreDraw" + res@tmYLFormat = "0@*+^sg" + res@tmYLPrecision = 2 + res@tmXBOn=False + + txres = True ; text mods desired + txres@txFontHeightF = 0.018 + txres@txAngleF = 90 ; text angle + txres@txJust = "TopCenter" ; puts text on top of bars + + res@gsnCenterString="heure=19" + +; plot ozone production + res@tiYAxisString = "ozone budget (ppbv/h)" ; string to use as the Y-Axis title + x=ispan(0,dimsizes(o3_bud_t1)-1,1) + plot_o3_t1 = gsn_csm_xy(wks,x,o3_bud_t1(:),res) + text=gsn_add_text(wks,plot_o3_t1,tostring(o3_list_t1(:)),x,o3_bud_t1(:),txres) ; add label + draw(plot_o3_t1) + frame(wks) + +; plot ozone production + res@tiYAxisString = "carbon monoxide budget (ppbv/h)" ; string to use as the Y-Axis title + x:=ispan(0,dimsizes(co_bud_t1)-1,1) + plot_co_t1 = gsn_csm_xy(wks,x,co_bud_t1(:),res) + text:=gsn_add_text(wks,plot_co_t1,tostring(co_list_t1(:)),x,co_bud_t1(:),txres) ; add label + draw(plot_co_t1) + frame(wks) + + res@gsnCenterString="heure=20" + +; plot ozone production + res@tiYAxisString = "ozone budget (ppbv/h)" ; string to use as the Y-Axis title + x:=ispan(0,dimsizes(o3_bud_t2)-1,1) + plot_o3_t2 = gsn_csm_xy(wks,x,o3_bud_t2(:),res) + text:=gsn_add_text(wks,plot_o3_t1,tostring(o3_list_t2(:)),x,o3_bud_t2(:),txres) ; add label + draw(plot_o3_t2) + frame(wks) + +; plot ozone production + res@tiYAxisString = "carbon monoxide budget (ppbv/h)" ; string to use as the Y-Axis title + x:=ispan(0,dimsizes(co_bud_t2)-1,1) + plot_co_t2 = gsn_csm_xy(wks,x,co_bud_t2(:),res) + text:=gsn_add_text(wks,plot_co_t2,tostring(co_list_t2(:)),x,co_bud_t2(:),txres) ; add label + draw(plot_co_t2) + frame(wks) + + + +;;;;;;;;;;;;;;;;;;;;;;;; + +end + + + -- GitLab