Subset test sep26

C/projects/CARDAMOM_GENERAL/CARDAMOM_NETCDF_DATA_STRUCTURE.c

@@ -126,5 +126,13 @@ TIMESERIES_DRIVER_STRUCT YIELD;
 //MCMCID
 MCMCID_STRUCT MCMCID;
 
+//Optional output subsets
+char **FLUXES_SUBSET_NAMES;
+int *FLUXES_SUBSET_INDICES;
+int FLUXES_SUBSET_COUNT;
+
+char **POOLS_SUBSET_NAMES;
+int *POOLS_SUBSET_INDICES;
+int POOLS_SUBSET_COUNT;
 
 }NETCDF_DATA;

C/projects/CARDAMOM_GENERAL/CARDAMOM_READ_NETCDF_DATA.c

@@ -389,6 +389,13 @@ double alpha = asin((sin(pi/180*DATA->LAT)*sin(pi/180*DA)+cos(pi/180*DATA->LAT)*
 printf("Done reading all data");
 
 
+DATA->FLUXES_SUBSET_NAMES = ncdf_read_string_array(ncid, "FLUXES_SUBSET", &DATA->FLUXES_SUBSET_COUNT);
+DATA->POOLS_SUBSET_NAMES = ncdf_read_string_array(ncid, "POOLS_SUBSET", &DATA->POOLS_SUBSET_COUNT);
+DATA->FLUXES_SUBSET_INDICES = NULL;
+DATA->POOLS_SUBSET_INDICES = NULL;
+
+printf("FLUXES_SUBSET_COUNT = %d\n", DATA->FLUXES_SUBSET_COUNT);
+printf("POOLS_SUBSET_COUNT = %d\n", DATA->POOLS_SUBSET_COUNT);
 
 
 MCMCID_STRUCT MCMCID;

C/projects/CARDAMOM_GENERAL/CARDAMOM_RUN_MODEL.c

@@ -28,6 +28,40 @@
 #define min(a,b) ({ __typeof__ (a) _a = (a);  __typeof__ (b) _b = (b);  _a < _b ? _a : _b; })
 
 
+int *build_subset_indices(char **abbreviations, int nabbrevs,
+                          char **subset_names, int nsubset, int *out_count) {
+    if (subset_names == NULL || nsubset == 0) {
+        *out_count = nabbrevs;
+        int *all_indices = calloc(nabbrevs, sizeof(int));
+        for (int i = 0; i < nabbrevs; i++) {
+            all_indices[i] = i;
+        }
+        return all_indices;
+    }
+
+    int *indices = calloc(nsubset, sizeof(int));
+    int found_count = 0;
+
+    for (int s = 0; s < nsubset; s++) {
+        int found = 0;
+        for (int i = 0; i < nabbrevs; i++) {
+            if (abbreviations[i] != NULL &&
+                strcmp(abbreviations[i], subset_names[s]) == 0) {
+                indices[found_count++] = i;
+                found = 1;
+                break;
+            }
+        }
+        if (!found) {
+            printf("Warning: Subset name '%s' not found in abbreviations\n", subset_names[s]);
+        }
+    }
+
+    *out_count = found_count;
+    return indices;
+}
+
+
 //This scans the string and removes all instances of the string toFind, and replaces them with the single char toReplace.
 void str_inplace_replace(char * str, const char * toFind, const char toReplace){
   //Yeah this implementation is N^2... but we have a small fixed max N, so don't @ me.
@@ -146,10 +180,31 @@ int sampleDimID,  timePoolsDimID,timeFluxesDimID, probIdxDimID,edcIdxDimID, noLi
 FAILONERROR(nc_def_dim(ncid,"Sample",N,&sampleDimID));
 
 
+struct FLUX_META_STRUCT fluxInfo_pre = ((DALEC *)CARDADATA.MODEL)->FLUX_META;
+struct POOLS_META_STRUCT poolsInfo_pre = ((DALEC *)CARDADATA.MODEL)->POOLS_META;
+
+CARDADATA.ncdf_data.FLUXES_SUBSET_INDICES = build_subset_indices(
+    fluxInfo_pre.ABBREVIATION, CARDADATA.nofluxes,
+    CARDADATA.ncdf_data.FLUXES_SUBSET_NAMES,
+    CARDADATA.ncdf_data.FLUXES_SUBSET_COUNT,
+    &CARDADATA.ncdf_data.FLUXES_SUBSET_COUNT);
+
+CARDADATA.ncdf_data.POOLS_SUBSET_INDICES = build_subset_indices(
+    poolsInfo_pre.ABBREVIATION, CARDADATA.nopools,
+    CARDADATA.ncdf_data.POOLS_SUBSET_NAMES,
+    CARDADATA.ncdf_data.POOLS_SUBSET_COUNT,
+    &CARDADATA.ncdf_data.POOLS_SUBSET_COUNT);
+
+int output_flux_count = CARDADATA.ncdf_data.FLUXES_SUBSET_COUNT;
+int output_pool_count = CARDADATA.ncdf_data.POOLS_SUBSET_COUNT;
+
+printf("Output flux count: %d (of %d total)\n", output_flux_count, CARDADATA.nofluxes);
+printf("Output pool count: %d (of %d total)\n", output_pool_count, CARDADATA.nopools);
+
 int poolDimID;
-FAILONERROR(nc_def_dim(ncid,"Pool",CARDADATA.nopools,&poolDimID ));
+FAILONERROR(nc_def_dim(ncid,"Pool",output_pool_count,&poolDimID ));
 int fluxDimID;
-FAILONERROR(nc_def_dim(ncid,"Flux",CARDADATA.nofluxes,&fluxDimID ));
+FAILONERROR(nc_def_dim(ncid,"Flux",output_flux_count,&fluxDimID ));
 int noParsDimID;
 FAILONERROR(nc_def_dim(ncid,"Parameter",CARDADATA.nopars,&noParsDimID ));
 
@@ -183,7 +238,9 @@ FAILONERROR(nc_def_var(	ncid,"FLUXES" , NC_DOUBLE, 3, fluxes_dems, &(fluxesVarID
 
 //Create each flux's mapping as an attribute
 struct FLUX_META_STRUCT fluxInfo = ((DALEC *)CARDADATA.MODEL)->FLUX_META;
-for(int i = 0; i < CARDADATA.nofluxes; i++){
+
+for(int s = 0; s < output_flux_count; s++){
+  int i = CARDADATA.ncdf_data.FLUXES_SUBSET_INDICES[s];
   const char* ncVarAbbreviation =(const char *) calloc(sizeof(char), METADATA_MAX_LEN );//WARNING: DO NOT FREE THIS ARRAY! Netcdf libs require a const char*, so whatever is inside the string should not change or be freed!
 
   if (fluxInfo.ABBREVIATION != NULL && fluxInfo.ABBREVIATION[i] != NULL){
@@ -195,7 +252,7 @@ for(int i = 0; i < CARDADATA.nofluxes; i++){
 
   }
   //FAILONERROR(nc_def_var(	ncid,ncVarAbbreviation , NC_DOUBLE, 2, fluxes_dems, &(fluxesVarID[i]) ));
-  WARNONERROR(nc_put_att_int	(	ncid,fluxesVarID,ncVarAbbreviation,NC_INT,1,&i));
+  WARNONERROR(nc_put_att_int	(	ncid,fluxesVarID,ncVarAbbreviation,NC_INT,1,&s));
 }
 //metadata vars
 FAILONERROR(nc_def_var(	ncid,"FLUX_NAMES" , NC_CHAR, 2, fluxes_meta_dems, &(fluxesNameVarID) ));
@@ -204,16 +261,18 @@ FAILONERROR(nc_def_var(	ncid,"FLUX_UNITS" , NC_CHAR, 2, fluxes_meta_dems, &(flux
 
 
 //POOLS DEFINITION
-//Create each pool variable as its own var inside 
+//Create each pool variable as its own var inside
 int poolsVarID,poolsNameVarID,poolsDescriptionVarID, poolsUnitVarID;
 struct POOLS_META_STRUCT poolsInfo = ((DALEC *)CARDADATA.MODEL)->POOLS_META;
+
 int pools_dems[] = {sampleDimID,timePoolsDimID, poolDimID}; //poolsDimId was last in the order
 int pools_meta_dems[] = {poolDimID, chidDimID};
 
 FAILONERROR(nc_def_var(	ncid,"POOLS" , NC_DOUBLE, 3, pools_dems, &(poolsVarID) ));
 
 
-for(int i = 0; i < CARDADATA.nopools; i++){
+for(int s = 0; s < output_pool_count; s++){
+  int i = CARDADATA.ncdf_data.POOLS_SUBSET_INDICES[s];
   const char* ncVarAbbreviation =(const char *) calloc(sizeof(char), METADATA_MAX_LEN );//WARNING: DO NOT FREE THIS ARRAY! Netcdf libs require a const char*, so whatever is inside the string should not change or be freed!
   if (poolsInfo.ABBREVIATION != NULL && poolsInfo.ABBREVIATION[i] != NULL ){
     snprintf( (char *) ncVarAbbreviation,METADATA_MAX_LEN-1,"POOL-%s", poolsInfo.ABBREVIATION[i] );//Write to it once, overriding the const qualifier so it is set
@@ -223,7 +282,7 @@ for(int i = 0; i < CARDADATA.nopools; i++){
     printf("ERROR in %s at %d: pool ID %d has no defined ABBREVIATION in it's POOLS_META. Add it to your DALEC_####_NC_INFO.c file! This pool will be called %s until you do!\n", __FILE__, __LINE__,i,ncVarAbbreviation);
 
   }
-  WARNONERROR(nc_put_att_int	(	ncid,poolsVarID,ncVarAbbreviation,NC_INT,1,&i));
+  WARNONERROR(nc_put_att_int	(	ncid,poolsVarID,ncVarAbbreviation,NC_INT,1,&s));
 
 
   /*if (poolsInfo.NAME != NULL && poolsInfo.NAME[i] != NULL){
@@ -299,34 +358,36 @@ nc_enddef(ncid);
 
 
 //Insert Fluxes metadata
-for(int i = 0; i < CARDADATA.nofluxes; i++){
+for(int s = 0; s < output_flux_count; s++){
+  int i = CARDADATA.ncdf_data.FLUXES_SUBSET_INDICES[s];
   if (fluxInfo.NAME != NULL && fluxInfo.NAME[i] != NULL){
     //"Name"
-    WARNONERROR(nc_put_vara_text	(	ncid,fluxesNameVarID,(const size_t[]){i,0},(const size_t[]){1,min(METADATA_MAX_LEN-1,strlen(fluxInfo.NAME[i]))},(const char *)fluxInfo.NAME[i]));
+    WARNONERROR(nc_put_vara_text	(	ncid,fluxesNameVarID,(const size_t[]){s,0},(const size_t[]){1,min(METADATA_MAX_LEN-1,strlen(fluxInfo.NAME[i]))},(const char *)fluxInfo.NAME[i]));
   }
   if (fluxInfo.DESCRIPTION != NULL && fluxInfo.DESCRIPTION[i] != NULL){
     //"Description"
-    WARNONERROR(nc_put_vara_text	(	ncid,fluxesDescriptionVarID,(const size_t[]){i,0},(const size_t[]){1,min(METADATA_MAX_LEN-1,strlen(fluxInfo.DESCRIPTION[i]))},(const char *)fluxInfo.DESCRIPTION[i]));
+    WARNONERROR(nc_put_vara_text	(	ncid,fluxesDescriptionVarID,(const size_t[]){s,0},(const size_t[]){1,min(METADATA_MAX_LEN-1,strlen(fluxInfo.DESCRIPTION[i]))},(const char *)fluxInfo.DESCRIPTION[i]));
   }
   if (fluxInfo.UNITS != NULL && fluxInfo.UNITS[i] != NULL){
     //"Units"
-    WARNONERROR(nc_put_vara_text	(	ncid,fluxesUnitVarID,(const size_t[]){i,0},(const size_t[]){1,min(METADATA_MAX_LEN-1,strlen(fluxInfo.UNITS[i]))},(const char *)fluxInfo.UNITS[i]));
+    WARNONERROR(nc_put_vara_text	(	ncid,fluxesUnitVarID,(const size_t[]){s,0},(const size_t[]){1,min(METADATA_MAX_LEN-1,strlen(fluxInfo.UNITS[i]))},(const char *)fluxInfo.UNITS[i]));
   }
 }
 
 //Insert Pools metadata
-for(int i = 0; i < CARDADATA.nopools; i++){
+for(int s = 0; s < output_pool_count; s++){
+  int i = CARDADATA.ncdf_data.POOLS_SUBSET_INDICES[s];
   if (poolsInfo.NAME != NULL && poolsInfo.NAME[i] != NULL){
     //"Name"
-    WARNONERROR(nc_put_vara_text	(	ncid,poolsNameVarID,(const size_t[]){i,0},(const size_t[]){1,min(METADATA_MAX_LEN-1,strlen(poolsInfo.NAME[i]))},(const char *)poolsInfo.NAME[i]));
+    WARNONERROR(nc_put_vara_text	(	ncid,poolsNameVarID,(const size_t[]){s,0},(const size_t[]){1,min(METADATA_MAX_LEN-1,strlen(poolsInfo.NAME[i]))},(const char *)poolsInfo.NAME[i]));
   }
   if (poolsInfo.DESCRIPTION != NULL && poolsInfo.DESCRIPTION[i] != NULL){
     //"Description"
-    WARNONERROR(nc_put_vara_text	(	ncid,poolsDescriptionVarID,(const size_t[]){i,0},(const size_t[]){1,min(METADATA_MAX_LEN-1,strlen(poolsInfo.DESCRIPTION[i]))},(const char *)poolsInfo.DESCRIPTION[i]));
+    WARNONERROR(nc_put_vara_text	(	ncid,poolsDescriptionVarID,(const size_t[]){s,0},(const size_t[]){1,min(METADATA_MAX_LEN-1,strlen(poolsInfo.DESCRIPTION[i]))},(const char *)poolsInfo.DESCRIPTION[i]));
   }
   if (poolsInfo.UNITS != NULL && poolsInfo.UNITS[i] != NULL){
     //"Units"
-    WARNONERROR(nc_put_vara_text	(	ncid,poolsUnitVarID,(const size_t[]){i,0},(const size_t[]){1,min(METADATA_MAX_LEN-1,strlen(poolsInfo.UNITS[i]))},(const char *)poolsInfo.UNITS[i]));
+    WARNONERROR(nc_put_vara_text	(	ncid,poolsUnitVarID,(const size_t[]){s,0},(const size_t[]){1,min(METADATA_MAX_LEN-1,strlen(poolsInfo.UNITS[i]))},(const char *)poolsInfo.UNITS[i]));
   }
 }
 
@@ -421,8 +482,25 @@ clock_t    end = clock();//End timer
 //(with N (Number of samples) being another dimension, applied to all vars)
 
 
-FAILONERROR(nc_put_vara_double(ncid,fluxesVarID,(const size_t []){n,0,0}, (const size_t[]){1,Ntimesteps,CARDADATA.nofluxes}, CARDADATA.M_FLUXES));
-FAILONERROR(nc_put_vara_double(ncid,poolsVarID,(const size_t []){n,0,0}, (const size_t[]){1,Ntimesteps+1,CARDADATA.nopools}, CARDADATA.M_POOLS));
+double *flux_subset = calloc(Ntimesteps * output_flux_count, sizeof(double));
+for(int t = 0; t < Ntimesteps; t++){
+  for(int s = 0; s < output_flux_count; s++){
+    int flux_idx = CARDADATA.ncdf_data.FLUXES_SUBSET_INDICES[s];
+    flux_subset[t * output_flux_count + s] = CARDADATA.M_FLUXES[t * CARDADATA.nofluxes + flux_idx];
+  }
+}
+FAILONERROR(nc_put_vara_double(ncid,fluxesVarID,(const size_t []){n,0,0}, (const size_t[]){1,Ntimesteps,output_flux_count}, flux_subset));
+free(flux_subset);
+
+double *pool_subset = calloc((Ntimesteps+1) * output_pool_count, sizeof(double));
+for(int t = 0; t < Ntimesteps+1; t++){
+  for(int s = 0; s < output_pool_count; s++){
+    int pool_idx = CARDADATA.ncdf_data.POOLS_SUBSET_INDICES[s];
+    pool_subset[t * output_pool_count + s] = CARDADATA.M_POOLS[t * CARDADATA.nopools + pool_idx];
+  }
+}
+FAILONERROR(nc_put_vara_double(ncid,poolsVarID,(const size_t []){n,0,0}, (const size_t[]){1,Ntimesteps+1,output_pool_count}, pool_subset));
+free(pool_subset);
 FAILONERROR(nc_put_vara_double(ncid,parsVarID,(const size_t[]){n,0}, (const size_t[]){1,CARDADATA.nopars}, pars));
 
 
@@ -454,6 +532,27 @@ FAILONERROR(nc_close(ncid));
 
 /*Step 6: Free memory*/
 /*exhaustive list of all malloc/calloc used fields*/
+
+if (CARDADATA.ncdf_data.FLUXES_SUBSET_NAMES != NULL) {
+    for (int i = 0; i < CARDADATA.ncdf_data.FLUXES_SUBSET_COUNT; i++) {
+        free(CARDADATA.ncdf_data.FLUXES_SUBSET_NAMES[i]);
+    }
+    free(CARDADATA.ncdf_data.FLUXES_SUBSET_NAMES);
+}
+if (CARDADATA.ncdf_data.FLUXES_SUBSET_INDICES != NULL) {
+    free(CARDADATA.ncdf_data.FLUXES_SUBSET_INDICES);
+}
+
+if (CARDADATA.ncdf_data.POOLS_SUBSET_NAMES != NULL) {
+    for (int i = 0; i < CARDADATA.ncdf_data.POOLS_SUBSET_COUNT; i++) {
+        free(CARDADATA.ncdf_data.POOLS_SUBSET_NAMES[i]);
+    }
+    free(CARDADATA.ncdf_data.POOLS_SUBSET_NAMES);
+}
+if (CARDADATA.ncdf_data.POOLS_SUBSET_INDICES != NULL) {
+    free(CARDADATA.ncdf_data.POOLS_SUBSET_INDICES);
+}
+
 free(pars);
 FREE_DATA_STRUCT(CARDADATA);
 

C/projects/CARDAMOM_GENERAL/NETCDF_AUXILLIARY_FUNCTIONS.c

@@ -6,7 +6,6 @@
 #define DEFAULT_DOUBLE_VAL -9999.0
 #define DEFAULT_INT_VAL -9999
 
-
 //NOTE ABOUT THIS MACRO:
 //If set to 1, netCDF methods will continue to run and return with default values if they fail to find the requested variable or attribute
 //if set to 0, they will instantly die on failing to find any variable or attribute
@@ -298,6 +297,65 @@ double ** ncdf_read_double_2D(int ncid, const char * varName, size_t * dimLen ){
 }
 
 
+/*
+ * Function:  ncdf_read_string_array
+ * --------------------
+ * Reads a comma-delimited string attribute from netCDF
+ *
+ *  ncid: netCDF file ID to pull the data from
+ *  attrName: This is the name of the attribute to read
+ *  count: pointer where the number of strings will be written
+ *
+ *  returns: array of string pointers, or NULL if attribute doesn't exist
+ *   Reads a global attribute like "GPP,rh_co2,ets" and splits by commas
+ */
+char **ncdf_read_string_array(int ncid, const char *attrName, int *count) {
+	int retval = 0;
+	size_t attr_len;
+
+	if ((retval = nc_inq_attlen(ncid, NC_GLOBAL, attrName, &attr_len))) {
+		if (retval == NC_ENOTATT && ALLOW_DEFAULTS) {
+			*count = 0;
+			return NULL;
+		}
+		ERR_ATTR_AND_CONTEXT(retval, attrName, "/", NC_GLOBAL);
+	}
+
+	char *attr_value = calloc(attr_len + 1, sizeof(char));
+	if ((retval = nc_get_att_text(ncid, NC_GLOBAL, attrName, attr_value))) {
+		free(attr_value);
+		ERR_ATTR_AND_CONTEXT(retval, attrName, "/", NC_GLOBAL);
+	}
+	attr_value[attr_len] = '\0';
+
+	int num_strings = 1;
+	for (size_t i = 0; i < attr_len; i++) {
+		if (attr_value[i] == ',') num_strings++;
+	}
+
+	char **strings = calloc(num_strings, sizeof(char *));
+	int string_idx = 0;
+	char *token = strtok(attr_value, ",");
+	while (token != NULL && string_idx < num_strings) {
+		while (*token == ' ') token++;
+
+		strings[string_idx] = calloc(100, sizeof(char));
+		strncpy(strings[string_idx], token, 99);
+		strings[string_idx][99] = '\0';
+
+		size_t len = strlen(strings[string_idx]);
+		while (len > 0 && strings[string_idx][len-1] == ' ') {
+			strings[string_idx][--len] = '\0';
+		}
+
+		string_idx++;
+		token = strtok(NULL, ",");
+	}
+
+	free(attr_value);
+	*count = string_idx;
+	return strings;
+}
 
 
 

C/projects/CARDAMOM_MODELS/DALEC/DALEC_1100/DALEC_1100_NC_INFO.c

@@ -67,7 +67,7 @@ void POPULATE_INFO_STRUCTS(DALEC * DALECmodel){
 
     // Carbon, Water, Energy Fluxes
     DALECmodel->FLUX_META.NAME[F.gpp]="Gross Primary productivity";
-    DALECmodel->FLUX_META.ABBREVIATION[F.gpp]="GPP";
+    DALECmodel->FLUX_META.ABBREVIATION[F.gpp]="gpp";
     DALECmodel->FLUX_META.UNITS[F.gpp]="gC/m2/day";
     DALECmodel->FLUX_META.DESCRIPTION[F.gpp]="GPP, doesn\"t include maintenance respiration";