Merge branch 'feature/memory_fix' into feature/markross_PrototypeMast…

…erVersion2.9

bin/CMakeLists.txt

@@ -1,7 +1,7 @@
 include_directories(${PROJECT_SOURCE_DIR}/inc)
 
 add_executable (sbnfit_make_covariance sbnfit_make_covariance.cxx)
-target_link_libraries (sbnfit_make_covariance SBNfitlib sbnfit_tinyxmllib EventWeight)
+target_link_libraries (sbnfit_make_covariance SBNfitlib sbnfit_tinyxmllib EventWeight SBNfit_DeNan)
 
 add_executable (sbnfit_plot_covariance sbnfit_plot_covariance.cxx)
 target_link_libraries (sbnfit_plot_covariance SBNfitlib sbnfit_tinyxmllib EventWeight)

bin/sbnfit_feldman_cousins.cxx

@@ -80,7 +80,7 @@ int main(int argc, char* argv[])
     std::string tag = "TEST";
     std::string mode_option;
     bool bool_stat_only = false;
-    int number = 2500;
+    int number = -1;
     int grid_pt = 0;
     double random_number_seed = -1;
 
@@ -161,8 +161,12 @@ int main(int argc, char* argv[])
 
     std::cout<<"Begining FeldmanCousins for tag: "<<tag<<std::endl;
 
-    NGrid mygrid;
+    NGrid mygrid(xml);
 
+
+    // update to able grid configuration in the xml 
+    // keep this as a record of old ranges used in previous sensitivity evaluation
+    /*
     if(tag == "NuMuDis"){	
        //grid for numu disappearance
        mygrid.AddDimension("m4", -1, 1.05, 0.5);//0.05 FULL
@@ -188,7 +192,7 @@ int main(int argc, char* argv[])
       //mygrid.AddDimension("ue4", -1.5, 0.05, 0.05); //0.1 full inject
       mygrid.AddFixedDimension("um4",0.0);         //0.05
     }
-
+    */
 
     //Print the grid interesting bits
     mygrid.Print();

bin/sbnfit_make_covariance.cxx

@@ -152,11 +152,6 @@ int main(int argc, char* argv[])
         return 1;
     }
 
