diff --git a/src/FitBase/MeasurementBase.cxx b/src/FitBase/MeasurementBase.cxx
index a77597e..235b470 100644
--- a/src/FitBase/MeasurementBase.cxx
+++ b/src/FitBase/MeasurementBase.cxx
@@ -1,561 +1,562 @@
 // Copyright 2016 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret
 
 /*******************************************************************************
 *    This file is part of NUISANCE.
 *
 *    NUISANCE is free software: you can redistribute it and/or modify
 *    it under the terms of the GNU General Public License as published by
 *    the Free Software Foundation, either version 3 of the License, or
 *    (at your option) any later version.
 *
 *    NUISANCE is distributed in the hope that it will be useful,
 *    but WITHOUT ANY WARRANTY; without even the implied warranty of
 *    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *    GNU General Public License for more details.
 *
 *    You should have received a copy of the GNU General Public License
 *    along with NUISANCE.  If not, see <http://www.gnu.org/licenses/>.
 *******************************************************************************/
 
 #include "MeasurementBase.h"
 
 /*
   Constructor/Destructors
 */
 
 //********************************************************************
 // 2nd Level Constructor (Inherits From MeasurementBase.h)
 MeasurementBase::MeasurementBase(void) {
   //********************************************************************
 
   fScaleFactor = 1.0;
   fMCFilled = false;
   fNoData = false;
   fInput = NULL;
 
   // Set the default values
   // After-wards this gets set in SetupMeasurement
   EnuMin = 0.;
   EnuMax = 1.E5;
 
   fMeasurementSpeciesType = kSingleSpeciesMeasurement;
   fEventVariables = NULL;
   fIsJoint = false;
 };
 
 void MeasurementBase::FinaliseMeasurement() {
 
   // Used to setup default data hists, covars, etc.
 
 
 
 }
 
 //********************************************************************
 // 2nd Level Destructor (Inherits From MeasurementBase.h)
 MeasurementBase::~MeasurementBase() {
   //********************************************************************
 
 };
 
 //********************************************************************
 double MeasurementBase::TotalIntegratedFlux(std::string intOpt, double low,
     double high) {
 //********************************************************************
 
   // Set Energy Limits
   if (low == -9999.9) low = this->EnuMin;
   if (high == -9999.9) high = this->EnuMax;
   return GetInput()->TotalIntegratedFlux(low, high, intOpt);
 };
 
 //********************************************************************
 double MeasurementBase::PredictedEventRate(std::string intOpt, double low,
     double high) {
   //********************************************************************
 
   // Set Energy Limits
   if (low == -9999.9) low = this->EnuMin;
   if (high == -9999.9) high = this->EnuMax;
 
   return GetInput()->PredictedEventRate(low, high, intOpt) * 1E-38;
 };
 
 //********************************************************************
 void MeasurementBase::SetupInputs(std::string inputfile) {
   //********************************************************************
 
 
   // Add this infile to the global manager
   if (FitPar::Config().GetParB("EventManager")) {
     fInput = FitBase::AddInput(fName, inputfile);
   } else {
     std::vector<std::string> file_descriptor =
       GeneralUtils::ParseToStr(inputfile, ":");
     if (file_descriptor.size() != 2) {
       ERR(FTL) << "File descriptor had no filetype declaration: \"" << inputfile
                << "\". expected \"FILETYPE:file.root\"" << std::endl;
       throw;
     }
     InputUtils::InputType inpType =
       InputUtils::ParseInputType(file_descriptor[0]);
 
     fInput = InputUtils::CreateInputHandler(fName, inpType, file_descriptor[1]);
   }
 
 
   fNEvents = fInput->GetNEvents();
 
   // Expect INPUTTYPE:FileLocation(s)
   std::vector<std::string> file_descriptor =
     GeneralUtils::ParseToStr(inputfile, ":");
   if (file_descriptor.size() != 2) {
     ERR(FTL) << "File descriptor had no filetype declaration: \"" << inputfile
              << "\". expected \"FILETYPE:file.root\"" << std::endl;
     throw;
   }
   fInputType = InputUtils::ParseInputType(file_descriptor[0]);
 
   fInputFileName = file_descriptor[1];
   if (EnuMin == 0 && EnuMax == 1.E5) {
     EnuMin = fInput->GetFluxHistogram()->GetBinLowEdge(1);
     EnuMax = fInput->GetFluxHistogram()->GetBinLowEdge(
                fInput->GetFluxHistogram()->GetNbinsX() + 1);
   }
 
   fFluxHist = fInput->GetFluxHistogram();
   fEventHist = fInput->GetEventHistogram();
 }
 
 //***********************************************
 int MeasurementBase::GetInputID() {
 //***********************************************
   return FitBase::GetInputID(fInputFileName);
 }
 
 //***********************************************
 SampleSettings MeasurementBase::LoadSampleSettings(nuiskey samplekey) {
 //***********************************************
   SampleSettings setting = SampleSettings(samplekey);
   fName = setting.GetS("name");
 
   // Used as an initial setup function incase we need to do anything here.
   LOG(SAM) << "Loading Sample : " << setting.GetName() << std::endl;
   if (!fIsJoint) SetupInputs( setting.GetS("input") );
 
   return setting;
 }
 
 //***********************************************
 SampleSettings MeasurementBase::LoadSampleSettings(std::string name, std::string input, std::string type) {
 //***********************************************
 
   nuiskey samplekey = Config::CreateKey("sample");
   samplekey.SetS("name",name);
   samplekey.SetS("input",input);
   samplekey.SetS("type",type);
 
   return LoadSampleSettings(samplekey);
 }
 
 void MeasurementBase::FinaliseSampleSettings() {
 
   EnuMin = fSettings.GetD("enu_min");
   EnuMax = fSettings.GetD("enu_max");
 
 }
 
 
 //***********************************************
 void MeasurementBase::Reconfigure() {
 //***********************************************
   
   LOG(REC) << " Reconfiguring sample " << fName << std::endl;
 
   // Reset Histograms
   ResetExtraHistograms();
   AutoResetExtraTH1();
   this->ResetAll();
 
   // FitEvent* cust_event = fInput->GetEventPointer();
   int fNEvents = fInput->GetNEvents();
   int countwidth = (fNEvents / 5);
 
 
   // MAIN EVENT LOOP
   FitEvent* cust_event = fInput->FirstNuisanceEvent();
   int i = 0;
   int npassed = 0;
   while(cust_event){
 
     cust_event->RWWeight = fRW->CalcWeight(cust_event);
     cust_event->Weight = cust_event->RWWeight * cust_event->InputWeight;
 
     Weight = cust_event->Weight;
 
     // Initialize
     fXVar = -999.9;
     fYVar = -999.9;
     fZVar = -999.9;
     Signal = false;
     Mode = cust_event->Mode;
 
     // Extract Measurement Variables
     this->FillEventVariables(cust_event);
     Signal = this->isSignal(cust_event);
     if (Signal) npassed++;
 
     GetBox()->SetX(fXVar);
     GetBox()->SetY(fYVar);
     GetBox()->SetZ(fZVar);
     GetBox()->SetMode(Mode);
     // GetBox()->fSignal = Signal;
 
     // Fill Histogram Values
     GetBox()->FillBoxFromEvent(cust_event);
     // this->FillExtraHistograms(GetBox(), Weight);
     this->FillHistogramsFromBox(GetBox(), Weight);
 
     // Print Out
     if (LOG_LEVEL(REC) && countwidth > 0 && !(i % countwidth)) {
       std::stringstream ss("");
       ss.unsetf(std::ios_base::fixed);
       ss << std::setw(7) << std::right << i << "/" << fNEvents << " events ("
          << std::setw(2) << double(i) / double(fNEvents) * 100. << std::left
          << std::setw(5) << "%) "
          << "[S,X,Y,Z,M,W] = [" << std::fixed << std::setprecision(2)
          << std::right << Signal << ", " << std::setw(5) << fXVar << ", "
          << std::setw(5) << fYVar << ", " << std::setw(5) << fYVar << ", "
          << std::setw(3) << (int)Mode << ", " << std::setw(5) << Weight << "] "
          << std::endl;
       LOG(SAM) << ss.str();
     }
 
     // iterate
     cust_event = fInput->NextNuisanceEvent();
     i++;
   }
 
   LOG(SAM) << npassed << "/" << fNEvents << " passed selection " << std::endl;
   if (npassed == 0) {
     LOG(SAM) << "WARNING: NO EVENTS PASSED SELECTION!" << std::endl;
   }
   LOG(REC) << std::setw(10) << std::right << NSignal << "/"
            << fNEvents << " events passed selection + binning after reweight"
            << std::endl;
 
   // Finalise Histograms
   fMCFilled = true;
   this->ConvertEventRates();
 }
 
 void MeasurementBase::FillHistogramsFromBox(MeasurementVariableBox* var, double weight) {
 
   fXVar  = var->GetX();
   fYVar  = var->GetY();
   fZVar  = var->GetZ();
   // Signal = var->fSignal;
   // Mode   = var->fMode;
   Weight = weight;
+  fEventVariables = var;
 
   FillHistograms();
   FillExtraHistograms(var, weight);
 
 }
 
 void MeasurementBase::FillHistograms(double weight){
   Weight = weight * GetBox()->GetSampleWeight();
   FillHistograms();
   FillExtraHistograms(GetBox(), Weight);
 }
 
 
 MeasurementVariableBox* MeasurementBase::FillVariableBox(FitEvent* event) {
 
   GetBox()->Reset();
   Mode = event->Mode;
   Weight = 1.0; //event->Weight;
 
   this->FillEventVariables(event);
   Signal = this->isSignal(event);
 
   GetBox()->FillBoxFromEvent(event);
   
   GetBox()->SetX(fXVar);
   GetBox()->SetY(fYVar);
   GetBox()->SetZ(fZVar);
   GetBox()->SetMode(event->Mode);
   GetBox()->SetSampleWeight(Weight);
   // GetBox()->fSignal = Signal;
 
   return GetBox();
 }
 
 MeasurementVariableBox* MeasurementBase::GetBox() {
   if (!fEventVariables) fEventVariables = CreateBox();
   return fEventVariables;
 }
 
 //***********************************************
 void MeasurementBase::ReconfigureFast() {
   //***********************************************
   this->Reconfigure();
 }
 
 //***********************************************
 void MeasurementBase::ConvertEventRates() {
   //***********************************************
 
   AutoScaleExtraTH1();
   ScaleExtraHistograms(GetBox());
   this->ScaleEvents();
 
   double normval = fRW->GetSampleNorm(this->fName);
   if (normval < 0.01 or normval > 10.0){
     ERR(WRN) << "Norm Value inside MeasurementBase::ConvertEventRates() looks off!" << std::endl;
     ERR(WRN) << "It could have become out of sync with the minimizer norm list." << std::endl;
     ERR(WRN) << "Setting it to 1.0" << std::endl;
     normval = 1.0;
   }
   AutoNormExtraTH1(normval);
   NormExtraHistograms(GetBox(), normval);
   this->ApplyNormScale(normval);
 
 }
 
 //***********************************************
 InputHandlerBase* MeasurementBase::GetInput() {
   //***********************************************
 
   if (!fInput) {
     ERR(FTL) << "MeasurementBase::fInput not set. Please submit your command "
              "line options and input cardfile with a bug report to: "
              "nuisance@projects.hepforge.org"
              << std::endl;
     throw;
   }
   return fInput;
 };
 
 //***********************************************
 void MeasurementBase::Renormalise() {
   //***********************************************
 
   // Called when the fitter has changed a measurements normalisation but not any
   // reweight dials
   // Means we don't have to call the time consuming reconfigure when this
   // happens.
   double norm = fRW->GetDialValue(this->fName + "_norm");
 
   if ((this->fCurrentNorm == 0.0 and norm != 0.0) or not fMCFilled) {
     this->ReconfigureFast();
     return;
   }
 
   if (this->fCurrentNorm == norm) return;
 
   this->ApplyNormScale(1.0 / this->fCurrentNorm);
   this->ApplyNormScale(norm);
 
   return;
 };
 
 //***********************************************
 void MeasurementBase::SetSignal(bool sig) {
   //***********************************************
   Signal = sig;
 }
 
 //***********************************************
 void MeasurementBase::SetSignal(FitEvent* evt) {
   //***********************************************
   Signal = this->isSignal(evt);
 }
 
 //***********************************************
 void MeasurementBase::SetWeight(double wght) {
   //***********************************************
   Weight = wght;
 }
 
 //***********************************************
 void MeasurementBase::SetMode(int md) {
   //***********************************************
   Mode = md;
 }
 
 //***********************************************
 std::vector<TH1*> MeasurementBase::GetFluxList() {
   //***********************************************
   return GetInput()->GetFluxList();
 }
 
 //***********************************************
 std::vector<TH1*> MeasurementBase::GetEventRateList() {
   //***********************************************
   return GetInput()->GetEventList();
 }
 
 //***********************************************
 std::vector<TH1*> MeasurementBase::GetXSecList() {
   //***********************************************
   return GetInput()->GetXSecList();
 }
 
 
 void MeasurementBase::ProcessExtraHistograms(int cmd,
     MeasurementVariableBox* vars,
     double weight) {
   // This should be overriden if we have extra histograms!!!
   // Add a flag to tell user this...
   return;
 }
 
 void MeasurementBase::FillExtraHistograms(MeasurementVariableBox* vars,
     double weight) {
   ProcessExtraHistograms(kCMD_Fill, vars, weight);
 }
 
 void MeasurementBase::ScaleExtraHistograms(MeasurementVariableBox* vars) {
   ProcessExtraHistograms(kCMD_Scale, vars, 1.0);
 }
 void MeasurementBase::ResetExtraHistograms() {
   ProcessExtraHistograms(kCMD_Reset, NULL, 1.0);
 }
 void MeasurementBase::NormExtraHistograms(MeasurementVariableBox* vars,
     double norm) {
   ProcessExtraHistograms(kCMD_Norm, vars, norm);
 }
 void MeasurementBase::WriteExtraHistograms() {
   ProcessExtraHistograms(kCMD_Write, NULL, 1.00);
 }
 
 void MeasurementBase::SetAutoProcessTH1(TH1* hist, int c1, int c2, int c3, int c4, int c5) {
   FakeStack* fake = new FakeStack(hist);
   SetAutoProcessTH1(fake, c1, c2, c3, c4, c5); // Need to add a destroy command!
 }
 
 void MeasurementBase::SetAutoProcess(TH1* hist, int c1, int c2, int c3, int c4, int c5) {
   FakeStack* fake = new FakeStack(hist);
   SetAutoProcessTH1(fake, c1, c2, c3, c4, c5); // Need to add a destroy command!                                                                                                                                                            
 }
 
 void MeasurementBase::SetAutoProcess(TGraph* g, int c1, int c2, int c3, int c4, int c5) {
   FakeStack* fake = new FakeStack(g);
   SetAutoProcessTH1(fake, c1, c2, c3, c4, c5); // Need to add a destroy command!                                                                                                                                                            
 }
 
 void MeasurementBase::SetAutoProcess(TF1* f, int c1, int c2, int c3, int c4, int c5) {
   FakeStack* fake = new FakeStack(f);
   SetAutoProcessTH1(fake, c1, c2, c3, c4, c5); // Need to add a destroy command!                                                                                                                                                            
 }
 
 void MeasurementBase::SetAutoProcess(StackBase* hist,  int c1, int c2, int c3, int c4, int c5){
   SetAutoProcessTH1(hist, c1, c2, c3, c4, c5);
 }
 
 void MeasurementBase::SetAutoProcessTH1(StackBase* hist,  int c1, int c2, int c3, int c4, int c5) {
 
   // Set Defaults
   // int ncommands = kCMD_extraplotflags;
   bool autoflags[5];
   autoflags[0] = false;
   autoflags[1] = false;
   autoflags[2] = false;
   autoflags[3] = false;
   autoflags[4] = false;
 
   int givenflags[5];
   givenflags[0] = c1;
   givenflags[1] = c2;
   givenflags[2] = c3;
   givenflags[3] = c4;
   givenflags[4] = c5;
   fExtraTH1s[hist] = std::vector<int>(5,0);
 
   // Setup a default one.
   if (c1 == -1 && c2 == -1 && c3 == -1 && c4 == -1 && c5 == -1){
     fExtraTH1s[hist][kCMD_Reset] = 1;
     fExtraTH1s[hist][kCMD_Scale] = 1;
     fExtraTH1s[hist][kCMD_Norm] = 1;
     fExtraTH1s[hist][kCMD_Write] = 1;
   }
 
   for (int i = 0; i < 5; i++) {
     switch (givenflags[i]) {
 
     // Skip over...
     case -1:
       break;
 
     case kCMD_Reset:
     case kCMD_Scale:
     case kCMD_Norm:
     case kCMD_Write:
       fExtraTH1s[hist][givenflags[i]] = 1;
       break;
 
     case kCMD_Fill:
       ERR(FTL) << "Can't auto fill yet!" << std::endl;
       autoflags[givenflags[i]] = 1;
       break;
 
     default:
       break;
     }
   }
 
   // LOG(SAM) << "AutoProcessing " << hist->GetName() << std::endl;
 };
 
 void MeasurementBase::AutoFillExtraTH1() {
   ERR(FTL) << "Can't auto fill yet! it's too inefficent!" << std::endl;
   return;
 }
 
 void MeasurementBase::AutoResetExtraTH1() {
 
   for (std::map<StackBase*, std::vector<int> >::iterator iter = fExtraTH1s.begin();
        iter != fExtraTH1s.end(); iter++) {
 
     if (!((*iter).second)[kCMD_Reset]) continue;
     (*iter).first->Reset();
   }
 };
 
 void MeasurementBase::AutoScaleExtraTH1() {
   for (std::map<StackBase*, std::vector<int> >::iterator iter = fExtraTH1s.begin();
        iter != fExtraTH1s.end(); iter++) {
 
     if (!((*iter).second)[kCMD_Scale]) continue;
     if (fIsNoWidth){
       (*iter).first->Scale(fScaleFactor);
     } else {
       (*iter).first->Scale(fScaleFactor, "width");
     }
   }
 };
 
 void MeasurementBase::AutoNormExtraTH1(double norm) {
   double sfactor = 0.0;
   if (norm != 0.0) sfactor = 1.0 / norm;
 
   for (std::map<StackBase*, std::vector<int> >::iterator iter = fExtraTH1s.begin();
        iter != fExtraTH1s.end(); iter++) {
 
     if (!((*iter).second)[kCMD_Norm]) continue;
     (*iter).first->Scale(sfactor);
   }
 };
 
 void MeasurementBase::AutoWriteExtraTH1() {
   for (std::map<StackBase*, std::vector<int> >::iterator iter = fExtraTH1s.begin();
        iter != fExtraTH1s.end(); iter++) {
 
     if (!(((*iter).second)[kCMD_Write])) continue;
     (*iter).first->Write();
   }
 };
 
 
 
