From 1beef09759d0453a520b24966d63e869413998d0 Mon Sep 17 00:00:00 2001 From: David Jean du Preez Date: Wed, 23 Feb 2022 16:24:18 +0100 Subject: [PATCH] additonal chem species --- README.txt | 8 +- .../augsburg_validation_summer_10.conf | 7 +- dynamic/.palm_wrf_utils.py.swp | Bin 16384 -> 0 bytes .../palm_dynamic_config.cpython-38.pyc | Bin 1122 -> 0 bytes .../palm_dynamic_output.cpython-38.pyc | Bin 12200 -> 0 bytes .../__pycache__/palm_wrf_utils.cpython-38.pyc | Bin 16169 -> 0 bytes dynamic/palm_dynamic_defaults.py | 3 - dynamic/palm_dynamic_output.py | 768 +++++++----------- dynamic/palm_wrf_utils.py | 56 +- 9 files changed, 310 insertions(+), 532 deletions(-) delete mode 100644 dynamic/.palm_wrf_utils.py.swp delete mode 100644 dynamic/__pycache__/palm_dynamic_config.cpython-38.pyc delete mode 100644 dynamic/__pycache__/palm_dynamic_output.cpython-38.pyc delete mode 100644 dynamic/__pycache__/palm_wrf_utils.cpython-38.pyc diff --git a/README.txt b/README.txt index 93250b7..e2e9722 100644 --- a/README.txt +++ b/README.txt @@ -113,7 +113,11 @@ origin_x, origin_y origin x and y of the domain origin_z origin of the domain in the vertical direction -# Changes made to wrf + CAMx scripts +# Changes from wrf + CAMx scripts made to wrf_chem_for_palm -Proj future warning resolved -Chemical species are read from wrf-chem files and interpolated at the same time as the dynamic variables --interpolated files are saved to a directory directory than the wrf-chem data files +-interpolated files are saved to a different directory than the wrf-chem data files + +# Latest update (Feb 2022) +-any chemical specie from wrf-chem can be selected in the configuration file and will be included in the dynamic driver +-extensive changes were made to palm_dynamic_output and palm_dynamic_config to accomdate the changes above diff --git a/configurations/augsburg_validation_summer_10.conf b/configurations/augsburg_validation_summer_10.conf index 7fd7b66..d31a3b0 100644 --- a/configurations/augsburg_validation_summer_10.conf +++ b/configurations/augsburg_validation_summer_10.conf @@ -8,7 +8,7 @@ scenario = 'validation_summer' nested_domain = False # dynamic driver output -dynamic_driver_file = "/cfs/home/d/u/dupreeda/MBEES/PALM/palm_model_system-v21.10/JOBS/augs9/INPUT/augs9_dynamic" +dynamic_driver_file = "/cfs/home/d/u/dupreeda/MBEES/PALM/wrf_chem_for_palm/augs9_dynamic" # import grid parameters for dynamic driver from static driver grid_from_static = True # static driver input @@ -36,8 +36,9 @@ vinterp_terrain_smoothing = None # interpolated files interp_dir_name = '/cfs/home/d/u/dupreeda/MBEES/PALM/wrf_chem_data/interp' -# PALM species -species_names = ['NO', 'NO2', 'O3', 'PM10', 'PM25'] +# WRF-chem species need for chemical reaction +wrfchem_spec = ['no','no2','o3','PM10'] +# possible species: no, no2, no3, pm10,PM2_5_DRY, o3, co, hno3, ho, h2o2 # radiation radiation_from_wrf = False diff --git a/dynamic/.palm_wrf_utils.py.swp b/dynamic/.palm_wrf_utils.py.swp deleted file mode 100644 index e661e7fcf27403eaa77707dfb4d022a6aefef2c3..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 16384 zcmeHNO^oD56|S%(AzAXHz%`<>t(D!gZFkS?h77~@ZYDF6?2?_C%}h@=(Pk;R?P|As z{fF(I>1G81X%R>fhbUaQff5ctxN$)U1X9ES@sk9!5uUxvKeoCQ?jVQfW?Mq3cHlnxzI(Q=qhEMbSQ6M7FYi&1}X+B1}X+B1}X+B1}X+B1}X+B1}X+B z2L3M@Ff2lzg6FDOFXR1xzW)E?U4*;>ybL@I90mS#Cm}BY&jVisJ`a2xaDh`m8#oFy zfxq5C$j^c2fFY0o8^B}0LEs)>AMnff67n_RGr&b41|r}wpavWOZU%n6pO9YxF9FX2 z=YToj^6iBD9e4?N0T=@#;3RMmxC!{hZG?Oi_yXVnePA8Xfc?Pjz-#Xz9KL9QOtH3GX1Hk252>BWCD)1HH)4+qk0pM>p6Y^u=N5G4K4J-l&fj8bw z$m_uGfNufM0w&M}?g!ove19Jy-vgcldca-4oxp1tr*8w#03QWT0%GhG;GI0EJS^I# zeyI1yVUVbUfi`GR$I+ZR$rT(28u&Rx_4|v$L8=Bp%-x>u8i^bDED6}aQE}0rxGY-` zx>+xt9XZz7)nRUwedYDd*0%CNzy? zFU`eHu2{4gMvkAT3Xh^7qWCg}4jk>6>n>2w7|{!@Pf&~74Qg!L^ula4PKW=2@`>zDnXgWVV_65K}2P72!AH+%R*tdLs^()vkVSRUXrvK`6C@t$69JC zKEG;Eb)*Rw>P>4jrJL8wU6^(ZE|1aE9!{n6UAh5FCVVWK7fwyuSSE&rLZ?QTdC&1V z?HO@G3)YOpL9lo{=?7L^4y$kn=9t`1c!X1p@=<~dHlCSCO~TJoYgB5DCbeo*j0g=P z$98;b20>&wzLD^l;!-+P_AF>Jk|a`Gu4Fr3Isc)4;X=c{W zo~Gp#z-84JoXpng47J7u$5KIXqWXf4OreZ}8nv7_F?^HLE>Cti2OlZLF%5dik8zxF z6$VJ5G2EuyviZcg$ZU=dU_e{)vbnz1IiansZ*)$uM>b9_v*nYkD{Gr41$V-T!oIAw zFb~${wa=_;k7Si-dQCmXlrCzUeTGY!$bRobxp#Dqo- zc_I$nG%yh%UfH^cjU?Ey#DCoGL_=P7BUFucWLVBHZqROiMnrjj zi2TrMXEttj7T1<#xRG2z{Kvd1lvdXlmsuK?A_~i-lIzXM#yT5YE2kDab9?E`PV3B0 z>P*gV5G0eqX{0$-;_nqE4?aBUY!tp8oY0V0JC85KZBl6B{ zJWMSk8G6)nViXafFOIY#cik{CNGIgIP;)_{LvaI(L&n>&j)1$~$6BABE2YMBCEQ1;g)aFSVF}M#RFJ%EI z@-JhK?<9;}d7D(0ij-I8T}&kLXat5Zd1=qVgO@@~S7ZQ`r8PRmoh|4NxS5FCF2f^K z`Vi^!9b@FcQ1Ov(*N8bwxEDfpxYg2WZEd|nAqj>F?X+4XvY07^b9(auNPAJ>rGvA5 zgeKhV`$58~y1^{kq$UgN&=<<{Lf4AM)Nmqd8!>XXFS1uc)sHmK+s#Mo#D?C7^!Zs@ z^P2G-$>tl%K6v;6eA3*~Tck#pWUUycf0aEinU&<7r(OmVUb| zKV)X;vgg>;@xmZV=!{gs{VE*mLlBOeNG3`zvIa6KGANJ20O7-PorxkgE1nSo8>3dJ zp^V7ZvL;N!H7We|LUBgYLWZZI823zuY>k@m6)mWlJ}ZPwk3-kNQ?QQ6o~`3Sx<1u3 z7E_to2+F19hmYi!vFX4u-*P;|=AvtPn#;44kY6DqhZr4bg?SO4qWk%IfbdWhMa)9= z5Pr3#6DzBADn1E!k}>rHYv|&60AVP;CP&6jza6JmYTDa+)DRY% zX5@qzKhZz3o5|zY@ojN?$xY3sg@}+n?9<7e;aujPfPw`BG&((!W#J~*3x4AugXQ_U$EZ)Gw>?#UEm^c9C#CJ``-e;0e%R4 z8Tb;A0Gq&Z;23Z>@Dr@>KL>mkuz};iy#N6&VNL%H;6>mWpbI<*h&BIZ*zz*)Q{X$m z*MW;b7gz&?tse%eqhg?9pkkn6pkkn6pkkn6pkm-Z$w0nT@9Cx+#9S?RZW3MYH}@vW zu*_OoTpjIN7f)#Ppdr^%`m&K=L7b>YWXHN`c%!T#r$>vDa$mGGpx6r#wxusIn1$tT zqPlz=+Y{VPjI%h!N=_( z$qTXHie)CVLls*KVSmCqT!@lBT1ZvGKI@OW*r^bj-mtVtzI%Df-%1x#nbHC_*NVs`>Wk8t)2ZEK+uJu|u3*oY9Up9ep(&jCWe>|uF1aCir#UIU{j zK@?M`z=;F9Tewp%@o`;I_GH6Yc38&IpX75W}1tVOa-kAHjJWGSX5T4 zNG2IkqLw+MgqKhEM^h%L&Auh8*4WkFR2s#;Hdny{+b;`rNU;v?9NgKh>3}+O=N1F*&hQd$qqsT(TYR-X$JCo4>hGbs zd+UI0Z{{r71~h;eF5uzJZSjz<-Qh)!{>+(U^&Pa)%-KV9eGj#|4cc$_P5lTmayhJm{eIYP3cXUO?gg6QcMq;{bekrY?Rf7R>SJ1 zdHk_EkUUe234%bM!aC8*(OvlMM|LwyXH!(7Y}y$20O+P%6Jtt+XR$Lcn>o{cBwjqm#Sj9 zHo=%_8_4wYT1r;xkp#nJ9XwK%OI~SZ{7S-P8ShH)Qrh6ASHai;)Z?mk_8>JFuxpkZ zPgJZlo@5z?=-v1YJQH;#vDZj;I{9BAZHQ!{Z{;nuc)oqTMZO{V1oHkli}R6+Ct}KC z8rLzcE6Etm;^PlTN9XbB`%lO57dcAv37e)PAyZr2c+tyvrPAFEpw?!s6Kl!#A%%Jl z47ty7~;}V+y9804mHA~>X4ZkK z?tVHvoKfKK>VN*ty7MoJ@*h|c|1_vvgRgrT3ZfVaQHfSo4V7IrLt|InNU&?tNWxVw zrz&Y9ttuozl7CW+0ir!tW>c*nGf$RFZnjV=+u6KhTb}J^tt!c`IVB{$WjjkYs#c|x z&l0E9w4H3#_VT2_HI^*bCfQmwTX$^yw;Y>fJjbfKg_=_V5t=GhJ=>|*%3vGq7Cy}{ z+ZCqop33|NAt?7}Rewj_DpyPpo5i*KPUHrEzrb}6zBl3P{&}yWco6S9RaHEd>hIL{ zPNESVqarDiCIiKEdz=iC4B0`3$S@fpqrlnaX`4D7qzN)Kt<_)fz2gWjG1Yhz$( z;+^_V1DlA;0B!_u8Ng)#mjT>v!0iOw?f{n{^l7DifF2+-MENHT{5U`qvJa^#rJaRHKcd+n zeJXr?R}P6N=r8Dj_F;N(bMym5j(swr zko|}YIMAD%1HE*Z_bx~S|Q~ zqCH3FHg}R7+0REKNDm3pc{(pi4}fnq z^MUjNy&y@C3DO0+AW4sd^lYGekzQod85R-J=g9d^bh85T61^lKkrt6B1mtJ*GXaUT zi2OuAUZ$4?B+}=|g-+`y1>_ZaML;4QBFEkjyID3L7XXhd!&P#z{W2(QlB|uyaY}}wXB*fw%3$y4XdK8!_JAqPMyz z7LV-hZF*Zm$liWIzmO2Jw_nmPC4}tl9eSsWAU9*Y{fd6oMMZl1HT_ybgx)R&dAUpP zN{GHRK(+?KrkB0xQ$54xxxZ$rdG`cOiI-i9OL8~Tle2)zvv zkLV)_5qkS&;KyV7xQieYlDBsP)D!xoi|X+2nDY zdK-@07xYCJLB5u}4QJU)`m&2+v^f9eF>FvnKqD=p!zjXxeoMa<&`68uFrpCsioO!i zNQ-DWy4UozfJQn%le;mIen-C(G>{hU%aNMY6wpYEXgSjF>GuK}X%Q_)YEer-BQ2uk zNSEl6fJQn%lY23d=4oEgKw7jfM@lFW&`67DIZ~V20vc%%Ek|0Q1p$q;h?XNQ(xQMy zIzW^AF_A9QWkCb!@RtU619GG#S`yGmi)cC06}lpzkrvT%q;KdO0gbeXmLo0GvVcZ9 zK$8bCkydC$&_G(WFGpIXRRN8(h?XO*(VBopT13l{)@fZpBQ2ukNLT5qfJRzGKZM`= zZ-Ot;GFwHp=0l7+<(O^!l@5V<#=j+pLzZFb^}g8PpL{9;t2?A z!U-tTdc@8ofBXcr0H>HH)BQu|9yrfV(k72YEb(;9`UTFXFdnbc2K0+Oi>BEqxJDbj znp=_l97<-U6sB>#uSPIa8S+LyBtO-nRo*_+xEZYh*g=S>gB3IuBQ)gM&m+CHJyMGF zc3%xK(l7cUBP0DXTBA>-cOo?8dC4Pvxjj;f)aa`rM*3|(WMrhTqBZ(N`Z_{GUR8Lc z-))bSA~pMJh>?EZ4;dM$6|K=H(xnIuc}?Vz=C?;mk&?a|Vx)FIWMrg;XpKIR79%v| z6`e=AyggEiwA5EajC7?RGBVOP(HeasEk|g`zehaM%JxVp(rRA~G16K;WMriEXpKIR zu109c`zRJ^dyPu`-niQ9w7Z)UqG)e)`B6T%TNL(GAjw4kI<()?HF|JDVRu*8LuyK) z8}tz8!hP0ONCW>ZJzQXbr&j{(^Qj?nh{JdxPQN zE;sILA8e0=cHka%XRq)uLSut|6TxlKM-g1`L-3sY(%m&bj?j?z*xd6c+ao2@(~jvkiT8V*U*LXi zhP>cBe!mxYYFpj?QMxP(lVKK!0#6cO-*4p}I<>s*x}|Ed_n99)T!cr7 zHR(yBueN-xwdN;VZmm@I)pcK6Z}{qluWdB_M2nOvAX>xA07D>^SZD0@2CHna%0?3| ztx+iKs#}7Gk+nj0{gTAkxKzaFi9h}GzyGy0w8)=iqF-9=2`KFkw%==ll;S{XAZ9;eEVL+N;9iCFVKttEo zJ!7owb^v9&;F;jW1UAue-`FWYoKkVwGkaKWt>G>As$sxAm}<=dj4>)$@@n;tQYSt` zgML>VJw61PVNNk1GJLh!1HzPr3^9F7tCPbJ!5HYIufYJ&vBS06ypg~THInrTbf7UT z_}VDWYTc~#h?*$c1E!p5bV2nK)_TcpjS6!C=U6Lu&vM-?DY;8l*{bI4OWC)lFD(?? z*~)y*%^tU!#q5>r+cTV%b7xy459)Rm-xc7m85jdN8%@U|#vl}?yJnWEt;F+OWsWn- z$S{UGGLw`FvwQty-GP@6e7$T}A=tW~vf$N(qU|ROU!t>*eE$F9_xwx5RA8gL2l z(~Fh|Z#8&+GGDH_wm-_>idYO+7k`8^o>@-GS}NPV-YV67%`suIvK_FfyN!xJ*bKz3 zKM2OEdCT)0mPJ45mTkN4Cj(d9Wvg!ciHdDi{j>!ySXf2d*@Zd9b;UW3iT%YN#1|vj zI_eMc_bGzfuHfwqlfRw8xc(rf%Qam)?`Jy8q6^EYS+<*Y+4YA*%5*)e2#YNYNMy2Q zBfw|`7G$zvWXiRwnYZyJ4xlYBp#2fdHgBj73%@^#bvD6W6B&l?RxE$GYP;~>he>J` z3%Y0k++jLCQt)1h?Ha>CxuuqEvMQVVY|&|<9-4aR`gUSRBFPC>F<{FeaHxTNYL?MzGp6 z_C#tlV)$|RHDj!6hvf(sNQqO4aiUr?rn*q!*9}&e7!_D@rXy9@Iy4g_!o}$*lZ6X< z>;j*LDTpH~XkS>->K@F4f)}{lD`0DUq=YT?krKAR2ifkm2GeS5^Q+C6=0OnCJV;`i z2T4rxAi?JO;1Jwu@YeY->*4dpBLuwd2mx<;H|vldQ;{tPac{8`h8-b@X*mpfOUp4a zuf@c?79@O2%WFZtBLp!qud!@z!B=ZBoehGR&IU=aW$A37CJpeLmMD z33nh;#65@m$k$87qkY93$g z^tDFO*Psm@R_4Db&JwioUypM6f}2~eRqPzeHFBf@8+4mkxqCNm-FlpRc>V6Z9BzH^ zl~LTPoA{P!F4*(r*ej#q1_<^B^@ewXe!RoeB0GfO7nIO&kA5A~^xP;G>Ud6QG{~WarbZ zB{Hh^dkU=o!K)tPhrkxMunm<+X#;>v=n3$hYix7X{*NGqh@h9kZ1ksu*#2MW1ml5? z1wEtBsol?nc0)U<^3Ul+MoYnX+%hJG5qVNO0owSTi=y&6O&(I!qkI*XCy6w*-GBDf5R1N0Xm^z`2=s5dEfErha zz}6J-#p^!ZAuy|F2mQ1Or;NO5Iv1cGoVBWrN`1pQi;Rm<_#>US%o|>*?D8`f>a*2I bDDsFNxLm1`M%liKr)qZy3XF$;hP3|$MHKx1 diff --git a/dynamic/__pycache__/palm_wrf_utils.cpython-38.pyc b/dynamic/__pycache__/palm_wrf_utils.cpython-38.pyc deleted file mode 100644 index 740171a90b171750c5328ab6d77e2e90f94865cf..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 16169 zcmbt*X^8q9RF zTD#(q1}Ga^S_;{g6d}#XwrtuChhycHWrrNWa47Uk@(+cyUYCAYid}m}du^{QJ1kQ& z;V-kQdkz50wT(e#<-78|dM`6yzI@rYySh>e{(Woq*R$XFi;D7XN^JfYN8&sl_j65A zgkmW|6hzT(%ro{APT~#=JLYx$* zTIscu>xy;CIaNFj8o%Wv)?OlfUA1PMm&NH-6>au9GnmI&XLjusXJ(ypUKTHL$xQK# zm>K1D5ZVWtz~2YjfeW@vy5L-J&W$M6Mdy-pv3NP;h0Drc8nmuBFF046i_R5rho&kYY4`7f z)>)=?;QyP}KXzu;uC6PXqf^D#XjP&aeT`@I+EX(b)fU>kk~!Zv|MtzR>BC0>k6*fW z<$5N)EM`pdM{FG#CE>Oom zRe@x*pu14<99gTDvYum4iJ-q!&D&e^GHPHhRSKR9L_Lxx3jg$tI~TsY_k;iRrC78i=Dnj1>}LOUL5sU{^yJ!o=GY zF??`5=)c0jIrIF20ez9kt)QeJj0CC^>G9`fB~_078$B z#~sACKa4G~3VA=#)I6oA3f+imiXBAk`hEj1Q@8&#EeVM-(nMn>up(Hs!>chOblkBV|x&n!drf3Owd1Yw<_s zBck8WJrAs^>&BP0X1tm3A2^+xUlA87~NFP!$uY9b|fg9p^$IE)2B*w@B z=#|{U$I55dCvpq!#A3A!-7~Q?A(m>=aYT0F`o$|(<|l|hCzj=com+IuHdf89WlQCW zh+`)xVN)@&wmV*{2Z=MKZ27Lpo;!>o{0SbVaqKC!7`MT-1<$s}&y=dUY{@+bBy48I z!Q7>4CEqY7kKylD#te*_tQKV#Xa(xsAnsL{^X|;aAYQ53!dY-Fy_CyYiCZ^sUwPxg zt*JIPZFOFM^TK8O(v>+I1y*wY)`hvt_8V``Ws;Jn$uiNGJr%yD!`EbBqIY|eyOV>c zC3g@QqTaR3cksN^7e$slLK&rGYv1M><$h-3c~+ zq!Cl@R2Hz_l7#1DEr`|XG;=|+EqNIu(DBUx)k>hQ1Zw@`SfGg&JoUgROxi1fQJl2v zK_@MTBRxAW3nJ(W*Ti;fi;i?;56H`20IEidWqVb-SeQo+bDp9mR87&IB(%ol53}~0 zi99t6aZuel3&tbUS4jXLX~W>~M-ab9`tYV%@wN6W$en18uU3nJT9A8@kb4O1Be3Z^ ziUsT0bj6I(Pwk$R=EwTJ1JnjtCnXeXMjKFZvMwO< z8!%N+5Bk=CURKsjU-R{K7!u7`Gw!KtF+b)TBMNEa5oKL#CcL;X))IbVExE1$rq+>5 zW&_%G_=azC+hjAvZIeF43%VEL7PU2>*_K8e+cM(#^G7Ta!^#gDwI3JBVdc4_Pd#e{ zyfbvn1;!m{JV|>(VBUeod#!QEkrJCOZM>z@@pHw(Y`UCXNxwVUOqY+Pvk&s=cc*}k zH;!&pEEbNRE0s1X3Mao&6^U^xM*Y~gtvGjL-)03U$Ix%%#0H()$6YU*&qH^JAKZDS zGe|zj%0l)oYy!wfEilWZftde8Ng!DX%*A@G>IG&^!uknfWNMIB!v_UM*~wNSVUkm* z9P&3yJk z9f=_2>{EN8_Xbt$4%^M3@sPw^CKf=%VC*@z9i(hKGz5Wn+V=gWY$+;<+qS6YY+If{ zZ}L0=9?k|Ja zTmL$N(5}FW;sAJ@_C{Y_SDPwh8Zc59Pz$!^=}-<3Y~TiAzE;%kY3__?zDKoVKGlZN zM#`ay^C+o@CRD{*T%d#uF|I|_sv zcv?}%wxbpep+DBJg{gx06JBD2I-@Pd|V1Q2q!3=04LY*Sxy?QOclo=H-=d z6ZRGXGA?29JImBuzCl@JrvZrhN+v^=Q(+blFR~Fo>`1gzCz$80n!Ze zECF77J$u)M-Fshp@NpOd-T*Ka4Q12(kQ@h44A}2&K_mP%WVn>8KteTGj`Y{5^gfXl zMPG%5XAELe2I+LI&B?2%nb8AHifEY5N)4w` ztp-(>uLXLZ5O2~l?`CDSynU=G&b6Vv++PPFmv=BI9Z75Z{NWMHg@bDk4fS65-JDwNFqf8dyS=D@bi>LO@ z8|}YB0!m@CMNX-^Espn0mNvJ%5P%n2doVIwoP z)q=9b_GfO(eRkW}wp!el4uzTcTwva?Z@qbg^^v7tI&l@?B!Lsx2)p*NCQsvA1KrN9 zz^SD!1X_71&})mHrRPp80GuRnVv(@L?Q`1=cSg>xE!vCqyRe*o1C$~zB^xpPJENX? zBEN!?EwfD=+MezAkodot?H_JNX}j4ztIwW4)a2PR=VTotiRQ|hd)v~jx!aQVQZSrP zMNI-f24I=+FWC#>>#2F^YR!vV$d|LfjruO>uh5w8R*$Ghp2)8vy~U?RTCP1M)ak~Q zgbTbGE`a&CoLO>nNr+Q9_BQB9Fe5YSv+H;&c0_GW7ubUb6%Th0tPo;w#D+{0Xg4kf z>aF}AFMj*cH~#cH=P(n@FkR#bBlE-du$e{MTCKY5z*>^R&egCb@{}#(OSbckNsYbV zMM)%(msR-?sVyu*bJ5oNe;vOivd&TN$uLRyFw0 zEw*+|4>m8&3T)g+c%oSYEzZNyck_IO)(7f7BYzckr>fk)=*4kF1rD(@=K}rOYvGBM z?SoHYLp*^Wj%QBI^S5G}vK5>X$uxjvS*tk}!CHyuQrSVRQVBda&(qz{W}jNiHEm~zXNNI^fBfx zg|~rs2i_*$opeVd^t|5eg0h&oEZ{ZxrjUADsdTEiOHh?}l!|so z@$@yQNdZSgvxjkw@m^0a!VVNM5f_Qadb7{(TQgyO!SR5dI)xP*YxYx~-;Z(h6k~or z-HcJ`{h##W$b)-(RYkp#<^b30@&~Yn-2}PTucB2#bU-~?yzBYpkfK3< zP;{b3x0l>WHsTH+BfOXRNR=o3cyp+-$4iMGI8u`97z<{iv)BdKOAQ6W|~r{HS}}~e~-WOIVGcL`GTKrt(zai`1koQY)S9;)BYgx4`7u>cqMieNBsRo z%mC(X6kCJJsZkm&=~Se~#27}1yHNkd?em0_$bRHtBz525`ijo3d#eAGPT?eca1MgP#M-aHIR*-baNV(+@% z%=m}>%+}tb^)T=?#9!%1Kad6Fql zwkc0Ca-9T$}nirhYD> zelDW^DpNn(rvB;=qQ175q7nSw=-n5-`TTwk_Gu%~8q0y%aI1w;&WGSUAlbrWJcxlcwpRu>lUZnwgY!w(BW2>yhXcGmENM&OSEt?R(GUp zlNTS>?=?OT1p19@7x|FfxHVs{R=veSB_9SfkYSdEAD1p0$m?5MLePTiq^k?tV;zyxg;>z72G&;y((+m2^{R{Rv(9{9h?(W%}m;>Ztef@ zyVG2DqM2!QwHmmKg#}oO*Kc18wClIS^Ua5t8hIC>G4R@rH2WUh^xb+I=A=`a4KFq5 z0w9k^QD(mONf4jEaplsr+t-7BzRAvd-p3Z#`sBb@=yWHyar>L8Dua%kQ-%6QU^~BRHswZ%fdVOj#|1;k__(vV%ednyf?QLu2sYxq&DgmvVql^i;za2Y=M zK$aa@bwg+AkD?yDQMM#kQqUca@vZ2#EHFw1$gY_yLC`YigIJ+Lfi*!KvP_5-#{rXb z5)GT#+6}r|SSZ*?hQdxv?*@sgEaVFnyMc&{P@aPxI;4h<4Wi)&qmhciiI*kbc3X-d zFdCQ*jDMM4306sr0njy@97dIg-CqtLQH#6aCY1yJDnra#=Dm>OU|-Hp&K32r_+kw1d>ka3ys zRbm>vvQJ{TA4Q$eJB!pFeFnh-3H6*FgD$0N)h32>vn<9J>*O14juiaaPSQ>rFF^J1D{l16PuoheN z(1+ZZczq)fGhBmk>e<_eOYUVheX{vdXhFE}Y|(WmqV5d0fNO0X-_8po+4A6<1~W)t zd9tP2Vm90jI$d1A=N4VqVYvLfPEF{(EYoF0aiDgno?NSw6jCb{s4y8GSZ0m2RJ<0h ztrdgiB%Bp1hKv31eYd()!NuufuIj*RTcQ_khL!V*nvVI z;6*b?^hZ>#wM9&5iQ8gw5bZIwxMmUR04)hZY6!!WBY+ckS0>+FnMH){Ae_>mdn*L# z!21d>XM%1(;HiONY{`fBiSsF_lDssu*bz32TlOXbc#tNa7H;la7?*Uwf2aFN13{Db zRQYv3#h3yNBOqE}rQfW;s|GKK311acrM})n_^IC!wfmmm8MXT!w?j|&l#f&w9==U2 zksG%3V_ql1J5%HnT-AVsr9ymx0z)DE3p&I-*0^vXos&z#Az@yuB6y3|7On`3j5rJIOpClLc0U`+g-+${A5XKNgv##|pt}}_GnTQU#xZaMaP}hY z5BP)RUKKsiojd&A&^5t2kp!Uk9Z=%w>v1 z>q9+^qh8QOzZgKA^jPJH*S|KfxehMC2gM-jDq&h1{ZQeR>SbGZ}sclaZSEF6Xk z+O;LM+Z)2V!yU@5^}EeI-Y{I3qyC<79zw|TGve>u3CFk?#h$g-AM^KyV_#LV(rI$O zwq~k0>hJRRz`66nrkvfBgP$>U!yU}dpVr)=?rv5NV{`H(*b)bLaJ0f=E>98mvjl#Q z0GTJ@&O!{#;dt`L344#gBLdu(_MOmBIRCXz`24%&`i<(4`X0c;qC;`^2)ISa?IZjS z;4xVy+T6$bMtrH2fD9KrC@u;=F4xHBfKu2nSZp#_5Cx*Si#ToW2Vr>t0AYT0vV41t`E zi$Hu>Of?Ut3sRF-ycSxYbZP^kfpXp6&PjBJcWImXw=~j#z)=GKjXKho@+C`ukT3Bf z!g29GsHn|z@=qwq%=UZ0evg>#A`k`9U;ALzbZdKX#_2WFxGvm%3tPFk7E`guP;OL@a_%7#B&-LMifCc~ zigXp7HlX2oXF_gRi%uV8m@zHR-N*s0PhfkaPqd-^0>6LZ9fE$d0s|~OyTH7OuwkVL zOC5SO<03Yt;FZ+7r?FspRdo3Id9O$yAhQd(8!fx2E&nuv_ET8y(Rg6GMHm9wLX7ojat^->Fram&7)}FLD}+}o z(;FsJB;vn#Uj9$49gDX7&xHMVfJ~pH@}Li9PkzB|nLQB5P{v8JWsv~MrZs@hY>?+4 zRp6&9h;VY{=TzL5N%%fQ)8zQxL4D=LmIR*H4!hSji8DVJls(bW7 z$U`UnNI*$E>C}f}W138&$EQwwv^B95{4$U{PJ4JvJTYIvmGulb3#K}(Bz;wf46ucS zU|a;;(iJ+62&$iiT|#&acq+v2EIazIuK?8L#P9=7_gW}I0wOnvV>4fvROEm zD8f1DA;ecBt_yc%-s5;PMeG56lKg8JdAMxkwBsW*AARKONaA#6+SwE}ff{sGb81>MkL+g740S1sQ~ICjWoFVSL*p_QT<{Fi{zp~>$i$BdfN zPOB1u1`36mlOb>0ICgeWRGlo;DEP6V9baxFjvt3*3_n|=KMYP!he`f99-UJcol*@o zbB);XUJT+>T+|1T75in z$%X>zS;7+{x3^jkR_zTJgFUXG(~=sB(u_E!kJPrxDYEnjAlO5h=Zj|luOf!`L) z^K6zWi>(nhjlM?7e@Qhk+2s|h|I&r)D_eg_y2;(d($S0B__1O53yO0oq+FUx{0IVS ziT>b_1P&J=R92^!&_h`P)Zmo-QtzdN!GFm_9Peb}E6UYWS7IR1le&<2A<>J#a5y&M dlw>Ms4~i+k8TB$CNg7FA|5h9+a)YMS{{tpCL)8EP diff --git a/dynamic/palm_dynamic_defaults.py b/dynamic/palm_dynamic_defaults.py index 85e267c..456f7c1 100644 --- a/dynamic/palm_dynamic_defaults.py +++ b/dynamic/palm_dynamic_defaults.py @@ -49,6 +49,3 @@ radiation_from_wrf = False wrf_rad_file_mask = "auxhist6_*" # smoothing distance for radiation radiation_smoothing_distance = 10000.0 - -# PALM species names - version for phstatp2 mechanism -species_names = ['NO', 'NO2', 'O3', 'PM10', 'PM25'] diff --git a/dynamic/palm_dynamic_output.py b/dynamic/palm_dynamic_output.py index 026929d..443b260 100644 --- a/dynamic/palm_dynamic_output.py +++ b/dynamic/palm_dynamic_output.py @@ -16,518 +16,314 @@ import time import numpy as np import netCDF4 from palm_wrf_utils import palm_wrf_gw +from palm_dynamic_config import * def palm_dynamic_output(wrf_files, interp_files, dynamic_driver_file, times_sec, dimensions, z_levels, z_levels_stag, ztop, z_soil_levels, dx, dy, lon_center, lat_center, rad_times_proc, rad_values_proc, sma, nested_domain): - print('Processing interpolated files to dynamic driver') + print('\nProcessing interpolated files to dynamic driver') # dimension of the time coordinate dimtimes = len(times_sec) # other coordinates dimnames = ['z', 'zw', 'zsoil', 'x','xu', 'y', 'yv'] # z height agl in m, zw staggered, zsoil 4 lev from wrf dimsize_names = ['zdim' , 'zwdim', 'zsoildim', 'xdim', 'xudim', 'ydim', 'yvdim'] # palm: zw = z - 1 - x = np.arange(dx, dimensions['xdim']*dx+dx, dx) # clean this - print('dimension of x:', len(x)) - y = np.arange(dy, dimensions['ydim']*dy+dy, dy) # clean this - print('dimension of y:', len(y)) + x = np.arange(dx, dimensions['xdim']*dx+dx, dx) + y = np.arange(dy, dimensions['ydim']*dy+dy, dy) # fill values fillvalue_float = float(-9999.0) - # check out file and remove for sure try: os.remove(dynamic_driver_file) except: pass - # create netcdf out file - print('Driver file:', dynamic_driver_file) + # boundaries + boundary = ['left','right','south', 'north','top'] + # atmospheric variables + atmos_var = ['pt','qv','u','v','w','soil_m','soil_t'] + all_variables = atmos_var + wrfchem_spec + + # prepare influx/outflux area sizes + zstag_all = np.r_[0., z_levels_stag, ztop] + zwidths = zstag_all[1:] - zstag_all[:-1] + areas_xb = np.zeros((len(z_levels), 1)) + areas_xb[:,0] = zwidths * dy + areas_yb = np.zeros((len(z_levels), 1)) + areas_yb[:,0] = zwidths * dx + areas_zb = dx*dy + area_boundaries = (areas_xb.sum()*dimensions['ydim']*2 + + areas_yb.sum()*dimensions['xdim']*2 + + areas_zb*dimensions['xdim']*dimensions['ydim']) + + # create netcdf output file + print('Creating Dynamic Driver:', dynamic_driver_file) outfile = netCDF4.Dataset(dynamic_driver_file, "w", format="NETCDF4" ) - try: - # time dimension and variable - outfile.createDimension('time', dimtimes) - _val_times = outfile.createVariable('time',"f4", ("time")) - # other dimensions and corresponding variables - for _dim in zip(dimnames, dimsize_names): #range (len(dimnames)): - print(_dim[0],_dim[1], dimensions[_dim[1]]) - outfile.createDimension(_dim[0], dimensions[_dim[1]]) - _val_z_levels = outfile.createVariable('z',"f4", ("z")) - _val_z_levels_stag = outfile.createVariable('zw',"f4", ("zw")) - _val_z_soil_levels = outfile.createVariable('zsoil',"f4", ("zsoil")) - _val_y = outfile.createVariable('y',"f4", ("y")) - _val_x = outfile.createVariable('x',"f4", ("x")) - - # prepare influx/outflux area sizes - zstag_all = np.r_[0., z_levels_stag, ztop] - zwidths = zstag_all[1:] - zstag_all[:-1] - print('zwidths', zwidths) - areas_xb = np.zeros((len(z_levels), 1)) - areas_xb[:,0] = zwidths * dy - areas_yb = np.zeros((len(z_levels), 1)) - areas_yb[:,0] = zwidths * dx - areas_zb = dx*dy - area_boundaries = (areas_xb.sum()*dimensions['ydim']*2 - + areas_yb.sum()*dimensions['xdim']*2 - + areas_zb*dimensions['xdim']*dimensions['ydim']) - - # write values for coordinates - _val_times[:] = times_sec[:] - _val_z_levels[:] = z_levels[:] - _val_z_levels_stag[:] = z_levels_stag[:] - _val_z_soil_levels[:] = z_soil_levels[:] - _val_y[:] = y[:] - _val_x[:] = x[:] - - # initialization of the variables and setting of init_* variables - print("Processing initialization from file", interp_files[0]) - # open corresponding - infile = netCDF4.Dataset(interp_files[0], "r", format="NETCDF4") - try: - # open variables in the input file - init_atmosphere_pt = infile.variables['init_atmosphere_pt'] - init_atmosphere_qv = infile.variables['init_atmosphere_qv'] - init_atmosphere_u = infile.variables['init_atmosphere_u'] - init_atmosphere_v = infile.variables['init_atmosphere_v'] - init_atmosphere_w = infile.variables['init_atmosphere_w'] - init_soil_m = infile.variables['init_soil_m'] - init_soil_t = infile.variables['init_soil_t'] - # chemistry - init_atmosphere_no = infile.variables['init_atmosphere_no'] - init_atmosphere_no2 = infile.variables['init_atmosphere_no2'] - init_atmosphere_o3 = infile.variables['init_atmosphere_o3'] - init_atmosphere_pm10 = infile.variables['init_atmosphere_pm10'] - init_atmosphere_pm2_5 = infile.variables['init_atmosphere_pm2_5'] - - # create netcdf structure - _val_init_atmosphere_pt = outfile.createVariable('init_atmosphere_pt', "f4", ("z", "y", "x"), - fill_value=fillvalue_float) - _val_init_atmosphere_pt.setncattr('lod', 2) - _val_init_atmosphere_qv = outfile.createVariable('init_atmosphere_qv', "f4", ("z", "y", "x"), - fill_value=fillvalue_float) - _val_init_atmosphere_qv.setncattr('lod', 2) - _val_init_atmosphere_u = outfile.createVariable('init_atmosphere_u', "f4", ("z", "y", "xu"), - fill_value=fillvalue_float) - _val_init_atmosphere_u.setncattr('lod', 2) - _val_init_atmosphere_v = outfile.createVariable('init_atmosphere_v', "f4", ("z", "yv", "x"), - fill_value=fillvalue_float) - _val_init_atmosphere_v.setncattr('lod', 2) - _val_init_atmosphere_w = outfile.createVariable('init_atmosphere_w', "f4", ("zw", "y", "x"), - fill_value=fillvalue_float) - _val_init_atmosphere_w.setncattr('lod', 2) - _val_init_soil_t = outfile.createVariable('init_soil_t', "f4", ("zsoil", "y", "x"), - fill_value=fillvalue_float) - _val_init_soil_t.setncattr('lod', 2) - _val_init_soil_m = outfile.createVariable('init_soil_m', "f4", ("zsoil", "y", "x"), - fill_value=fillvalue_float) - _val_init_soil_m.setncattr('lod', 2) - - # chemistry - _val_init_atmosphere_no = outfile.createVariable('init_atmosphere_no', 'f4', ('z',), - fill_value=fillvalue_float) - _val_init_atmosphere_no.setncattr('lod', 1) - _val_init_atmosphere_no.setncattr('units','ppm') - - _val_init_atmosphere_no2 = outfile.createVariable('init_atmosphere_no2', 'f4', ('z',), - fill_value=fillvalue_float) - _val_init_atmosphere_no2.setncattr('lod', 1) - _val_init_atmosphere_no2.setncattr('units','ppm') - - _val_init_atmosphere_o3 = outfile.createVariable('init_atmosphere_o3', 'f4', ('z',), - fill_value=fillvalue_float) - _val_init_atmosphere_o3.setncattr('lod', 1) - _val_init_atmosphere_o3.setncattr('units','ppm') - - _val_init_atmosphere_pm10 = outfile.createVariable('init_atmosphere_pm10', 'f4', ('z',), - fill_value=fillvalue_float) - _val_init_atmosphere_pm10.setncattr('lod', 1) - _val_init_atmosphere_pm10.setncattr('units','kg/m3') - - _val_init_atmosphere_pm2_5 = outfile.createVariable('init_atmosphere_pm2_5', 'f4', ('z',), - fill_value=fillvalue_float) - _val_init_atmosphere_pm2_5.setncattr('lod', 1) - _val_init_atmosphere_pm2_5.setncattr('units','kg/m3') - - - # time dependent variables + # create dimensions (time and coordinates) + outfile.createDimension('time', dimtimes) + for _dim in zip(dimnames, dimsize_names): + outfile.createDimension(_dim[0], dimensions[_dim[1]]) + # create variables (time and coordinates) + _val_times = outfile.createVariable('time',"f4", ("time")) + _val_times[:] = times_sec[:] + _val_z_levels = outfile.createVariable('z',"f4", ("z")) + _val_z_levels[:] = z_levels[:] + _val_z_levels_stag = outfile.createVariable('zw',"f4", ("zw")) + _val_z_levels_stag[:] = z_levels_stag[:] + _val_z_soil_levels = outfile.createVariable('zsoil',"f4", ("zsoil")) + _val_z_soil_levels[:] = z_soil_levels[:] + _val_y = outfile.createVariable('y',"f4", ("y")) + _val_y[:] = y[:] + _val_x = outfile.createVariable('x',"f4", ("x")) + _val_x[:] = x[:] + + # create all other variables in outout file + def add_interpDim(all_variables): + print('Adding dimensions to Dynamic Driver') + # surface pressure + _val_surface_forcing_surface_pressure = outfile.createVariable('surface_forcing_surface_pressure', "f4",("time")) + # geostrophic wind + _val_ls_forcing_ug = outfile.createVariable('ls_forcing_ug', "f4", ("time", "z"),fill_value=fillvalue_float) + _val_ls_forcing_vg = outfile.createVariable('ls_forcing_vg', "f4", ("time", "z"),fill_value=fillvalue_float) + # create variables in outfile for all variables + for var in all_variables: + if (var == 'pt' or var == 'qv'): + _val_init_var = outfile.createVariable('init_atmosphere_'+ var, "f4", ("z", "y", "x"),fill_value=fillvalue_float) + elif (var == 'soil_m' or var == 'soil_t'): + _val_init_var = outfile.createVariable('init_'+ var,"f4", ("zsoil", "y", "x"),fill_value=fillvalue_float) + elif var == 'u': + _val_init_var = outfile.createVariable('init_atmosphere_'+ var, "f4", ("z", "y", "xu"),fill_value=fillvalue_float) + elif var == 'v': + _val_init_var = outfile.createVariable('init_atmosphere_'+ var, "f4", ("z", "yv", "x"),fill_value=fillvalue_float) + elif var == 'w': + _val_init_var = outfile.createVariable('init_atmosphere_'+ var, "f4", ("zw", "y", "x"),fill_value=fillvalue_float) + else: + _val_init_var = outfile.createVariable('init_atmosphere_'+ var, "f4",('z',),fill_value=fillvalue_float) + # add attributes + if var in atmos_var: + _val_init_var.setncattr('lod', 2) + else: + _val_init_var.setncattr('lod', 1) + if (var == 'PM10' or var == 'PM2_5_DRY'): + _val_init_var.setncattr('units','kg/m3') + else: + _val_init_var.setncattr('units','ppm') + # create time dependent variables if not nested_domain: - # SURFACE PRESSURE - _val_surface_forcing_surface_pressure = outfile.createVariable('surface_forcing_surface_pressure', "f4", - ("time")) - # BOUNDARY - vertical slices from left, right, south, north, top - varname = 'pt' - _val_ls_forcing_pt_left = outfile.createVariable('ls_forcing_left_' + varname, "f4", ("time", "z", "y"), - fill_value=fillvalue_float) - _val_ls_forcing_pt_left.setncattr('lod', 2) - _val_ls_forcing_pt_right = outfile.createVariable('ls_forcing_right_' + varname, "f4", ("time", "z", "y"), - fill_value=fillvalue_float) - _val_ls_forcing_pt_right.setncattr('lod', 2) - _val_ls_forcing_pt_south = outfile.createVariable('ls_forcing_south_' + varname, "f4", ("time", "z", "x"), - fill_value=fillvalue_float) - _val_ls_forcing_pt_south.setncattr('lod', 2) - _val_ls_forcing_pt_north = outfile.createVariable('ls_forcing_north_' + varname, "f4", ("time", "z", "x"), - fill_value=fillvalue_float) - _val_ls_forcing_pt_north.setncattr('lod', 2) - _val_ls_forcing_pt_top = outfile.createVariable('ls_forcing_top_' + varname, "f4", ("time", "y", "x"), - fill_value=fillvalue_float) - _val_ls_forcing_pt_top.setncattr('lod', 2) - - varname = 'qv' - _val_ls_forcing_qv_left = outfile.createVariable('ls_forcing_left_' + varname, "f4", ("time", "z", "y"), - fill_value=fillvalue_float) - _val_ls_forcing_qv_left.setncattr('lod', 2) - _val_ls_forcing_qv_right = outfile.createVariable('ls_forcing_right_' + varname, "f4", ("time", "z", "y"), - fill_value=fillvalue_float) - _val_ls_forcing_qv_right.setncattr('lod', 2) - _val_ls_forcing_qv_south = outfile.createVariable('ls_forcing_south_' + varname, "f4", ("time", "z", "x"), - fill_value=fillvalue_float) - _val_ls_forcing_qv_south.setncattr('lod', 2) - _val_ls_forcing_qv_north = outfile.createVariable('ls_forcing_north_' + varname, "f4", ("time", "z", "x"), - fill_value=fillvalue_float) - _val_ls_forcing_qv_north.setncattr('lod', 2) - _val_ls_forcing_qv_top = outfile.createVariable('ls_forcing_top_' + varname, "f4", ("time", "y", "x"), - fill_value=fillvalue_float) - _val_ls_forcing_qv_top.setncattr('lod', 2) - - varname = 'u' - _val_ls_forcing_u_left = outfile.createVariable('ls_forcing_left_' + varname, "f4", ("time", "z", "y"), - fill_value=fillvalue_float) - _val_ls_forcing_u_left.setncattr('lod', 2) - _val_ls_forcing_u_right = outfile.createVariable('ls_forcing_right_' + varname, "f4", ("time", "z", "y"), - fill_value=fillvalue_float) - _val_ls_forcing_u_right.setncattr('lod', 2) - _val_ls_forcing_u_south = outfile.createVariable('ls_forcing_south_' + varname, "f4", ("time", "z", "xu"), - fill_value=fillvalue_float) - _val_ls_forcing_u_south.setncattr('lod', 2) - _val_ls_forcing_u_north = outfile.createVariable('ls_forcing_north_' + varname, "f4", ("time", "z", "xu"), - fill_value=fillvalue_float) - _val_ls_forcing_u_north.setncattr('lod', 2) - _val_ls_forcing_u_top = outfile.createVariable('ls_forcing_top_' + varname, "f4", ("time", "y", "xu"), - fill_value=fillvalue_float) - _val_ls_forcing_u_top.setncattr('lod', 2) - - varname = 'v' - _val_ls_forcing_v_left = outfile.createVariable('ls_forcing_left_' + varname, "f4", ("time", "z", "yv"), - fill_value=fillvalue_float) - _val_ls_forcing_v_left.setncattr('lod', 2) - _val_ls_forcing_v_right = outfile.createVariable('ls_forcing_right_' + varname, "f4", ("time", "z", "yv"), - fill_value=fillvalue_float) - _val_ls_forcing_v_right.setncattr('lod', 2) - _val_ls_forcing_v_south = outfile.createVariable('ls_forcing_south_' + varname, "f4", ("time", "z", "x"), - fill_value=fillvalue_float) - _val_ls_forcing_v_south.setncattr('lod', 2) - _val_ls_forcing_v_north = outfile.createVariable('ls_forcing_north_' + varname, "f4", ("time", "z", "x"), - fill_value=fillvalue_float) - _val_ls_forcing_v_north.setncattr('lod', 2) - _val_ls_forcing_v_top = outfile.createVariable('ls_forcing_top_' + varname, "f4", ("time", "yv", "x"), - fill_value=fillvalue_float) - _val_ls_forcing_v_top.setncattr('lod', 2) - - varname = 'w' - _val_ls_forcing_w_left = outfile.createVariable('ls_forcing_left_' + varname, "f4", ("time", "zw", "y"), - fill_value=fillvalue_float) - _val_ls_forcing_w_left.setncattr('lod', 2) - _val_ls_forcing_w_right = outfile.createVariable('ls_forcing_right_' + varname, "f4", ("time", "zw", "y"), - fill_value=fillvalue_float) - _val_ls_forcing_w_right.setncattr('lod', 2) - _val_ls_forcing_w_south = outfile.createVariable('ls_forcing_south_' + varname, "f4", ("time", "zw", "x"), - fill_value=fillvalue_float) - _val_ls_forcing_w_south.setncattr('lod', 2) - _val_ls_forcing_w_north = outfile.createVariable('ls_forcing_north_' + varname, "f4", ("time", "zw", "x"), - fill_value=fillvalue_float) - _val_ls_forcing_w_north.setncattr('lod', 2) - _val_ls_forcing_w_top = outfile.createVariable('ls_forcing_top_' + varname, "f4", ("time", "y", "x"), - fill_value=fillvalue_float) - _val_ls_forcing_w_top.setncattr('lod', 2) - - # chemistry - varname = 'no' - _val_ls_forcing_no_left = outfile.createVariable('ls_forcing_left_'+varname,'f4', ('time','z','y'), - fill_value=fillvalue_float) - _val_ls_forcing_no_left.setncattr('lod', 2) - _val_ls_forcing_no_left.setncattr('units','ppm') - _val_ls_forcing_no_right = outfile.createVariable('ls_forcing_right_'+varname,'f4', ('time','z','y'), - fill_value=fillvalue_float) - _val_ls_forcing_no_right.setncattr('lod', 2) - _val_ls_forcing_no_right.setncattr('units','ppm') - _val_ls_forcing_no_south = outfile.createVariable('ls_forcing_south_'+varname,'f4', ('time','z','x'), - fill_value=fillvalue_float) - _val_ls_forcing_no_south.setncattr('lod', 2) - _val_ls_forcing_no_south.setncattr('units','ppm') - _val_ls_forcing_no_north = outfile.createVariable('ls_forcing_north_'+varname,'f4', ('time','z','x'), - fill_value=fillvalue_float) - _val_ls_forcing_no_north.setncattr('lod', 2) - _val_ls_forcing_no_north.setncattr('units','ppm') - _val_ls_forcing_no_top = outfile.createVariable('ls_forcing_top_'+varname,'f4', ('time','y','x'), - fill_value=fillvalue_float) - _val_ls_forcing_no_top.setncattr('lod', 2) - _val_ls_forcing_no_top.setncattr('units','ppm') - - varname = 'no2' - _val_ls_forcing_no2_left = outfile.createVariable('ls_forcing_left_'+varname,'f4', ('time','z','y'), - fill_value=fillvalue_float) - _val_ls_forcing_no2_left.setncattr('lod', 2) - _val_ls_forcing_no2_left.setncattr('units','ppm') - _val_ls_forcing_no2_right = outfile.createVariable('ls_forcing_right_'+varname,'f4', ('time','z','y'), - fill_value=fillvalue_float) - _val_ls_forcing_no2_right.setncattr('lod', 2) - _val_ls_forcing_no2_right.setncattr('units','ppm') - _val_ls_forcing_no2_south = outfile.createVariable('ls_forcing_south_'+varname,'f4', ('time','z','x'), - fill_value=fillvalue_float) - _val_ls_forcing_no2_south.setncattr('lod', 2) - _val_ls_forcing_no2_south.setncattr('units','ppm') - _val_ls_forcing_no2_north = outfile.createVariable('ls_forcing_north_'+varname,'f4', ('time','z','x'), - fill_value=fillvalue_float) - _val_ls_forcing_no2_north.setncattr('lod', 2) - _val_ls_forcing_no2_north.setncattr('units','ppm') - _val_ls_forcing_no2_top = outfile.createVariable('ls_forcing_top_'+varname,'f4', ('time','y','x'), - fill_value=fillvalue_float) - _val_ls_forcing_no2_top.setncattr('lod', 2) - _val_ls_forcing_no2_top.setncattr('units','ppm') - - varname = 'o3' - _val_ls_forcing_o3_left = outfile.createVariable('ls_forcing_left_'+varname,'f4', ('time','z','y'), - fill_value=fillvalue_float) - _val_ls_forcing_o3_left.setncattr('lod', 2) - _val_ls_forcing_o3_left.setncattr('units','ppm') - _val_ls_forcing_o3_right = outfile.createVariable('ls_forcing_right_'+varname,'f4', ('time','z','y'), - fill_value=fillvalue_float) - _val_ls_forcing_o3_right.setncattr('lod', 2) - _val_ls_forcing_o3_right.setncattr('units','ppm') - _val_ls_forcing_o3_south = outfile.createVariable('ls_forcing_south_'+varname,'f4', ('time','z','x'), - fill_value=fillvalue_float) - _val_ls_forcing_o3_south.setncattr('lod', 2) - _val_ls_forcing_o3_south.setncattr('units','ppm') - _val_ls_forcing_o3_north = outfile.createVariable('ls_forcing_north_'+varname,'f4', ('time','z','x'), - fill_value=fillvalue_float) - _val_ls_forcing_o3_north.setncattr('lod', 2) - _val_ls_forcing_o3_north.setncattr('units','ppm') - _val_ls_forcing_o3_top = outfile.createVariable('ls_forcing_top_'+varname,'f4', ('time','y','x'), - fill_value=fillvalue_float) - _val_ls_forcing_o3_top.setncattr('lod', 2) - _val_ls_forcing_o3_top.setncattr('units','ppm') - - varname = 'pm10' - _val_ls_forcing_pm10_left = outfile.createVariable('ls_forcing_left_'+varname,'f4', ('time','z','y'), - fill_value=fillvalue_float) - _val_ls_forcing_pm10_left.setncattr('lod', 2) - _val_ls_forcing_pm10_left.setncattr('units','kg/m3') - _val_ls_forcing_pm10_right = outfile.createVariable('ls_forcing_right_'+varname,'f4', ('time','z','y'), - fill_value=fillvalue_float) - _val_ls_forcing_pm10_right.setncattr('lod', 2) - _val_ls_forcing_pm10_right.setncattr('units','kg/m3') - _val_ls_forcing_pm10_south = outfile.createVariable('ls_forcing_south_'+varname,'f4', ('time','z','x'), - fill_value=fillvalue_float) - _val_ls_forcing_pm10_south.setncattr('lod', 2) - _val_ls_forcing_pm10_south.setncattr('units','kg/m3') - _val_ls_forcing_pm10_north = outfile.createVariable('ls_forcing_north_'+varname,'f4', ('time','z','x'), - fill_value=fillvalue_float) - _val_ls_forcing_pm10_north.setncattr('lod', 2) - _val_ls_forcing_pm10_north.setncattr('units','kg/m3') - _val_ls_forcing_pm10_top = outfile.createVariable('ls_forcing_top_'+varname,'f4', ('time','y','x'), - fill_value=fillvalue_float) - _val_ls_forcing_pm10_top.setncattr('lod', 2) - _val_ls_forcing_pm10_top.setncattr('units','kg/m3') - - varname = 'pm2_5' - _val_ls_forcing_pm2_5_left = outfile.createVariable('ls_forcing_left_'+varname,'f4', ('time','z','y'), - fill_value=fillvalue_float) - _val_ls_forcing_pm2_5_left.setncattr('lod', 2) - _val_ls_forcing_pm2_5_left.setncattr('units','kg/m3') - _val_ls_forcing_pm2_5_right = outfile.createVariable('ls_forcing_right_'+varname,'f4', ('time','z','y'), - fill_value=fillvalue_float) - _val_ls_forcing_pm2_5_right.setncattr('lod', 2) - _val_ls_forcing_pm2_5_right.setncattr('units','kg/m3') - _val_ls_forcing_pm2_5_south = outfile.createVariable('ls_forcing_south_'+varname,'f4', ('time','z','x'), - fill_value=fillvalue_float) - _val_ls_forcing_pm2_5_south.setncattr('lod', 2) - _val_ls_forcing_pm2_5_south.setncattr('units','kg/m3') - _val_ls_forcing_pm2_5_north = outfile.createVariable('ls_forcing_north_'+varname,'f4', ('time','z','x'), - fill_value=fillvalue_float) - _val_ls_forcing_pm2_5_north.setncattr('lod', 2) - _val_ls_forcing_pm2_5_north.setncattr('units','kg/m3') - _val_ls_forcing_pm2_5_top = outfile.createVariable('ls_forcing_top_'+varname,'f4', ('time','y','x'), - fill_value=fillvalue_float) - _val_ls_forcing_pm2_5_top.setncattr('lod', 2) - _val_ls_forcing_pm2_5_top.setncattr('units','kg/m3') - - # geostrophic wind - _val_ls_forcing_ug = outfile.createVariable('ls_forcing_ug', "f4", ("time", "z"), - fill_value=fillvalue_float) - _val_ls_forcing_vg = outfile.createVariable('ls_forcing_vg', "f4", ("time", "z"), - fill_value=fillvalue_float) - - time.sleep(1) - - # write values for initialization variables - _val_init_atmosphere_pt[:, :, :] = init_atmosphere_pt[0, :, :, :] - _val_init_atmosphere_qv[:, :, :] = init_atmosphere_qv[0, :, :, :] - _val_init_atmosphere_u[:, :, :] = init_atmosphere_u[0, :, :, 1:] - _val_init_atmosphere_v[:, :, :] = init_atmosphere_v[0, :, 1:, :] - _val_init_atmosphere_w[:, :, :] = init_atmosphere_w[0, :, :, :] - _val_init_soil_t[:, :, :] = init_soil_t[0, :, :, :] - for k in range(0,_val_init_soil_m.shape[0]): - # adjust soil moisture according soil_moisture_adjust field (if exists) - _val_init_soil_m[k, :, :] = init_soil_m[0, k, :, :] * sma[:, :] - # chemistry - _val_init_atmosphere_no[:] = init_atmosphere_no[0,:,:,:].mean(axis=(1,2)) - _val_init_atmosphere_no2[:] = init_atmosphere_no2[0,:,:,:].mean(axis=(1,2)) - _val_init_atmosphere_o3[:] = init_atmosphere_o3[0,:,:,:].mean(axis=(1,2)) - _val_init_atmosphere_pm10[:] = init_atmosphere_pm10[0,:,:,:].mean(axis=(1,2)) - _val_init_atmosphere_pm2_5[:] = init_atmosphere_pm2_5[0,:,:,:].mean(axis=(1,2)) + if (var == 'soil_m' or var == 'soil_t'): + continue + else: + for side in boundary: + if (side == 'left' or side == 'right'): + if var == 'u': + _val_ls_forcing = outfile.createVariable('ls_forcing_'+ side +'_'+var, "f4", ("time", "z", "y"),fill_value=fillvalue_float) + elif var == 'v': + _val_ls_forcing = outfile.createVariable('ls_forcing_'+ side +'_'+var, "f4", ("time", "z", "yv"),fill_value=fillvalue_float) + elif var == 'w': + _val_ls_forcing = outfile.createVariable('ls_forcing_'+ side +'_'+var, "f4", ("time", "zw", "y"),fill_value=fillvalue_float) + else: + _val_ls_forcing = outfile.createVariable('ls_forcing_'+ side +'_'+var, "f4", ("time", "z", "y"),fill_value=fillvalue_float) + + elif (side == 'south' or side == 'north'): + if var == 'u': + _val_ls_forcing = outfile.createVariable('ls_forcing_'+ side +'_'+var, "f4", ("time", "z", "xu"),fill_value=fillvalue_float) + elif var == 'v': + _val_ls_forcing = outfile.createVariable('ls_forcing_'+ side +'_'+var, "f4", ("time", "z", "x"),fill_value=fillvalue_float) + elif var == 'w': + _val_ls_forcing = outfile.createVariable('ls_forcing_'+ side +'_'+var, "f4", ("time", "zw", "x"),fill_value=fillvalue_float) + else: + _val_ls_forcing = outfile.createVariable('ls_forcing_'+ side +'_'+var, "f4", ("time", "z", "x"),fill_value=fillvalue_float) + else: + if var == 'u': + _val_ls_forcing = outfile.createVariable('ls_forcing_'+ side +'_'+var, "f4", ("time", "y", "xu"),fill_value=fillvalue_float) + elif var == 'v': + _val_ls_forcing = outfile.createVariable('ls_forcing_'+ side +'_'+var, "f4", ("time", "yv", "x"),fill_value=fillvalue_float) + elif var == 'w': + _val_ls_forcing = outfile.createVariable('ls_forcing_'+ side +'_'+var, "f4", ("time", "y", "x"),fill_value=fillvalue_float) + else: + _val_ls_forcing = outfile.createVariable('ls_forcing_'+ side +'_'+var, "f4", ("time", "y", "x"),fill_value=fillvalue_float) + # atrributes + _val_ls_forcing.setncattr('lod', 2) + if var not in atmos_var: + if (var == 'PM10' or var == 'PM2_5_DRY'): + _val_ls_forcing.setncattr('units', 'kg/m3') + else: + _val_ls_forcing.setncattr('units', 'ppm') + # close output file + outfile.close() + # create all other variables in outout file + add_interpDim(all_variables) + + # read interpolated files and write values to outfile + def add_interpValues(all_variables): + print('\nAdding initializing variable values to Dynamic Driver') + infile = netCDF4.Dataset(interp_files[0], "r", format="NETCDF4") + outfile = netCDF4.Dataset(dynamic_driver_file, "r+", format="NETCDF4") + # initialization variables + for var in all_variables: + if (var == 'soil_m' or var == 'soil_t'): + init_var = infile.variables['init_'+var] + _val_init_var = outfile.variables['init_'+var] + else: + init_var = infile.variables['init_atmosphere_'+ var] + _val_init_var = outfile.variables['init_atmosphere_'+ var] + # write values + if (var == 'pt' or var == 'qv' or var == 'w' or var == 'soil_t'): + _val_init_var[:, :, :] = init_var[0, :, :, :] + elif var == 'u': + _val_init_var[:, :, :] = init_var[0, :, :, 1:] + elif var == 'v': + _val_init_var[:, :, :] = init_var[0, :, 1:, :] + elif var == 'soil_m': + for k in range(0,_val_init_var.shape[0]): + _val_init_var[k, :, :] = init_var[0, k, :, :] * sma[:, :] + if var in wrfchem_spec: + _val_init_var[:] = init_var[0,:,:,:].mean(axis=(1,2)) + infile.close() + outfile.close() - finally: - # close interpolated file - infile.close() - # + # time dependent variables (all time steps) if not nested_domain: - # cycle over all included time steps + print('\nAdding time dependent variables value to Dynamic Driver') + outfile = netCDF4.Dataset(dynamic_driver_file, "r+", format="NETCDF4") for ts in range(0, len(interp_files)): - print("Processing file",interp_files[ts]) - # open corresponding interpolated file + # geostrophic wind + print('Open wrf file: '+wrf_files[ts]) + nc_wrf = netCDF4.Dataset(wrf_files[ts], 'r') + ug, vg = palm_wrf_gw(nc_wrf, lon_center, lat_center, z_levels) + _val_ls_forcing_ug = outfile.variables['ls_forcing_ug'] + _val_ls_forcing_vg = outfile.variables['ls_forcing_vg'] + _val_ls_forcing_ug = ug + _val_ls_forcing_vg = vg + nc_wrf.close() + + print("Processing interpolated file: ",interp_files[ts]) infile = netCDF4.Dataset(interp_files[ts], "r", format="NETCDF4") - try: - # open variables in the input file - init_atmosphere_pt = infile.variables['init_atmosphere_pt'] - init_atmosphere_qv = infile.variables['init_atmosphere_qv'] - init_atmosphere_u = infile.variables['init_atmosphere_u'] - init_atmosphere_v = infile.variables['init_atmosphere_v'] - init_atmosphere_w = infile.variables['init_atmosphere_w'] - surface_forcing_surface_pressure = infile.variables['surface_forcing_surface_pressure'] - # chemistry - init_atmosphere_no = infile.variables['init_atmosphere_no'] - init_atmosphere_no2 = infile.variables['init_atmosphere_no2'] - init_atmosphere_o3 = infile.variables['init_atmosphere_o3'] - init_atmosphere_pm10 = infile.variables['init_atmosphere_pm10'] - init_atmosphere_pm2_5 = infile.variables['init_atmosphere_pm2_5'] - - # write values for time dependent values - # surface pressure - _val_surface_forcing_surface_pressure[ts] = np.average(surface_forcing_surface_pressure[:,:,:], axis = (1,2))[0] - # boundary conditions - _val_ls_forcing_pt_left[ts, :, :] = init_atmosphere_pt[0, :, :, 0] - _val_ls_forcing_pt_right[ts, :, :] = init_atmosphere_pt[0, :, :, dimensions['xdim'] - 1] - _val_ls_forcing_pt_south[ts, :, :] = init_atmosphere_pt[0, :, 0, :] - _val_ls_forcing_pt_north[ts, :, :] = init_atmosphere_pt[0, :, dimensions['ydim'] - 1, :] - _val_ls_forcing_pt_top[ts, :, :] = init_atmosphere_pt[0, dimensions['zdim'] - 1, :, :] - - _val_ls_forcing_qv_left[ts, :, :] = init_atmosphere_qv[0, :, :, 0] - _val_ls_forcing_qv_right[ts, :, :] = init_atmosphere_qv[0, :, :, dimensions['xdim'] - 1] - _val_ls_forcing_qv_south[ts, :, :] = init_atmosphere_qv[0, :, 0, :] - _val_ls_forcing_qv_north[ts, :, :] = init_atmosphere_qv[0, :, dimensions['ydim'] - 1, :] - _val_ls_forcing_qv_top[ts, :, :] = init_atmosphere_qv[0, dimensions['zdim'] - 1, :, :] - - _val_ls_forcing_no_left[ts, :, :] = init_atmosphere_no[0, :, :, 0] - _val_ls_forcing_no_right[ts, :, :] = init_atmosphere_no[0, :, :, dimensions['xdim'] - 1] - _val_ls_forcing_no_south[ts, :, :] = init_atmosphere_no[0, :, 0, :] - _val_ls_forcing_no_north[ts, :, :] = init_atmosphere_no[0, :, dimensions['ydim'] - 1, :] - _val_ls_forcing_no_top[ts, :, :] = init_atmosphere_no[0, dimensions['zdim'] - 1, :, :] - - # chemistry - _val_ls_forcing_no2_left[ts, :, :] = init_atmosphere_no2[0, :, :, 0] - _val_ls_forcing_no2_right[ts, :, :] = init_atmosphere_no2[0, :, :, dimensions['xdim'] - 1] - _val_ls_forcing_no2_south[ts, :, :] = init_atmosphere_no2[0, :, 0, :] - _val_ls_forcing_no2_north[ts, :, :] = init_atmosphere_no2[0, :, dimensions['ydim'] - 1, :] - _val_ls_forcing_no2_top[ts, :, :] = init_atmosphere_no2[0, dimensions['zdim'] - 1, :, :] - - _val_ls_forcing_o3_left[ts, :, :] = init_atmosphere_o3[0, :, :, 0] - _val_ls_forcing_o3_right[ts, :, :] = init_atmosphere_o3[0, :, :, dimensions['xdim'] - 1] - _val_ls_forcing_o3_south[ts, :, :] = init_atmosphere_o3[0, :, 0, :] - _val_ls_forcing_o3_north[ts, :, :] = init_atmosphere_o3[0, :, dimensions['ydim'] - 1, :] - _val_ls_forcing_o3_top[ts, :, :] = init_atmosphere_o3[0, dimensions['zdim'] - 1, :, :] - - _val_ls_forcing_pm10_left[ts, :, :] = init_atmosphere_pm10[0, :, :, 0] - _val_ls_forcing_pm10_right[ts, :, :] = init_atmosphere_pm10[0, :, :, dimensions['xdim'] - 1] - _val_ls_forcing_pm10_south[ts, :, :] = init_atmosphere_pm10[0, :, 0, :] - _val_ls_forcing_pm10_north[ts, :, :] = init_atmosphere_pm10[0, :, dimensions['ydim'] - 1, :] - _val_ls_forcing_pm10_top[ts, :, :] = init_atmosphere_pm10[0, dimensions['zdim'] - 1, :, :] - - _val_ls_forcing_pm2_5_left[ts, :, :] = init_atmosphere_pm2_5[0, :, :, 0] - _val_ls_forcing_pm2_5_right[ts, :, :] = init_atmosphere_pm2_5[0, :, :, dimensions['xdim'] - 1] - _val_ls_forcing_pm2_5_south[ts, :, :] = init_atmosphere_pm2_5[0, :, 0, :] - _val_ls_forcing_pm2_5_north[ts, :, :] = init_atmosphere_pm2_5[0, :, dimensions['ydim'] - 1, :] - _val_ls_forcing_pm2_5_top[ts, :, :] = init_atmosphere_pm2_5[0, dimensions['zdim'] - 1, :, :] - - # Perform mass balancing - uxleft = init_atmosphere_u[0, :, :, 0] - uxright = init_atmosphere_u[0, :, :, dimensions['xdim'] - 1] - vysouth = init_atmosphere_v[0, :, 0, :] - vynorth = init_atmosphere_v[0, :, dimensions['ydim'] - 1, :] - wztop = init_atmosphere_w[0, dimensions['zwdim'] - 1, :, :] - mass_disbalance = ((uxleft * areas_xb).sum() + # surface pressure + surface_forcing_surface_pressure = infile.variables['surface_forcing_surface_pressure'] + _val_surface_forcing_surface_pressure = outfile.variables['surface_forcing_surface_pressure'] + _val_surface_forcing_surface_pressure[ts] = np.average(surface_forcing_surface_pressure[:,:,:], axis = (1,2))[0] + # all other variables except soil_t, soli_m, u, v, w + #['pt', 'qv', 'u', 'v', 'w', 'soil_m', 'soil_t', 'no', 'no2', 'o3', 'PM10', 'PM2_5_DRY'] + for var in all_variables: + if (var == 'soil_m' or var == 'soil_t' or var == 'u' or var == 'v' or var == 'w'): + continue + else: + init_var= infile.variables['init_atmosphere_'+ var] + for side in boundary: + _val_ls_forcing_var = outfile.variables['ls_forcing_'+ side +'_'+var] + if side == 'left': + _val_ls_forcing_var[ts, :, :] = init_var[0, :, :, 0] + elif side == 'right': + _val_ls_forcing_var[ts, :, :] = init_var[0, :, :, dimensions['xdim'] - 1] + elif side == 'south': + _val_ls_forcing_var[ts, :, :] = init_var[0, :, 0, :] + elif side == 'north': + _val_ls_forcing_var[ts, :, :] = init_var[0, :, dimensions['ydim'] - 1, :] + elif side == 'top': + _val_ls_forcing_var[ts, :, :] = init_var[0, dimensions['zdim'] - 1, :, :] + + # variables:u, v, w - mass balancing + init_atmosphere_u = infile.variables['init_atmosphere_u'] + init_atmosphere_v = infile.variables['init_atmosphere_v'] + init_atmosphere_w = infile.variables['init_atmosphere_w'] + + uxleft = init_atmosphere_u[0, :, :, 0] + uxright = init_atmosphere_u[0, :, :, dimensions['xdim'] - 1] + vysouth = init_atmosphere_v[0, :, 0, :] + vynorth = init_atmosphere_v[0, :, dimensions['ydim'] - 1, :] + wztop = init_atmosphere_w[0, dimensions['zwdim'] - 1, :, :] + mass_disbalance = ((uxleft * areas_xb).sum() - (uxright * areas_xb).sum() + (vysouth * areas_yb).sum() - (vynorth * areas_yb).sum() - (wztop * areas_zb).sum()) - mass_corr_v = mass_disbalance / area_boundaries - print('Mass disbalance: {0:8g} m3/s (avg = {1:8g} m/s)'.format( - mass_disbalance, mass_corr_v)) - uxleft -= mass_corr_v - uxright += mass_corr_v - vysouth -= mass_corr_v - vynorth += mass_corr_v - wztop += mass_corr_v - - # Verify mass balance (optional) - #mass_disbalance = ((uxleft * areas_xb).sum() - # - (uxright * areas_xb).sum() - # + (vysouth * areas_yb).sum() - # - (vynorth * areas_yb).sum() - # - (wztop * areas_zb).sum()) - #mass_corr_v = mass_disbalance / area_boundaries - #print('Mass balanced: {0:8g} m3/s (avg = {1:8g} m/s)'.format( - # mass_disbalance, mass_corr_v)) - - _val_ls_forcing_u_left[ts, :, :] = uxleft - _val_ls_forcing_u_right[ts, :, :] = uxright - _val_ls_forcing_u_south[ts, :, :] = init_atmosphere_u[0, :, 0, 1:] - _val_ls_forcing_u_north[ts, :, :] = init_atmosphere_u[0, :, dimensions['ydim'] - 1, 1:] - _val_ls_forcing_u_top[ts, :, :] = init_atmosphere_u[0, dimensions['zdim'] - 1, :, 1:] - - _val_ls_forcing_v_left[ts, :, :] = init_atmosphere_v[0, :, 1:, 0] - _val_ls_forcing_v_right[ts, :, :] = init_atmosphere_v[0, :, 1:, dimensions['xdim'] - 1] - _val_ls_forcing_v_south[ts, :, :] = vysouth - _val_ls_forcing_v_north[ts, :, :] = vynorth - _val_ls_forcing_v_top[ts, :, :] = init_atmosphere_v[0, dimensions['zdim'] - 1, 1:, :] - - _val_ls_forcing_w_left[ts, :, :] = init_atmosphere_w[0, :, :, 0] - _val_ls_forcing_w_right[ts, :, :] = init_atmosphere_w[0, :, :, dimensions['xdim'] - 1] - _val_ls_forcing_w_south[ts, :, :] = init_atmosphere_w[0, :, 0, :] - _val_ls_forcing_w_north[ts, :, :] = init_atmosphere_w[0, :, dimensions['ydim'] - 1, :] - _val_ls_forcing_w_top[ts, :, :] = wztop - - finally: - # close interpolated file - infile.close() + mass_corr_v = mass_disbalance / area_boundaries + print('Mass disbalance: {0:8g} m3/s (avg = {1:8g} m/s)'.format(mass_disbalance, mass_corr_v)) + uxleft -= mass_corr_v + uxright += mass_corr_v + vysouth -= mass_corr_v + vynorth += mass_corr_v + wztop += mass_corr_v + + _val_ls_forcing_u_left = outfile.variables['ls_forcing_left_u'] + _val_ls_forcing_u_right = outfile.variables['ls_forcing_right_u'] + _val_ls_forcing_u_south = outfile.variables['ls_forcing_south_u'] + _val_ls_forcing_u_north = outfile.variables['ls_forcing_north_u'] + _val_ls_forcing_u_top = outfile.variables['ls_forcing_top_u'] + + _val_ls_forcing_u_left[ts, :, :] = uxleft + _val_ls_forcing_u_right[ts, :, :] = uxright + _val_ls_forcing_u_south[ts, :, :] = init_atmosphere_u[0, :, 0, 1:] + _val_ls_forcing_u_north[ts, :, :] = init_atmosphere_u[0, :, dimensions['ydim'] - 1, 1:] + _val_ls_forcing_u_top[ts, :, :] = init_atmosphere_u[0, dimensions['zdim'] - 1, :, 1:] + + _val_ls_forcing_v_left = outfile.variables['ls_forcing_left_v'] + _val_ls_forcing_v_right = outfile.variables['ls_forcing_right_v'] + _val_ls_forcing_v_south = outfile.variables['ls_forcing_south_v'] + _val_ls_forcing_v_north = outfile.variables['ls_forcing_north_v'] + _val_ls_forcing_v_top = outfile.variables['ls_forcing_top_v'] + + _val_ls_forcing_v_left[ts, :, :] = init_atmosphere_v[0, :, 1:, 0] + _val_ls_forcing_v_right[ts, :, :] = init_atmosphere_v[0, :, 1:, dimensions['xdim'] - 1] + _val_ls_forcing_v_south[ts, :, :] = vysouth + _val_ls_forcing_v_north[ts, :, :] = vynorth + _val_ls_forcing_v_top[ts, :, :] = init_atmosphere_v[0, dimensions['zdim'] - 1, 1:, :] + + _val_ls_forcing_w_left = outfile.variables['ls_forcing_left_w'] + _val_ls_forcing_w_right = outfile.variables['ls_forcing_right_w'] + _val_ls_forcing_w_south = outfile.variables['ls_forcing_south_w'] + _val_ls_forcing_w_north = outfile.variables['ls_forcing_north_w'] + _val_ls_forcing_w_top = outfile.variables['ls_forcing_top_w'] + + _val_ls_forcing_w_left[ts, :, :] = init_atmosphere_w[0, :, :, 0] + _val_ls_forcing_w_right[ts, :, :] = init_atmosphere_w[0, :, :, dimensions['xdim'] - 1] + _val_ls_forcing_w_south[ts, :, :] = init_atmosphere_w[0, :, 0, :] + _val_ls_forcing_w_north[ts, :, :] = init_atmosphere_w[0, :, dimensions['ydim'] - 1, :] + _val_ls_forcing_w_top[ts, :, :] = wztop + + infile.close() + outfile.close() + + # read interpolated files and write values to outfile + add_interpValues(all_variables) + + # add radiation input + if len(rad_times_proc) > 0: + outfile = netCDF4.Dataset(dynamic_driver_file, "a", format="NETCDF4") + # radiation time dimension and variable + outfile.createDimension('time_rad', len(rad_times_proc)) + _val_times = outfile.createVariable('time_rad',"f4", ("time_rad")) + _val_times[:] = rad_times_proc[:] + # radiation variables + var = outfile.createVariable('rad_sw_in', "f4", ("time_rad"), fill_value=fillvalue_float) + var.setncattr('lod', 1) + var.units = 'W/m2' + var[:] = rad_values_proc[0][:] + var = outfile.createVariable('rad_lw_in', "f4", ("time_rad"), fill_value=fillvalue_float) + var.setncattr('lod', 1) + var.units = 'W/m2' + var[:] = rad_values_proc[1][:] + var = outfile.createVariable('rad_sw_in_dif', "f4", ("time_rad"), fill_value=fillvalue_float) + var.setncattr('lod', 1) + var.units = 'W/m2' + var[:] = rad_values_proc[2][:] - # write geostrophic wind - print('Open wrf file '+wrf_files[ts]) - nc_wrf = netCDF4.Dataset(wrf_files[ts], 'r') - try: - ug, vg = palm_wrf_gw(nc_wrf, lon_center, lat_center, z_levels) - _val_ls_forcing_ug[ts, :] = ug - _val_ls_forcing_vg[ts, :] = vg - finally: - nc_wrf.close() - - - if len(rad_times_proc) > 0: - # process radiation inputs - # radiation time dimension and variable - outfile.createDimension('time_rad', len(rad_times_proc)) - _val_times = outfile.createVariable('time_rad',"f4", ("time_rad")) - _val_times[:] = rad_times_proc[:] - # radiation variables - var = outfile.createVariable('rad_sw_in', "f4", ("time_rad"), fill_value=fillvalue_float) - var.setncattr('lod', 1) - var.units = 'W/m2' - var[:] = rad_values_proc[0][:] - var = outfile.createVariable('rad_lw_in', "f4", ("time_rad"), fill_value=fillvalue_float) - var.setncattr('lod', 1) - var.units = 'W/m2' - var[:] = rad_values_proc[1][:] - var = outfile.createVariable('rad_sw_in_dif', "f4", ("time_rad"), fill_value=fillvalue_float) - var.setncattr('lod', 1) - var.units = 'W/m2' - var[:] = rad_values_proc[2][:] - - finally: outfile.close() + + # correct chemical species names to PALM output + outfile = netCDF4.Dataset(dynamic_driver_file, "a", format="NETCDF4") + if 'PM2_5_DRY' in all_variables: + outfile.renameVariable('init_atmosphere_PM2_5_DRY','init_atmosphere_pm2_5') + outfile.renameVariable('ls_forcing_left_PM2_5_DRY', 'ls_forcing_left_pm2_5') + outfile.renameVariable('ls_forcing_right_PM2_5_DRY','ls_forcing_right_pm2_5') + outfile.renameVariable('ls_forcing_south_PM2_5_DRY','ls_forcing_south_pm2_5') + outfile.renameVariable('ls_forcing_north_PM2_5_DRY','ls_forcing_north_pm2_5') + outfile.renameVariable('ls_forcing_top_PM2_5_DRY','ls_forcing_top_pm2_5') + if 'PM10' in all_variables: + outfile.renameVariable('init_atmosphere_PM10','init_atmosphere_pm10') + outfile.renameVariable('ls_forcing_left_PM10', 'ls_forcing_left_pm10') + outfile.renameVariable('ls_forcing_right_PM10','ls_forcing_right_pm10') + outfile.renameVariable('ls_forcing_south_PM10','ls_forcing_south_pm10') + outfile.renameVariable('ls_forcing_north_PM10','ls_forcing_north_pm10') + outfile.renameVariable('ls_forcing_top_PM10','ls_forcing_top_pm10') + outfile.close() diff --git a/dynamic/palm_wrf_utils.py b/dynamic/palm_wrf_utils.py index 33f394a..96c95c7 100644 --- a/dynamic/palm_wrf_utils.py +++ b/dynamic/palm_wrf_utils.py @@ -364,44 +364,24 @@ def palm_wrf_vertical_interp(infile, outfile, wrffile, z_levels, z_levels_stag, # zsoil is taken from wrf - not need to define it # ======================== CHEMISTRY =================================================== # convert wrf-chem units to PALM units - #NO:ppmv to ppm - no_mixing_raw = nc_infile.variables['no'][0] - no_mixing_raw = np.r_[no_mixing_raw[0:1], no_mixing_raw] - no_mixing = interpolate_1d(z_levels, height, no_mixing_raw) - vdata = nc_outfile.createVariable('init_atmosphere_no',"f4", ("Time", "z","south_north","west_east")) - vdata[0,:,:,:] = no_mixing - - #NO2: ppmv to ppm - no2_mixing_raw = nc_infile.variables['no2'][0] - no2_mixing_raw = np.r_[no2_mixing_raw[0:1], no2_mixing_raw] - no2_mixing = interpolate_1d(z_levels, height, no2_mixing_raw) - vdata = nc_outfile.createVariable('init_atmosphere_no2',"f4", ("Time", "z","south_north","west_east")) - vdata[0,:,:,:] = no2_mixing - - #O3:ppmv to ppm - o3_mixing_raw = nc_infile.variables['o3'][0] - o3_mixing_raw = np.r_[o3_mixing_raw[0:1], o3_mixing_raw] - o3_mixing = interpolate_1d(z_levels, height, o3_mixing_raw) - vdata = nc_outfile.createVariable('init_atmosphere_o3',"f4", ("Time", "z","south_north","west_east")) - vdata[0,:,:,:] = o3_mixing - - #PM10: micrograms m-3 to kg/m3 - pm10_raw = nc_infile.variables['PM10'][0] - pm10_raw = np.r_[pm10_raw[0:1], pm10_raw] - pm10_mass = interpolate_1d(z_levels, height, (pm10_raw)*1e-9) - vdata = nc_outfile.createVariable('init_atmosphere_pm10',"f4", ("Time", "z","south_north","west_east")) - vdata[0,:,:,:] = pm10_mass - - #PM2_5_DRY: micrograms m-3 to kg/m3 - pm2_5_raw = nc_infile.variables['PM2_5_DRY'][0] - pm2_5_raw = np.r_[pm2_5_raw[0:1], pm2_5_raw] - pm2_5_mass = interpolate_1d(z_levels, height, (pm2_5_raw)*1e-9) - vdata = nc_outfile.createVariable('init_atmosphere_pm2_5',"f4", ("Time", "z","south_north","west_east")) - vdata[0,:,:,:] = pm2_5_mass - - nc_infile.close() - nc_wrf.close() - nc_outfile.close() + # convert ppmv to ppm for PALM except PM10 & PM2_5_DRY: micrograms m-3 to kg/m3 + + def chem_from_wrfchem(wrfchem_spec): + for spec in wrfchem_spec: + chem_data_raw = nc_infile.variables[spec][0] + chem_data_raw = np.r_[chem_data_raw[0:1], chem_data_raw] + if (spec == 'PM10' or spec == 'PM2_5_DRY'): + chem_data = interpolate_1d(z_levels, height, (chem_data_raw)*1e-9) + else: + chem_data = interpolate_1d(z_levels, height, chem_data_raw) + vdata = nc_outfile.createVariable('init_atmosphere_'+spec,"f4",("Time", "z","south_north","west_east")) + vdata[0,:,:,:] = chem_data + + nc_infile.close() + nc_wrf.close() + nc_outfile.close() + + chem_from_wrfchem(wrfchem_spec) def palm_wrf_gw(f, lon, lat, levels, tidx=0): '''Calculate geostrophic wind from WRF using metpy''' -- GitLab