-    //gSystem->Load("/uboone/app/users/markrl/Hive_v3.0/hellstroms_hive/hive/root_linkdefs/loc/SL7/denan_cxx.so");
-
-    //std::string dict_location = "../libio/libEventWeight.so";
-    //std::cout<<"Trying to load dictionary: "<<dict_location<<std::endl;
-    //gSystem->Load(  (dict_location).c_str());
 
     /*************************************************************
      *************************************************************
@@ -208,17 +203,19 @@ int main(int argc, char* argv[])
         //Will be outputted in the form: SBNfit_covariance_plots_TAG.root
         example_covar.PrintMatricies(tag,minimal_matrix_plots);
 
-        //Constraint will be patched in shortly: mark
-        std::ofstream ratio_con;
-        std::vector<double> average_ratio;
-        ratio_con.open("ratio_list_"+tag+".txt",std::ofstream::trunc);
         if(constrain_mode){  
+            //Constraint will be patched in shortly: mark
+            std::ofstream ratio_con;
+            std::vector<double> average_ratio;
+            ratio_con.open("ratio_list_"+tag+".txt",std::ofstream::trunc);
+
             for (int i=0;i<example_covar.variations.size();i++){
                 std::cout<<"var " <<example_covar.variations.at(i)<<std::endl;
                 average_ratio=example_covar.DoConstraint(0,1,tag,i);
                 //ratio_con<<i<<" " <<average_ratio<<std::endl;
                 ratio_con<<"var "<<i<<" name "<<example_covar.variations.at(i)<<" events "<<average_ratio[0]<<" uncon "<<average_ratio[1]<<" con "<<average_ratio[2]<<" ratio "<<average_ratio[3]<<std::endl;
             }
+	    ratio_con.close();
 
             example_covar.DoConstraint_test(0,1,tag);
         }  

bin/sbnfit_singlephoton.cxx

@@ -215,7 +215,7 @@ int main(int argc, char* argv[])
 
     std::cout<<"Begining Single Photon module"<<std::endl;
 
-    NGrid mygrid, poly_grid;
+    NGrid mygrid(xml), poly_grid;
 
     // flat fit to all subchannels
     //mygrid.AddConstrainedDimension("All", 0.5, 1.5, 0.01, 1.19);   //0.1 FULL
@@ -231,11 +231,10 @@ int main(int argc, char* argv[])
     //mygrid.AddConstrainedDimension("NCPi0NotCoh", 0.5, 1.25, 0.05, 1.0);   //coarser grid
     //mygrid.AddConstrainedDimension("NCPi0Coh", 0, 8, 0.2, 1.0); 
     //mygrid.AddFixedDimension("NCPi0NotCoh", 1.19);   //fixed
-    mygrid.AddDimension("NCDelta", 0, 6, 3.18);
+    //mygrid.AddDimension("NCDelta", 0, 6, 1);
     //mygrid.AddDimension("NCDeltaLEE", 0, 5, 0.01 );
     //mygrid.AddDimension("NCDeltaLEE", 0, 2.5, 0.005 );
 
-
     poly_grid.AddConstrainedDimension("NCPi0NotCoh", -2.6, 1.0, 0.2, 1);  //zoomed in first order
     //poly_grid.AddConstrainedDimension("NCPi0NotCoh", -4.0, 2.0, 0.4, -1.1);  //first order
     //poly_grid.AddFixedDimension("NCPi0NotCoh", -1.05); // second order 
@@ -318,7 +317,7 @@ int main(int argc, char* argv[])
 
 	sp.PublicDataPrintOut({"NCDelta"}, {0});
 	sp.PublicDataPrintOut({"NCDelta"}, {1});
-	sp.PublicDataPrintOut({"NCDelta"}, {3.18});
+	//sp.PublicDataPrintOut({"NCDelta"}, {3.18});
     }
     else{
         std::cout << "Mode input is not identified, please try with a valid input.." << std::endl;

inc/SBNchi.h

@@ -73,6 +73,15 @@ namespace sbn{
             //Initialise a stat_only one;
             SBNchi(SBNspec, bool is_stat_only);
             SBNchi(std::string);
+	    ~SBNchi(){
+	std::cout << "start free SBNchi " << std::endl;
+		 if(rangen_twister)   {delete rangen_twister; rangen_twister = nullptr;}
+        	 if(rangen_linear)    {delete rangen_linear; rangen_linear = nullptr;}
+        	 if(rangen_carry)     {delete rangen_carry; rangen_carry = nullptr;}
+        	 if(rangen)	      {delete rangen; rangen = nullptr;}
+  		 if(m_dist_normal)  {delete  m_dist_normal; m_dist_normal = nullptr;}
+	std::cout << "end free SBNchi " << std::endl;
+	    }
 
 	//This is the core spectra that you are comparing too. This is used to calculate covariance matrix and in a way is on the 'bottom' of the chi^2.
 	SBNspec core_spectrum;
@@ -105,13 +114,13 @@ namespace sbn{
 
             /***** Random Number Generation ****/
             std::random_device random_device_seed;
-            std::mt19937 *rangen_twister; //merseinne twister
-            std::minstd_rand * rangen_linear;
-            std::ranlux24_base * rangen_carry;
+            std::mt19937 *rangen_twister = nullptr; //merseinne twister
+            std::minstd_rand * rangen_linear = nullptr;
+            std::ranlux24_base * rangen_carry = nullptr;
             void InitRandomNumberSeeds();
             void InitRandomNumberSeeds(double);
-            TRandom3 * rangen;
-            std::normal_distribution<float>* m_dist_normal;
+            TRandom3 * rangen = nullptr;
+            std::normal_distribution<float>* m_dist_normal = nullptr;
 
             /*********************************** Member Functions ********************************/	
 

inc/SBNconfig.h

@@ -14,6 +14,7 @@
 #include <algorithm>
 #include <iomanip>
 #include <exception>
+#include <memory>
 
 #include "TH1D.h"
 #include "TFile.h"
