@@ -744,6 +744,10 @@ def call_ncrcat(cmd):
744744 # End if cam
745745 # End if has_lev
746746
747+ # Optional, add additional time varying variables to clm2.h1 files
748+ if ("clm" in hist_str ) and ("h1" in hist_str ):
749+ ncrcat_var_list = ncrcat_var_list + ",pfts1d_wtgcell,pfts1d_wtlunit,pfts1d_wtcol"
750+
747751 cmd = (
748752 ["ncrcat" , "-O" , "-4" , "-h" , "--no_cll_mth" , "-v" , ncrcat_var_list ]
749753 + hist_files
@@ -764,23 +768,12 @@ def call_ncrcat(cmd):
764768 "ncatted" , "-O" ,
765769 "-a" , "adf_user,global,a,c," + f"{ self .user } " ,
766770 "-a" , "hist_file_locs,global,a,c," + f"{ hist_locs_str } " ,
767- "-a" , "hist_file_list,global,a,c," + f"{ hist_files_str } " ,
771+ # This list is too long and fails
772+ # "-a", "hist_file_list,global,a,c," + f"{hist_files_str}",
768773 ts_outfil_str
769774 ]
770775
771- if "clm" in hist_str :
772- # Step 3b: Optional, add additional variables to clm2.h1 files
773- if "h1" in hist_str :
774- cmd_add_clm_h1_fields = [
775- "ncrcat" , "-A" , "-v" ,
776- "pfts1d_ixy,pfts1d_jxy,pfts1d_itype_veg,lat,lon" ,
777- hist_files ,
778- ts_outfil_str
779- ]
780- # add time varrying information to clm2.h1 fields
781- list_of_hist_commands .append (cmd_add_clm_h1_fields )
782-
783- # Step 3c: Create the ncatted command to remove the history attribute
776+ # Step 3c: Create the ncatted command to remove the history attribute
784777 cmd_remove_history = [
785778 "ncatted" , "-O" , "-h" ,
786779 "-a" , "history,global,d,," ,
@@ -858,7 +851,19 @@ def call_ncrcat(cmd):
858851 ts_ds ['area' ] = ds .area
859852 ts_ds ['landfrac' ] = ds .landfrac
860853 ts_ds ['landmask' ] = ds .landmask
861-
854+ # Optional, add additional variables to clm2.h1 files
855+ # Note: this is currently set up for PFT output
856+ if "h1" in hist_str :
857+ ds = xr .open_dataset (hist_files [0 ], decode_times = False )
858+ ts_ds ['pfts1d_ixy' ] = ds .pfts1d_ixy
859+ ts_ds ['pfts1d_jxy' ] = ds .pfts1d_jxy
860+ ts_ds ['pfts1d_gi' ] = ds .pfts1d_gi
861+ ts_ds ['pfts1d_li' ] = ds .pfts1d_li
862+ ts_ds ['pfts1d_ci' ] = ds .pfts1d_li
863+ ts_ds ['pfts1d_itype_veg' ] = ds .pfts1d_itype_veg
864+ ts_ds ['pfts1d_itype_col' ] = ds .pfts1d_itype_col
865+ ts_ds ['pfts1d_itype_lunit' ] = ds .pfts1d_itype_lunit
866+ ts_ds ['pfts1d_active' ] = ds .pfts1d_active
862867 ts_ds ['time' ] = time
863868 ts_ds .assign_coords (time = time )
864869 ts_ds_fixed = xr .decode_cf (ts_ds )
@@ -1353,6 +1358,10 @@ def derive_variables(self, res=None, hist_str=None, vars_to_derive=None, ts_dir=
13531358 # create new file name for derived variable
13541359 derived_file = constit_files [0 ].replace (constit_list [0 ], var )
13551360
1361+ clmPDB = 0.0112372 # ratio of 13C/12C in Pee Dee Belemnite (C isotope standard)
1362+ clm14C = 1e-12 # not accepted value of 1.176 x 10-12
1363+ min13C = - 40. # prevent wacky values when 12C stock or fluxes are very small
1364+ min14C = - 400. # arbitrary
13561365 # Check if clobber is true for file
13571366 if Path (derived_file ).is_file ():
13581367 if overwrite :
@@ -1370,6 +1379,51 @@ def derive_variables(self, res=None, hist_str=None, vars_to_derive=None, ts_dir=
13701379 der_val = 100 * ds ["FSR" ]/ ds ["FSDS" ].where (ds ["FSDS" ]> 0 )
13711380 elif var == "RNET" :
13721381 der_val = ds ["FSA" ]- ds ["FIRA" ]
1382+ elif var == "WUE" :
1383+ der_val = ds ["GPP" ]/ ds ["FCTR" ].where (ds ["FCTR" ]> 0 )
1384+
1385+ # ----------------------------------------------------------------------------------
1386+ # Isotope-specific derived variables
1387+ # ----------------------------------------------------------------------------------
1388+ # NOTE: del13C : valid range = -40 to 0 per mil PDB
1389+ # formulas similar to Jain et al 1996 Tellus B 48B: 583-600
1390+ # as applied in land_diags /glade/campaign/cgd/tss/people/oleson/FROM_LMWG/diag/lnd_diag4.2/code/shared/lnd_func.ncl
1391+ # TODO, this would be nice to avoid repeating for all isotopes enables variables
1392+ # TODO, check for accuracy of equations, neither as as in Jain et al 1996...
1393+ # Should they just be ((ratio sample - ratio standard) / ratio standard) * 1000 ?
1394+
1395+ elif var == "C13_GPP_pm" :
1396+ der_val = (((ds ["C13_GPP" ]/ ds ["GPP" ].where (ds ["GPP" ]> 0 )) / clmPDB ) - 1. ) * 1e3
1397+ der_val = der_val .where (der_val > min13C )
1398+
1399+ elif var == "C14_GPP_pm" :
1400+ der_val = (((ds ["C14_GPP" ]/ ds ["GPP" ].where (ds ["GPP" ]> 0 )) / clm14C ) - 1. ) * 1e3
1401+ der_val = der_val .where (der_val > min14C )
1402+
1403+ elif var == "C13_NPP_pm" :
1404+ der_val = (((ds ["C13_NPP" ]/ ds ["NPP" ].where (ds ["NPP" ]> 0 )) / clmPDB ) - 1. ) * 1e3
1405+ der_val = der_val .where (der_val > min13C )
1406+
1407+ elif var == "C14_NPP_pm" :
1408+ der_val = (((ds ["C14_NPP" ]/ ds ["NPP" ].where (ds ["NPP" ]> 0 )) / clm14C ) - 1. ) * 1e3
1409+ der_val = der_val .where (der_val > min14C )
1410+
1411+ elif var == "C13_TOTVEGC_pm" :
1412+ der_val = (((ds ["C13_TOTVEGC" ]/ ds ["TOTVEGC" ].where (ds ["TOTVEGC" ]> 0 )) / clmPDB ) - 1. ) * 1e3
1413+ der_val = der_val .where (der_val > min13C )
1414+
1415+ elif var == "C14_TOTVEGC_pm" :
1416+ der_val = (((ds ["C14_TOTVEGC" ]/ ds ["TOTVEGC" ].where (ds ["TOTVEGC" ]> 0 )) / clm14C ) - 1. ) * 1e3
1417+ der_val = der_val .where (der_val > min14C )
1418+
1419+ elif var == "C13_TOTECOSYSC_pm" :
1420+ der_val = (((ds ["C13_TOTECOSYSC" ]/ ds ["TOTECOSYSC" ].where (ds ["TOTECOSYSC" ]> 0 )) / clmPDB ) - 1. ) * 1e3
1421+ der_val = der_val .where (der_val > min13C )
1422+
1423+ elif var == "C14_TOTECOSYSC_pm" :
1424+ der_val = (((ds ["C14_TOTECOSYSC" ]/ ds ["TOTECOSYSC" ].where (ds ["TOTECOSYSC" ]> 0 )) / clm14C ) - 1. ) * 1e3
1425+ #der_val = der_val.where(der_val > min14C)
1426+
13731427 else :
13741428 # Loop through all constituents and sum
13751429 der_val = 0
0 commit comments