diff --git a/src/MINERvA/MINERvAVariableBoxes.h b/src/MINERvA/MINERvAVariableBoxes.h
index 56288f8..06ed9ed 100644
--- a/src/MINERvA/MINERvAVariableBoxes.h
+++ b/src/MINERvA/MINERvAVariableBoxes.h
@@ -1,30 +1,57 @@
 #ifndef MINERvA_VARIABLES_BOX_H
 #define MINERvA_VARIABLES_BOX_H
 
 #include "MeasurementVariableBox.h"
 #include "MeasurementVariableBox1D.h"
 #include "MeasurementVariableBox2D.h"
 
 /*!
  *  \addtogroup FitBase
  *  @{
  */
 
 /// Custom box used to also save All Pion Tpi for each event.
 class NTpiVariableBox1D : public MeasurementVariableBox1D {
 public:
-	inline NTpiVariableBox1D() { Reset(); };
+	inline NTpiVariableBox1D() { };
 	inline void Reset() {  fTpiVect.clear(); } 
+
+    inline MeasurementVariableBox* CloneSignalBox(){
+        NTpiVariableBox1D* box = new NTpiVariableBox1D(); 
+        box->fX = this->fX;
+        box->fSampleWeight = this->fSampleWeight;
+
+        box->fTpiVect.clear();
+        for (int i = 0; i < this->fTpiVect.size(); i++){
+            box->fTpiVect.push_back( this->fTpiVect[i] );
+        }
+        return box;
+    }
+    inline void Print(){
+        std::cout << "Box Print Size : " << this->fTpiVect.size() << std::endl;
+    }
+
 	std::vector<double> fTpiVect;
 };
 
 /// Custom box used to also save All Pion Tpi for each event.
 class NthpiVariableBox1D : public MeasurementVariableBox1D {
  public:
-  inline NthpiVariableBox1D() { Reset(); };
+  inline NthpiVariableBox1D() {  };
   inline void Reset() {  fthpiVect.clear(); }
+  inline MeasurementVariableBox* CloneSignalBox(){
+        NthpiVariableBox1D* box = new NthpiVariableBox1D(); 
+        box->fX = this->fX;
+        box->fSampleWeight = this->fSampleWeight;
+
+        box->fthpiVect.clear();
+        for (int i = 0; i < this->fthpiVect.size(); i++){
+            box->fthpiVect.push_back( this->fthpiVect[i] );
+        }
+        return box;
+    }
   std::vector<double> fthpiVect;
 };
 
 
 #endif