@@ -69,6 +70,10 @@ class SBNconfig {
 	SBNconfig(std::string);
 	SBNconfig(){};
 	SBNconfig(std::vector<std::string>, std::vector<std::string>, std::vector<std::string>, std::vector<std::vector<std::string>>, std::vector<std::vector<double>>);
+
+	//destructor
+	~SBNconfig();
+
 	//This is going to be a manual Setup thing
     int LoadFromXML(const char* filedata,bool,bool);
 
@@ -97,10 +102,10 @@ class SBNconfig {
 
 	//vectors of length num_channels
 	std::vector<int> num_subchannels; 
-	int* a_num_subchannels;
+	int* a_num_subchannels = nullptr;
 	std::vector<int> num_subchannels_xml; 
 	std::vector<int> num_bins;
-	int* a_num_bins;
+	int* a_num_bins = nullptr;
 
 	std::string xmlname;	
 
@@ -200,7 +205,7 @@ class SBNconfig {
 	std::vector<double> pot;
 
 	std::vector<std::vector<std::string>> parameter_names;	//obsolete code
-	std::vector<std::vector<BranchVariable*>> branch_variables;
+	std::vector<std::vector<std::shared_ptr<BranchVariable>>> branch_variables;
 
  	/**********************created for single photon****************************/
         //systematics root files provided correspond to

inc/SBNcovariance.h

@@ -60,7 +60,7 @@ namespace sbn{
     SBNcovariance(std::string xmlname);
     SBNcovariance(std::string xmlname, bool);
     SBNcovariance(std::string xmlname, std::string);
-
+    ~SBNcovariance();
 
     void ShapeOnlyProcessing();
     int FormCovarianceMatrix(std::string tag);
@@ -94,9 +94,9 @@ namespace sbn{
     TMatrixD input_frac_covariance;
 
     //for plotting, hsa been superseeded by PrintMatricies a bit
-    TH2D * hist_frac_cov;
-    TH2D * hist_full_cor;
-    TH2D * hist_full_cov;
+    TH2D * hist_frac_cov = nullptr;
+    TH2D * hist_full_cor = nullptr;
+    TH2D * hist_full_cov = nullptr;
 
     //Some checks on montecarlos
     double tolerence_positivesemi;
@@ -135,9 +135,7 @@ namespace sbn{
     std::map<std::string,bool> m_variations_to_use;
 
     // In testing 2D fits with 4D covariance..
-    SBNspec template_spec;	
-    SBNspec spec_central_value2;	
-    THnSparseD *frac_covariance_4d;
+    THnSparseD *frac_covariance_4d = nullptr;
     TMatrixD full_covariance2D;
     TMatrixD frac_covariance2D;
     TMatrixD full_correlation2D;
@@ -194,6 +192,8 @@ namespace sbn{
 	     }
 	   }
        }
+     delete m_check1; m_check1 = nullptr;
+     delete m_check2; m_check2 = nullptr;
      return true;
      }
 

inc/SBNspec.h

@@ -71,7 +71,7 @@ namespace sbn{
 			SBNspec(std::vector<double> input_full_vec, std::vector<double> input_full_err, std::string whichxml);
 			SBNspec(std::vector<double> input_full_vec, std::vector<double> input_full_err, std::string whichxml, bool isverbose);
             		SBNspec(std::vector<double> input_full_vec, std::vector<double> input_full_err, std::string whichxml, int universe, bool isverbose);
-
+			~SBNspec(){std::cout << "free SBNspec " << std::endl;}
 
             std::map<int,std::vector<int>> GetCollapsedChannelIndicies();
             std::vector<TH1D> GetBlankChannelHists();

inc/branch_variable.h

@@ -50,6 +50,20 @@ struct BranchVariable{
 	int SetOscillate(bool inbool){ oscillate = inbool;};
 	bool GetOscillate(){ return oscillate;};
 
+       ~BranchVariable(){
+	    if(branch_formula){
+		  std::cout << " BranchVariable || delete reco formula " << branch_formula << std::endl;
+		  delete branch_formula; branch_formula = nullptr;
+	    }
+	    if(branch_true_value_formula){
+ 		std::cout << " BranchVariable || delete true value formula " << branch_true_value_formula << std::endl;
+		delete branch_true_value_formula; branch_true_value_formula = nullptr;
+  	    }
+	    if(branch_true_L_formula){
+		std::cout << "BranchVariable || Delete true L formula " << branch_true_L_formula << std::endl;
+		delete branch_true_L_formula; branch_true_L_formula = nullptr;
+	    }
+       }
 };
 
 /*struct BranchVariable_i: public BranchVariable{

inc/ngrid.h

@@ -8,9 +8,14 @@
 #include <sstream>
 #include <map>
 #include <time.h>
+#include <stdexcept>
+#include <cstdlib>
 
 #include "TMatrixT.h"
 
+//--- local header ----
+#include "tinyxml.h"
+
 struct SBNfitResult{
 
     size_t bf_pt;
@@ -101,14 +106,62 @@ struct NGridDimension{
 };
 
 struct NGrid{
-    int f_num_dimensions;
-    int f_num_total_points;
+    int f_num_dimensions = 0;
+    int f_num_total_points = 1;
 
     std::vector<NGridDimension> f_dimensions;
 
-    NGrid(){
-        f_num_dimensions=0;
-        f_num_total_points=1;
+    //constructors 
+    NGrid(){}
+
+    // given provided configuration xml, grab the grids in it
+    NGrid(std::string xml){
+
+	std::cout << "Load Grid using configuration xml: " << xml << std::endl;
+
+	TiXmlDocument doc(xml.c_str());
+	bool loadOkay = doc.LoadFile();
+        if(!loadOkay){
+	    throw std::runtime_error("\tERROR: fail to load configuration XML, generally means broken .xml brackets or attribute syntax.");
+        }
+        TiXmlHandle hDoc(&doc);
+	TiXmlElement *pGrid;
+	pGrid = doc.FirstChildElement("GridDimension");
+
+	while(pGrid){
+	    const char* dim_name = pGrid->Attribute("name");
+	    const char* dim_fixval = pGrid->Attribute("fix_val");
+	    const char* dim_lowedge = pGrid->Attribute("min");
+	    const char* dim_upedge = pGrid->Attribute("max");
+	    const char* dim_stepsize = pGrid->Attribute("step");
+	    char* pEnd; 
+
+            if(dim_name==NULL){
+                std::cerr << "\tERROR! Grid dimension need a name! Please kindly define a name for it"<<std::endl;
+		print_config_rule();
+		throw std::runtime_error("Invalid Syntax to setup Grid");
+            }	
+	    std::cout << "--> On Grid dimension " << dim_name;
+
+	    if(dim_fixval){
+		std::cout << " with fixed value : " << dim_fixval << std::endl;
+	 	this->AddFixedDimension(std::string(dim_name), strtod(dim_fixval, &pEnd));	
+	    }else{
+
+		if(dim_lowedge == NULL || dim_upedge == NULL || dim_stepsize == NULL){
+		    print_config_rule();
+		    throw std::runtime_error("Invalid Syntax to setup Grid");
+		}
+
+		std::cout << " with min: " << dim_lowedge << " max: " << dim_upedge << " and step size: " << dim_stepsize << std::endl;
+		this->AddDimension(std::string(dim_name), strtod(dim_lowedge, &pEnd),  strtod(dim_upedge, &pEnd),  strtod(dim_stepsize, &pEnd));
+	    }
+
+            pGrid = pGrid->NextSiblingElement("GridDimension");
+        }
+
+	std::cout << "Finish reading grids from xml... Total of " << f_num_dimensions << " dimensions."  << std::endl;
+
     }
 
     void AddDimension(std::string name, double min, double max, double step ){
@@ -225,6 +278,15 @@ struct NGrid{
 
     }
 
+    void print_config_rule(){
+	std::cerr << "\n"<< std::string(30, '_') <<  std::endl;
+	std::cerr << "\tRULEs for xml-based grid setup" << std::endl;
+	std::cerr << "\tPlease specify the name, lower, upper range for the grid, as well as the stepsize. Please follow style below: " << std::endl; 
+	std::cerr << "\t\tOption 1: name=\"xx\" min=\"xx\" max=\"xx\" step=\"xx\" to setup evenly distributed grid" << std::endl;
+	std::cerr << "\t\tOption 2: name=\"xx\" fix_val=\"xx\" to setup fixed value for this dimension" << std::endl;
+	std::cerr << "\n"<< std::string(30, '_') <<  std::endl;
+ 	return;	
+    }
 };
 
 #endif

libio/CMakeLists.txt

@@ -6,9 +6,15 @@ unset(CMAKE_CXX_FLAGS)
 set(CMAKE_CXX_FLAGS " ${CMAKE_SO_CXX_FLAGS}")
 message("@libio CMAKE_CXX_FLAGS: ${CMAKE_CXX_FLAGS}")
 
+#--- create eventweight library
 ROOT_GENERATE_DICTIONARY(G__EventWeight eventweight.h LINKDEF eventweight_LinkDef.h)
 
 #---Create a shared library with geneated dictionary
 add_library(EventWeight SHARED eventweight.cxx G__EventWeight.cxx)
 target_link_libraries(EventWeight ${ROOT_LIBRARIES})
 
+
+#--- create DeNan library ---
+ROOT_GENERATE_DICTIONARY(DeNan_Dict denan.h LINKDEF denan_LinkDef.h)
+add_library(SBNfit_DeNan SHARED denan.cxx DeNan_Dict.cxx)
+target_link_libraries(SBNfit_DeNan ${ROOT_LIBRARIES})

libio/denan.cxx

@@ -0,0 +1,10 @@
+#include "denan.h"
+#include "denan_LinkDef.h"
+
+double DeNan(double x, double replace){
+    if(x!=x || isinf(x))return replace;
+    return x;
+}
+
+
+

libio/denan.h

@@ -0,0 +1,9 @@
+#ifndef DENAN_H
+#define DENAN_H
+
+#include <cmath>
+
+double DeNan(double x, double replace);
+
+
+#endif