diff --git a/src/MINERvA/MINERvA_CCNpip_XSec_1DTpi_nu.cxx b/src/MINERvA/MINERvA_CCNpip_XSec_1DTpi_nu.cxx
index 51fa133..74f9d9a 100644
--- a/src/MINERvA/MINERvA_CCNpip_XSec_1DTpi_nu.cxx
+++ b/src/MINERvA/MINERvA_CCNpip_XSec_1DTpi_nu.cxx
@@ -1,267 +1,267 @@
 // Copyright 2016 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret
 
 /*******************************************************************************
 *    This file is part of NUISANCE.
 *
 *    NUISANCE is free software: you can redistribute it and/or modify
 *    it under the terms of the GNU General Public License as published by
 *    the Free Software Foundation, either version 3 of the License, or
 *    (at your option) any later version.
 *
 *    NUISANCE is distributed in the hope that it will be useful,
 *    but WITHOUT ANY WARRANTY; without even the implied warranty of
 *    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *    GNU General Public License for more details.
 *
 *    You should have received a copy of the GNU General Public License
 *    along with NUISANCE.  If not, see <http://www.gnu.org/licenses/>.
 *******************************************************************************/
 
 #include "MINERvA_SignalDef.h"
 #include "MINERvA_CCNpip_XSec_1DTpi_nu.h"
 
 //********************************************************************
 MINERvA_CCNpip_XSec_1DTpi_nu::MINERvA_CCNpip_XSec_1DTpi_nu(nuiskey samplekey) {
 //********************************************************************
 
   // Sample overview ---------------------------------------------------
   std::string descrip = "MINERvA_CCNpip_XSec_1DTpi_nu sample. \n" \
                         "Target: CH \n" \
                         "Flux: MINERvA Forward Horn Current nue + nuebar \n" \
                         "Signal: Any event with 1 electron, any nucleons, and no other FS particles \n";
 
   // Setup common settings
   fSettings = LoadSampleSettings(samplekey);
   fSettings.SetDescription(descrip);
   fSettings.SetXTitle("T_{#pi} (MeV)");
   fSettings.SetYTitle("(1/T#Phi) dN_{#pi}/dT_{#pi} (cm^{2}/MeV/nucleon)");
   fSettings.SetAllowedTypes("FIX,FREE,SHAPE/DIAG,FULL/NORM/MASK", "FIX/FULL");
   fSettings.SetEnuRange(1.5, 10.0);
   fSettings.DefineAllowedTargets("C,H");
   fSettings.DefineAllowedSpecies("numu");
 
   fFullPhaseSpace = !fSettings.Found("name", "_20deg");
   fFluxCorrection = fSettings.Found("name", "fluxcorr");
   fUpdatedData    = !fSettings.Found("name", "2015");
   fSettings.SetTitle("MINERvA_CCNpip_XSec_1DTpi_nu");
 
   FinaliseSampleSettings();
 
   // Scaling Setup ---------------------------------------------------
   // ScaleFactor automatically setup for DiffXSec/cm2/Nucleon
   fScaleFactor = GetEventHistogram()->Integral("width") * double(1E-38) / double(fNEvents) / TotalIntegratedFlux("width");
 
   // Plot Setup -------------------------------------------------------
   // Full Phase Space
   if (fFullPhaseSpace) {
     // 2016 release
     if (fUpdatedData) {
       SetDataFromTextFile(GeneralUtils::GetTopLevelDir() + "/data/MINERvA/CCNpip/2016/nu-ccNpi+-xsec-pion-kinetic-energy.csv");
 
       // MINERvA has the error quoted as a percentage of the cross-section
       // Need to make this into an absolute error before we go from correlation
       // matrix -> covariance matrix since it depends on the error in the ith
       // bin
       for (int i = 0; i < fDataHist->GetNbinsX() + 1; i++) {
         fDataHist->SetBinError(i + 1, fDataHist->GetBinContent(i + 1) * (fDataHist->GetBinError(i + 1) / 100.));
       }
 
       // This is a correlation matrix, not covariance matrix, so needs to be
       // converted
       SetCorrelationFromTextFile(GeneralUtils::GetTopLevelDir() + "/data/MINERvA/CCNpip/2016/nu-ccNpi+-correlation-pion-kinetic-energy.csv");
 
 
       // 2015 release
     } else {
       // If we're doing shape only
       if (fIsShape) {
         SetDataFromTextFile(GeneralUtils::GetTopLevelDir() + "/data/MINERvA/CCNpip/2015/MINERvA_CCNpi_Tpi_shape.txt");
         SetCorrelationFromTextFile(GeneralUtils::GetTopLevelDir() + "/data/MINERvA/CCNpip/2015/MINERvA_CCNpi_Tpi_shape_cov.txt");
 
         // If we're doing full cross-section
       } else {
         SetDataFromTextFile(GeneralUtils::GetTopLevelDir() + "/data/MINERvA/CCNpip/2015/MINERvA_CCNpi_Tpi.txt");
         SetCorrelationFromTextFile( GeneralUtils::GetTopLevelDir() +  "/data/MINERvA/CCNpip/2015/MINERvA_CCNpi_Tpi_cov.txt");
       }
 
     }
 
     // Restricted Phase Space
   } else {
     // Only 2015 data released restricted muon phase space cross-section
     // unfortunately
     if (fUpdatedData) {
       ERR(FTL) << fName << " has no updated 2016 data for restricted phase space! Using 2015 data." << std::endl;
       throw;
     }
 
     // If we're using the shape only data
     if (fIsShape) {
       SetDataFromTextFile(GeneralUtils::GetTopLevelDir() + "/data/MINERvA/CCNpip/2015/MINERvA_CCNpi_Tpi_20deg_shape.txt");
       SetCorrelationFromTextFile(GeneralUtils::GetTopLevelDir() + "/data/MINERvA/CCNpip/2015/MINERvA_CCNpi_Tpi_20deg_shape_cov.txt");
 
       // Or total cross-section
     } else {
       SetDataFromTextFile(GeneralUtils::GetTopLevelDir() + "/data/MINERvA/CCNpip/2015/MINERvA_CCNpi_Tpi_20deg.txt");
       SetCorrelationFromTextFile(GeneralUtils::GetTopLevelDir() + "/data/MINERvA/CCNpip/2015/MINERvA_CCNpi_Tpi_20deg_cov.txt");
     }
   }
 
   // Scale the MINERvA data to account for the flux difference
   // Adjust MINERvA data to flux correction; roughly a 11% normalisation increase in data
   // Please change when MINERvA releases new data!
   if (fFluxCorrection) {
     for (int i = 0; i < fDataHist->GetNbinsX() + 1; i++) {
       fDataHist->SetBinContent(i + 1, fDataHist->GetBinContent(i + 1) * 1.11);
     }
   }
 
   // Make some auxillary helper plots
   onePions  = (TH1D*)(fDataHist->Clone());
   onePions->SetNameTitle((fName + "_1pions").c_str(), (fName + "_1pions" + fPlotTitles).c_str());
   SetAutoProcessTH1(onePions, kCMD_Reset, kCMD_Scale, kCMD_Norm);
 
   twoPions  = (TH1D*)(fDataHist->Clone());
   twoPions->SetNameTitle((fName + "_2pions").c_str(), (fName + "_2pions;" + fPlotTitles).c_str());
   SetAutoProcessTH1(twoPions, kCMD_Reset, kCMD_Scale, kCMD_Norm);
 
   threePions = (TH1D*)(fDataHist->Clone());
   threePions->SetNameTitle((fName + "_3pions").c_str(), (fName + "_3pions" + fPlotTitles).c_str());
   SetAutoProcessTH1(threePions, kCMD_Reset, kCMD_Scale, kCMD_Norm);
 
   morePions = (TH1D*)(fDataHist->Clone());
   morePions->SetNameTitle((fName + "_4pions").c_str(), (fName + "_4pions" + fPlotTitles).c_str());
   SetAutoProcessTH1(morePions, kCMD_Reset, kCMD_Scale, kCMD_Norm);
 
   // Final setup  ---------------------------------------------------
   FinaliseMeasurement();
 
 };
 
 //********************************************************************
 // Here we have to fill for every pion we find in the event
 void MINERvA_CCNpip_XSec_1DTpi_nu::FillEventVariables(FitEvent *event) {
   //********************************************************************
 
   if (event->NumFSParticle(211) == 0 && event->NumFSParticle(-211) == 0) return;
   if (event->NumFSParticle(13) == 0) return;
 
   // Need to make this use event boxes
 
   // Clear out the vectors
   GetPionBox()->Reset();
   TLorentzVector Pnu  = event->GetNeutrinoIn()->fP;
   TLorentzVector Pmu  = event->GetHMFSParticle(13)->fP;
 
   // Loop over the particle stack
   for (unsigned int j = 2; j < event->Npart(); ++j) {
     
     // Only include alive particles
     if (event->GetParticleState(j) != kFinalState) continue;
 
     int PID = (event->PartInfo(j))->fPID;
     // Pick up the charged pions in the event
     if (abs(PID) == 211) {
         double ppi = FitUtils::T(event->PartInfo(j)->fP) * 1000.;
 	GetPionBox()->fTpiVect.push_back(ppi);
     }
   }
   
 
   fXVar = 0;
 
   return;
 };
 
 //********************************************************************
 // The last bool refers to if we're using restricted phase space or not
 bool MINERvA_CCNpip_XSec_1DTpi_nu::isSignal(FitEvent *event) {
 //********************************************************************
   // Last false refers to that this is NOT the restricted MINERvA phase space,
   // in which only forward-going muons are accepted
   return SignalDef::isCCNpip_MINERvA(event, EnuMin, EnuMax, !fFullPhaseSpace, !fUpdatedData);
 }
 
 //********************************************************************
 // Need to override FillHistograms() here because we fill the histogram N_pion
 // times
 void MINERvA_CCNpip_XSec_1DTpi_nu::FillHistograms() {
 //********************************************************************
 
   if (Signal) {
-    unsigned int nPions = GetPionBox()->fTpiVect.size();
 
+    unsigned int nPions = GetPionBox()->fTpiVect.size();
     // Need to loop over all the pions in the sample
     for (size_t k = 0; k < nPions; ++k) {
       double tpi = GetPionBox()->fTpiVect[k];
       this->fMCHist->Fill(tpi, Weight);
       this->fMCFine->Fill(tpi, Weight);
       this->fMCStat->Fill(tpi, 1.0);
 
       if (nPions == 1) {
         onePions->Fill(tpi, Weight);
       } else if (nPions == 2) {
         twoPions->Fill(tpi, Weight);
       } else if (nPions == 3) {
         threePions->Fill(tpi, Weight);
       } else if (nPions > 3) {
         morePions->Fill(tpi, Weight);
       }
 
       if (fMCHist_Modes) fMCHist_Modes->Fill(Mode, tpi, Weight);
     }
   }
 }
 
 //********************************************************************
 void MINERvA_CCNpip_XSec_1DTpi_nu::ScaleEvents() {
   //********************************************************************
   Measurement1D::ScaleEvents();
 
   onePions->Scale(this->fScaleFactor, "width");
   twoPions->Scale(this->fScaleFactor, "width");
   threePions->Scale(this->fScaleFactor, "width");
   morePions->Scale(this->fScaleFactor, "width");
 
   return;
 }
 
 //********************************************************************
 void MINERvA_CCNpip_XSec_1DTpi_nu::Write(std::string drawOpts) {
 //********************************************************************
   Measurement1D::Write(drawOpts);
 
   // Make an auto processed pion stack
   // Draw the npions stack
   onePions->SetTitle("1#pi");
   onePions->SetLineColor(kBlack);
   // onePions->SetFillStyle(0);
   onePions->SetFillColor(onePions->GetLineColor());
 
   twoPions->SetTitle("2#pi");
   twoPions->SetLineColor(kRed);
   // twoPions->SetFillStyle(0);
   twoPions->SetFillColor(twoPions->GetLineColor());
 
   threePions->SetTitle("3#pi");
   threePions->SetLineColor(kGreen);
   // threePions->SetFillStyle(0);
   threePions->SetFillColor(threePions->GetLineColor());
 
   morePions->SetTitle(">3#pi");
   morePions->SetLineColor(kBlue);
   // morePions->SetFillStyle(0);
   morePions->SetFillColor(morePions->GetLineColor());
 
   THStack pionStack =
     THStack((fName + "_pionStack").c_str(), (fName + "_pionStack").c_str());
 
   pionStack.Add(onePions);
   pionStack.Add(twoPions);
   pionStack.Add(threePions);
   pionStack.Add(morePions);
 
   pionStack.Write();
 
   return;
 }