diff --git a/inc/LauCategoryFlavTag.hh b/inc/LauCategoryFlavTag.hh
index 33c1e40..3115928 100644
--- a/inc/LauCategoryFlavTag.hh
+++ b/inc/LauCategoryFlavTag.hh
@@ -1,185 +1,180 @@
 
 /*
 Copyright 2017 University of Warwick
 
 Licensed under the Apache License, Version 2.0 (the "License");
 you may not use this file except in compliance with the License.
 You may obtain a copy of the License at
 
     http://www.apache.org/licenses/LICENSE-2.0
 
 Unless required by applicable law or agreed to in writing, software
 distributed under the License is distributed on an "AS IS" BASIS,
 WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 See the License for the specific language governing permissions and
 limitations under the License.
 */
 
 /*
 Laura++ package authors:
 John Back
 Paul Harrison
 Thomas Latham
 */
 
 /*! \file LauCategoryFlavTag.hh
   \brief File containing declaration of LauCategoryFlavTag class.
  */
 
 /*! \class LauCategoryFlavTag
   \brief Class for defining the flavour tagging approach.
 
   Define the flavour tagging categories and all associated parameters
   to be passed to the relevant fit models.
  */
 
 #ifndef LAU_CATEGORYFLAVTAG
 #define LAU_CATEGORYFLAVTAG
 
 #include <vector>
 #include <map>
 #include <set>
-#include <iostream>
 
 #include "TString.h"
-#include "TStopwatch.h"
-#include "TSystem.h"
-#include "TH1D.h"
 
-#include "LauConstants.hh"
 #include "LauParameter.hh"
-#include "LauFitDataTree.hh"
-#include "LauAbsFitModel.hh"
-#include "LauDecayTimePdf.hh"
+
+class LauFitDataTree;
+
 
 class LauCategoryFlavTag {
 
 	public:
 		//! Constructor
 		/*!
 		  \param [in] modelB0bar DP model for the antiparticle
 		  \param [in] modelB0 DP model for the particle
 		  \param [in] useUntaggedEvents should the untagged sample be used or excluded?
 		  \param [in] tagVarName the variable name in the data tree that specifies the event tag
 		  \param [in] tagCatVarName the variable name in the data tree that specifies the event tagging category
 		 */
 		LauCategoryFlavTag(const std::vector<LauParameter*>& params, const Bool_t useUntaggedEvents = kTRUE, const TString& tagVarName = "tagFlv", const TString& tagCatVarName = "tagCat");
 
 		//! Destructor
 		virtual ~LauCategoryFlavTag();
 
 		void initialise();
 
 		//! Add a tagging category to the list of valid categories
 		/*!
 		  NB category 0 is always valid and corresponds to untagged events.
 		  Whether untagged events are used in the fit or note is controlled by a constructor argument.
 
 		  \param [in] tagCat the tagging category ID
 		 */
 		void addValidTagCat(const Int_t tagCat);
 
 		//! Add several tagging categories to the list of valid categories
 		/*!
 		  NB category 0 is always valid and corresponds to untagged events.
 		  Whether untagged events are used in the fit or note is controlled by a constructor argument.
 
 		  \param [in] tagCats the list of tagging category IDs
 		 */
 		void addValidTagCats(const std::vector<Int_t>& tagCats);
 
 		//! Check the validity of the given tagging category
 		/*!
 		  \param [in] tagCat the tagging category ID
 		 */
 		Bool_t validTagCat(const Int_t tagCat) const;
 
 		//! Return the set of valid tagging categories
 		std::set<Int_t> getValidTagCats();
 
 		//! Change the dilutions, delta dilutions and tagCatFrac for signal if needed
 		/*!
 		  \param [in] tagCat the tagging category to adjust
 		  \param [in] tagCatFrac the tagging category fraction
 		  \param [in] dilution the tagging category average dilution = (1 - 2 * avg_mistag_fraction)
 		  \param [in] deltaDilution the tagging category dilution difference TODO - check sign convention
 		  \param [in] fixTCFrac whether to fix or float the tagging category fraction
 		 */
 		void setSignalTagCatPars(const Int_t tagCat, const Double_t tagCatFrac, const Double_t dilution, const Double_t deltaDilution,
 				const Bool_t fixTCFrac = kTRUE);
 
 		//! Read in the input fit data variables, e.g. m13Sq and m23Sq
 		void cacheInputFitVars(LauFitDataTree* inputFitData);
 
 		typedef std::map< Int_t, LauParameter> LauTagCatParamMap;
 
 		LauTagCatParamMap getDilution(){return dilution_;};
 		LauTagCatParamMap getDeltaDilution(){return deltaDilution_;};
 		LauTagCatParamMap& getSignalTagCatFrac(){return signalTagCatFrac_;};
 
 		//! Update the fraction of the first tagging category for all species
 		void setFirstTagCatFractions();
 
 		Int_t getEvtTagFlvVals(UInt_t iEvt){return evtTagFlvVals_[iEvt];};
 		Int_t getEvtTagCatVals(UInt_t iEvt){return evtTagCatVals_[iEvt];};
 		Int_t getCurEvtTagCat(){return curEvtTagCat_;};
 
 		LauParameter* findParameter(const TString& parName);
 		
 	protected:
 		//! Check the signal tagging category fractions
 		void checkSignalTagCatFractions();
 
 		//! Check that the background tagging category fractions are all present and sum to unity
 		/*!
 		    \param [in] theMap the container of tagcat fractions
 		*/
 		Bool_t checkTagCatFracMap(const LauTagCatParamMap& theMap) const;
 
 		//! Calculates the fraction of the first tagging category based on the others
 		/*!
 		    \param [in] theMap the container of tagcat fractions
 		*/
 		void setFirstTagCatFrac(LauTagCatParamMap& theMap);
 
 	private:
 		//! Whether or not to use untagged events
                 const Bool_t useUntaggedEvents_;
 
 		//! Signal tagging category fractions
 		LauTagCatParamMap signalTagCatFrac_;
 
 		//! Flavour tagging variable name
 		TString tagVarName_;
 
 		//! Tagging category variable name
 		TString tagCatVarName_;
 
 		//! The allowed tagging categories
 		std::set<Int_t> validTagCats_;
 
 		//! Flavour tag for current event
 		Int_t curEvtTagFlv_;
 
 		//! Tagging category for current event
 		Int_t curEvtTagCat_;
 
 		//! Vector to store event flavour tags
 		std::vector<Int_t> evtTagFlvVals_;
 
 		//! Vector to store event tagging categories
 		std::vector<Int_t> evtTagCatVals_;
 
 		//! Average dilutions
 		LauTagCatParamMap dilution_;
 
 		//! Dilution differences
 		LauTagCatParamMap deltaDilution_;
 
 		//! Input parameters
 		std::vector<LauParameter*> params_;
 
 		ClassDef(LauCategoryFlavTag,0) // Flaour tagging set up
 };
 
 #endif
diff --git a/inc/LauFlavTag.hh b/inc/LauFlavTag.hh
index 902f874..47d6197 100644
--- a/inc/LauFlavTag.hh
+++ b/inc/LauFlavTag.hh
@@ -1,236 +1,229 @@
 
 /*
 Copyright 2017 University of Warwick
 
 Licensed under the Apache License, Version 2.0 (the "License");
 you may not use this file except in compliance with the License.
 You may obtain a copy of the License at
 
     http://www.apache.org/licenses/LICENSE-2.0
 
 Unless required by applicable law or agreed to in writing, software
 distributed under the License is distributed on an "AS IS" BASIS,
 WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 See the License for the specific language governing permissions and
 limitations under the License.
 */
 
 /*
 Laura++ package authors:
 John Back
 Paul Harrison
 Thomas Latham
 */
 
 /*! \file LauFlavTag.hh
   \brief File containing declaration of LauFlavTag class.
  */
 
 /*! \class LauFlavTag
   \brief Class for defining the flavour tagging approach.
 
   Define the flavour tagging categories and all associated parameters
   to be passed to the relevant fit models.
  */
 
 #ifndef LAU_FLAVTAG
 #define LAU_FLAVTAG
 
-// TODO - audit these includes, there seem to be a number that are not necessary
-
 #include <vector>
 #include <map>
-#include <set>
-#include <iostream>
 
 #include "TString.h"
-#include "TStopwatch.h"
-#include "TSystem.h"
 
-#include "LauConstants.hh"
 #include "LauParameter.hh"
-#include "LauFitDataTree.hh"
-#include "LauAbsFitModel.hh"
-#include "LauDecayTimePdf.hh"
-#include "LauAbsPdf.hh"
+
+class LauAbsPdf;
+class LauFitDataTree;
+
 
 class LauFlavTag final {
 
 	public:
 		//! Constructor
 		/*!
 		    \param [in] useAveDelta use average and delta variables for tagging calibration and efficiency
 		    \param [in] useEtaPrime use eta prime rather the eta as the mistag throughout
 		*/    
 		LauFlavTag(const Bool_t useAveDelta = kFALSE, const Bool_t useEtaPrime = kFALSE);
 
 		//! Initialise
 		// TODO is this needed? Commented for the moment (here and where called in LauTimeDepFitModel)
 		//void initialise();
 
 		// TODO - need to decide which functions need to be public (interface) and which should be private (implementation details)
 
 		//! Change the dilutions, delta dilutions and tagCatFrac for signal if needed
 		/*!
 		  \param [in] name the name of the tagger
 		  \param [in] tagVarName the tagging variable name of the tagger in the ntuple
 		  \param [in] mistagVarName the associated mistag variable name of the same tagger in the ntuple
 		  \param [in] etapdf the mistag distribution for the tagger
 		  \param [in] tagEff tagging efficiency - (particle, antiparticle) or (average, delta) depending on useAveDelta_ flag
 		  \param [in] calib_p0 calibration parameter p0 - (particle, antiparticle) or (average, delta) depending on useAveDelta_ flag
 		  \param [in] calib_p1 calibration parameter p1 - (particle, antiparticle) or (average, delta) depending on useAveDelta_ flag
 		 */
 		// Need to set remember the position in the vector using a map for later reference
 		//void addTagger(const TString& name, const TString& tagVarName, const TString& mistagVarName, LauAbsPdf* etapdf, 
 		//		const Double_t tagEff_b0=1.0, const Double_t calib_p0_b0=1.0, const Double_t calib_p1_b0=1.0, 
 		//		const Double_t tagEff_b0bar=-1.0, const Double_t calib_p0_b0bar=-1.0, const Double_t calib_p1_b0bar=-1.0);
 		void addTagger(const TString& name, const TString& tagVarName, const TString& mistagVarName, LauAbsPdf* etapdf, 
 				const std::pair <Double_t, Double_t> tagEff, const std::pair <Double_t, Double_t> calib_p0,
 				const std::pair <Double_t, Double_t> calib_p1);
 
 		//! Read in the input fit data variables, e.g. m13Sq and m23Sq
 		void cacheInputFitVars(LauFitDataTree* inputFitData);
 
 		void updateEventInfo(const ULong_t iEvt);
 
 		const std::vector<TString>& getTagVarNames() const {return tagVarName_;};
 		const std::vector<TString>& getMistagVarNames() const {return mistagVarName_;};
 		const TString& getTrueTagVarName() const {return trueTagVarName_;};
 		const TString& getDecayFlvVarName() const {return decayFlvVarName_;};
 
 		Int_t getCurEvtTrueTagFlv() const {return curEvtTrueTagFlv_;};
 		Int_t getCurEvtDecayFlv() const {return curEvtDecayFlv_;};
 		const std::vector<Int_t>& getCurEvtTagFlv() const {return curEvtTagFlv_;};
 		const std::vector<Double_t>& getCurEvtMistag() const {return curEvtMistag_;};
 
 		ULong_t getNTaggers() const {return tagVarName_.size();}
 
 		//! Get map of calibration parameters for each tagging category
 		std::vector<LauParameter*> getCalibP0B0(){return calib_p0_B0_;};
 		std::vector<LauParameter*> getCalibP0B0bar(){return calib_p0_B0bar_;};
 		std::vector<LauParameter*> getCalibP1B0(){return calib_p1_B0_;};
 		std::vector<LauParameter*> getCalibP1B0bar(){return calib_p1_B0bar_;};
 		
 		//! Get map of alternative calibration parameters for each tagging category
 		std::vector<LauParameter*> getCalibP0Ave(){return calib_p0_ave_;};
 		std::vector<LauParameter*> getCalibP0Delta(){return calib_p0_delta_;};
 		std::vector<LauParameter*> getCalibP1Ave(){return calib_p1_ave_;};
 		std::vector<LauParameter*> getCalibP1Delta(){return calib_p1_delta_;};
 
 		//! Get map of alternative calibration parameters for each tagging category
 		std::vector<LauParameter*> getTagEffAve(){return tagEff_ave_;};
 		std::vector<LauParameter*> getTagEffDelta(){return tagEff_delta_;};
 		std::vector<LauParameter*> getTagEffB0(){return tagEff_B0_;};
 		std::vector<LauParameter*> getTagEffB0bar(){return tagEff_B0bar_;};
 
 		const std::vector<Double_t>& getPerEvtAvgMistag() const {return perEvtAvgMistag_;};
 
 		Double_t getLittleOmega(const ULong_t position, const Int_t flag) const;
 		Double_t getCapitalOmega(const ULong_t position, const Int_t flag) const;
 		
 		Double_t getEtaGen(const ULong_t position);
 
 		//! Return the Boolean controlling if we use the alternative tagging calibration parameters
 		Bool_t getUseAveDelta() const {return useAveDelta_;};
 
 		void setTrueTagVarName(TString trueTagVarName);
 		void setDecayFlvVarName(TString decayFlvVarName);
 
 		//! Gaussian constraints for P0 parameters for a given tagger
 		/*!
 		    param [in] name name of the tagger
 		    param [in] constraint1 the (mean, sigma) for the particle or average parameter
 		    param [in] constraint2 the (mean, sigma) for the antiparticle or delta parameter
 		*/    
 		void addP0GaussianConstraints(const TString name, const std::pair <Double_t, Double_t> constraint1, const std::pair <Double_t, Double_t> constraint2);
 
 		//! Gaussian constraints for P1 parameters for a given tagger
 		/*!
 		    param [in] name name of the tagger
 		    param [in] constraint1 the (mean, sigma) for the particle or average parameter
 		    param [in] constraint2 the (mean, sigma) for the antiparticle or delta parameter
 		*/    
 		void addP1GaussianConstraints(const TString name, const std::pair <Double_t, Double_t> constraint1, const std::pair <Double_t, Double_t> constraint2);
 
 		//! Gaussian constraints for tagging efficiency parameters for a given tagger
 		/*!
 		    param [in] name name of the tagger
 		    param [in] constraint1 the (mean, sigma) for the particle or average parameter
 		    param [in] constraint2 the (mean, sigma) for the antiparticle or delta parameter
 		*/    
 		void addTagEffGaussianConstraints(const TString name, const std::pair <Double_t, Double_t> constraint1, const std::pair <Double_t, Double_t> constraint2);
 
 	private:
 		//! Map to link taggers to their vector position
 		std::map<TString, Int_t> taggerPosition_;
 
 		//! Flavour tagging variable name
 		std::vector<TString> tagVarName_;
 
 		//! Per event mistag variable name
 		std::vector<TString> mistagVarName_;
 
 		//! True tag variable name for normalisation decays
 		TString trueTagVarName_;
 
 		//! Decay flavour tag variable name for normalisation decays
 		TString decayFlvVarName_;
 
 		//! Vector of flavour tags for each event
 		std::vector< std::vector<Int_t> > evtTagFlv_;
 
 		//! Flavour tag for current event
 		std::vector<Int_t> curEvtTagFlv_;
 
 		//! Vector of mistags for each event
 		std::vector< std::vector<Double_t> > evtMistag_;
 
 		//! Per event mistag for current event
 		std::vector<Double_t> curEvtMistag_;
 
 		//! Vector of true tags for each event
 		std::vector< Int_t > evtTrueTagFlv_;
 
 		//! Vector of decay tags for each event
 		std::vector< Int_t > evtDecayFlv_;
 
 		//! True tag from normalisation mode for current event
 		Int_t curEvtTrueTagFlv_{0};
 
 		//! True tag from normalisation mode for current event
 		Int_t curEvtDecayFlv_{0};
 
 		//! Per-event average mistag value (eta hat)
 		std::vector<Double_t> perEvtAvgMistag_;
 
 		//! Calibration parameters
 		std::vector<LauParameter*> calib_p0_B0_;
 		std::vector<LauParameter*> calib_p0_B0bar_;
 		std::vector<LauParameter*> calib_p1_B0_;
 		std::vector<LauParameter*> calib_p1_B0bar_;
 
 		//! Alternative calibration parameters
 		std::vector<LauParameter*> calib_p0_ave_;
 		std::vector<LauParameter*> calib_p0_delta_;
 		std::vector<LauParameter*> calib_p1_ave_;
 		std::vector<LauParameter*> calib_p1_delta_;
 
 		//! Flag to use alternative calibration parameters
 		Bool_t useAveDelta_; 
 
 		//! Flag to use eta prime not eta for the mistag 
 		Bool_t useEtaPrime_; 
 
 		//! Tagging efficiency parameters
 		std::vector<LauParameter*> tagEff_B0_;
 		std::vector<LauParameter*> tagEff_B0bar_;
 		std::vector<LauParameter*> tagEff_ave_;
 		std::vector<LauParameter*> tagEff_delta_;
 
 		//! Eta PDFs
 		std::vector<LauAbsPdf*> etaPdfs_;
 
 		ClassDef(LauFlavTag,0) // Flavour tagging set up
 };
 
 #endif
diff --git a/src/LauCategoryFlavTag.cc b/src/LauCategoryFlavTag.cc
index ec17565..bb1fd32 100644
--- a/src/LauCategoryFlavTag.cc
+++ b/src/LauCategoryFlavTag.cc
@@ -1,301 +1,280 @@
 
 /*
 Copyright 2017 University of Warwick
 
 Licensed under the Apache License, Version 2.0 (the "License");
 you may not use this file except in compliance with the License.
 You may obtain a copy of the License at
 
     http://www.apache.org/licenses/LICENSE-2.0
 
 Unless required by applicable law or agreed to in writing, software
 distributed under the License is distributed on an "AS IS" BASIS,
 WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 See the License for the specific language governing permissions and
 limitations under the License.
 */
 
 /*
 Laura++ package authors:
 John Back
 Paul Harrison
 Thomas Latham
 */
 
 /*! \file LauCategoryFlavTag.cc
   \brief File containing implementation of LauCategoryFlavTag class.
  */
 
 #include <iostream>
-#include <iomanip>
-#include <fstream>
-#include <map>
+#include <set>
 #include <vector>
 
-#include "TFile.h"
-#include "TMinuit.h"
-#include "TRandom.h"
+#include "TString.h"
 #include "TSystem.h"
-#include "TVirtualFitter.h"
-#include "TH1D.h"
-
-#include "LauAbsBkgndDPModel.hh"
-#include "LauAbsCoeffSet.hh"
-#include "LauAbsPdf.hh"
-#include "LauAsymmCalc.hh"
-#include "LauComplex.hh"
-#include "LauConstants.hh"
-#include "LauDPPartialIntegralInfo.hh"
-#include "LauDaughters.hh"
-#include "LauDecayTimePdf.hh"
-#include "LauFitNtuple.hh"
-#include "LauGenNtuple.hh"
-#include "LauIsobarDynamics.hh"
-#include "LauKinematics.hh"
-#include "LauPrint.hh"
-#include "LauRandom.hh"
-#include "LauScfMap.hh"
+
 #include "LauCategoryFlavTag.hh"
-#include "Lau1DHistPdf.hh"
+#include "LauFitDataTree.hh"
+#include "LauParameter.hh"
 
 ClassImp(LauCategoryFlavTag)
 
 	LauCategoryFlavTag::LauCategoryFlavTag(const std::vector<LauParameter*>& params, const Bool_t useUntaggedEvents, const TString& tagVarName, const TString& tagCatVarName) :
 		useUntaggedEvents_(useUntaggedEvents),
 		signalTagCatFrac_(),
 		tagVarName_(tagVarName),
 		tagCatVarName_(tagCatVarName),
 		validTagCats_(),
 		curEvtTagFlv_(0),
 		curEvtTagCat_(0),
 		evtTagFlvVals_(0),
 		evtTagCatVals_(0),
 		dilution_(),
 		deltaDilution_(),
 		params_(params)
 {
 	// Add the untagged category as a valid category
 	this->addValidTagCat(0);
 
 	// Set the fraction, average dilution and dilution difference for the untagged category
 	this->setSignalTagCatPars(0, 1.0, 1.0, 0.0, kTRUE);
 }
 
 LauCategoryFlavTag::~LauCategoryFlavTag()
 {
 }
 
 void LauCategoryFlavTag::initialise()
 {
 	this->checkSignalTagCatFractions();
 }
 
 void LauCategoryFlavTag::addValidTagCats(const std::vector<Int_t>& tagCats)
 {
 	for (std::vector<Int_t>::const_iterator iter = tagCats.begin(); iter != tagCats.end(); ++iter) {
 		this->addValidTagCat(*iter);
 	}
 }
 
 void LauCategoryFlavTag::addValidTagCat(const Int_t tagCat)
 {
 	validTagCats_.insert(tagCat);
 }
 
 void LauCategoryFlavTag::setSignalTagCatPars(const Int_t tagCat, const Double_t tagCatFrac, const Double_t dilution, const Double_t deltaDilution, const Bool_t fixTCFrac)
 {
 	if (!this->validTagCat(tagCat)) {
 		std::cerr<<"ERROR in LauCategoryFlavTag::setSignalTagCatPars : Tagging category \""<<tagCat<<"\" not valid, not changing the parameters."<<std::endl;
 		return;
 	}
 
 	TString tagCatFracName("signalTagCatFrac");
 	tagCatFracName += tagCat;
 
 	signalTagCatFrac_[tagCat].name(tagCatFracName);
 	signalTagCatFrac_[tagCat].range(0.0,1.0);
 	signalTagCatFrac_[tagCat].value(tagCatFrac); signalTagCatFrac_[tagCat].initValue(tagCatFrac); signalTagCatFrac_[tagCat].genValue(tagCatFrac);
 	signalTagCatFrac_[tagCat].fixed(fixTCFrac);
 
 	TString dilutionName("dilution");
 	dilutionName += tagCat;
 
 	dilution_[tagCat].name(dilutionName);
 	dilution_[tagCat].range(0.0,1.0);
 	dilution_[tagCat].value(dilution); dilution_[tagCat].initValue(dilution); dilution_[tagCat].genValue(dilution);
 	dilution_[tagCat].fixed(kTRUE);
 
 	TString deltaDilutionName("deltaDilution");
 	deltaDilutionName += tagCat;
 
 	deltaDilution_[tagCat].name(deltaDilutionName);
 	deltaDilution_[tagCat].range(-2.0,2.0);
 	deltaDilution_[tagCat].value(deltaDilution); deltaDilution_[tagCat].initValue(deltaDilution); deltaDilution_[tagCat].genValue(deltaDilution);
 	deltaDilution_[tagCat].fixed(kTRUE);
 }
 
 void LauCategoryFlavTag::checkSignalTagCatFractions()
 {
 	Double_t totalTaggedFrac(0.0);
 	for (LauTagCatParamMap::const_iterator iter=signalTagCatFrac_.begin(); iter!=signalTagCatFrac_.end(); ++iter) {
 		if (iter->first != 0) {
 			const LauParameter& par = iter->second;
 			totalTaggedFrac += par.value();
 		}
 	}
 
 	if ( ((totalTaggedFrac < (1.0-1.0e-8))&&!useUntaggedEvents_) || (totalTaggedFrac > (1.0+1.0e-8)) ) {
 		std::cerr<<"WARNING in LauCategoryFlavTag::checkSignalTagCatFractions : Tagging category fractions add up to "<<totalTaggedFrac<<", not 1.0; normalizing them."<<std::endl;
 		for (LauTagCatParamMap::iterator iter=signalTagCatFrac_.begin(); iter!=signalTagCatFrac_.end(); ++iter) {
 			LauParameter& par = iter->second;
 			Double_t newVal = par.value() / totalTaggedFrac;
 			par.value(newVal); par.initValue(newVal); par.genValue(newVal);
 		}
 	} else if (useUntaggedEvents_) {
 		Double_t tagCatFrac = 1.0 - totalTaggedFrac;
 		TString tagCatFracName("signalTagCatFrac0");
 		signalTagCatFrac_[0].name(tagCatFracName);
 		signalTagCatFrac_[0].range(0.0,1.0);
 		signalTagCatFrac_[0].value(tagCatFrac); signalTagCatFrac_[0].initValue(tagCatFrac); signalTagCatFrac_[0].genValue(tagCatFrac);
 		signalTagCatFrac_[0].fixed(kTRUE);
 
 		TString dilutionName("dilution0");
 		dilution_[0].name(dilutionName);
 		dilution_[0].range(0.0,1.0);
 		dilution_[0].value(0.0); dilution_[0].initValue(0.0); dilution_[0].genValue(0.0);
 		dilution_[0].fixed(kTRUE);
 
 		TString deltaDilutionName("deltaDilution0");
 		deltaDilution_[0].name(deltaDilutionName);
 		deltaDilution_[0].range(-2.0,2.0);
 		deltaDilution_[0].value(0.0); deltaDilution_[0].initValue(0.0); deltaDilution_[0].genValue(0.0);
 		deltaDilution_[0].fixed(kTRUE);
 	}
 
 	for (LauTagCatParamMap::const_iterator iter=dilution_.begin(); iter!=dilution_.end(); ++iter) {
 		std::cout<<"INFO in LauCategoryFlavTag::checkSignalTagCatFractions : Setting dilution for tagging category "<<(*iter).first<<" to "<<(*iter).second<<std::endl;
 	}
 	for (LauTagCatParamMap::const_iterator iter=deltaDilution_.begin(); iter!=deltaDilution_.end(); ++iter) {
 		std::cout<<"INFO in LauCategoryFlavTag::checkSignalTagCatFractions : Setting Delta(dilution) for tagging category "<<(*iter).first<<" to "<<(*iter).second<<std::endl;
 	}
 	for (LauTagCatParamMap::const_iterator iter=signalTagCatFrac_.begin(); iter!=signalTagCatFrac_.end(); ++iter) {
 		std::cout<<"INFO in LauCategoryFlavTag::checkSignalTagCatFractions : Setting (signal) fraction for tagging category "<<(*iter).first<<" to "<<(*iter).second<<std::endl;
 	}
 }
 
 void LauCategoryFlavTag::setFirstTagCatFractions()
 {
 	this->setFirstTagCatFrac(signalTagCatFrac_);
 	// TODO - add background maps
 }
 
 void LauCategoryFlavTag::setFirstTagCatFrac(LauTagCatParamMap& theMap)
 {
 	Double_t firstCatFrac = 1.0;
 	Int_t firstCat(0);
 	for (LauTagCatParamMap::iterator iter = theMap.begin(); iter != theMap.end(); ++iter) {
 		if (iter == theMap.begin()) {
 			firstCat = iter->first;
 			continue;
 		}
 		LauParameter& par = iter->second;
 		firstCatFrac -= par.unblindValue();
 	}
 	theMap[firstCat].value(firstCatFrac);
 }
 
 Bool_t LauCategoryFlavTag::validTagCat(Int_t tagCat) const
 {
 	return (validTagCats_.find(tagCat) != validTagCats_.end());
 }
 
 std::set<Int_t> LauCategoryFlavTag::getValidTagCats()
 {
 	return validTagCats_;
 }
 
 Bool_t LauCategoryFlavTag::checkTagCatFracMap(const LauTagCatParamMap& theMap) const
 {
 	// TODO - this is for checking the the background maps are OK (so it is unused at the moment)
 
 	// First check that there is an entry for each tagging category.
 	// NB an entry won't have been added if it isn't a valid category
 	// so don't need to check for that here.
 	if (theMap.size() != signalTagCatFrac_.size()) {
 		std::cerr<<"ERROR in LauCategoryFlavTag::checkTagCatFracMap : Not all tagging categories present."<<std::endl;
 		return kFALSE;
 	}
 
 	// Now check that the fractions sum up to unity.
 	Double_t tagCatFracSum(0.0);
 	for (LauTagCatParamMap::const_iterator iter = theMap.begin(); iter != theMap.end(); ++iter) {
 		const LauParameter& par = (*iter).second;
 		tagCatFracSum += par.unblindValue();
 	}
 	if ((tagCatFracSum - 1.0) > 1E-10) {
 		std::cerr<<"ERROR in LauCategoryFlavTag::checkTagCatFracMap : Tagging category event fractions do not sum to unity."<<std::endl;
 		return kFALSE;
 	}
 
 	// If we've got to here then all is OK.
 	return kTRUE;
 }
 
 void LauCategoryFlavTag::cacheInputFitVars(LauFitDataTree* inputFitData)
 {
 	// Fill event by event data
 	// Start by caching the tagging information
 	evtTagCatVals_.clear();
 	evtTagFlvVals_.clear();
 	if ( ! inputFitData->haveBranch( tagCatVarName_ ) ) {
 		std::cerr << "ERROR in LauCategoryFlavTag::cacheInputFitVars : Input data does not contain branch \"" << tagCatVarName_ << "\"." << std::endl;
 		gSystem->Exit(EXIT_FAILURE);
 	}
 	if ( ! inputFitData->haveBranch( tagVarName_ ) ) {
 		std::cerr << "ERROR in LauCategoryFlavTag::cacheInputFitVars : Input data does not contain branch \"" << tagVarName_ << "\"." << std::endl;
 		gSystem->Exit(EXIT_FAILURE);
 	}
 
 	UInt_t nEvents = inputFitData->nEvents();
 	evtTagCatVals_.reserve( nEvents );
 	evtTagFlvVals_.reserve( nEvents );
 
 	LauFitData::const_iterator fitdata_iter;
 	for (UInt_t iEvt = 0; iEvt < nEvents; iEvt++) {
 		const LauFitData& dataValues = inputFitData->getData(iEvt);
 
 		fitdata_iter = dataValues.find( tagCatVarName_ );
 		curEvtTagCat_ = static_cast<Int_t>( fitdata_iter->second );
 		if ( ! this->validTagCat( curEvtTagCat_ ) ) {
 			std::cerr << "WARNING in LauCategoryFlavTag::cacheInputFitVars : Invalid tagging category " << curEvtTagCat_ << " for event " << iEvt << ", setting it to untagged" << std::endl;
 			curEvtTagCat_ = 0;
 		}
 		evtTagCatVals_.push_back( curEvtTagCat_ );
 
 		fitdata_iter = dataValues.find( tagVarName_ );
 		curEvtTagFlv_ = static_cast<Int_t>( fitdata_iter->second );
 		if ( TMath::Abs( curEvtTagFlv_ ) != 1 && curEvtTagCat_!=0) {
 			if ( curEvtTagFlv_ > 0 ) {
 				std::cerr << "WARNING in LauCategoryFlavTag::cacheInputFitVars : Invalid tagging output " << curEvtTagFlv_ << " for event " << iEvt << ", setting it to +1" << std::endl;
 				curEvtTagFlv_ = +1;
 			} else {
 				std::cerr << "WARNING in LauCategoryFlavTag::cacheInputFitVars : Invalid tagging output " << curEvtTagFlv_ << " for event " << iEvt << ", setting it to -1" << std::endl;
 				curEvtTagFlv_ = -1;
 			}
 		}
 		evtTagFlvVals_.push_back( curEvtTagFlv_ );
 
 	}
 }
 
 LauParameter* LauCategoryFlavTag::findParameter(const TString& parName)
 {
 	for ( std::vector<LauParameter*>::iterator iter = params_.begin(); iter != params_.end(); ++iter ) {
 		if ((*iter)->name().Contains(parName)) {
 			return (*iter);
 		}
 	}
 	std::cerr << "ERROR in LauCategoryFlavTag::findParameter : Parameter \"" << parName << "\" not found." << std::endl;
 	return 0;
 }
 
diff --git a/src/LauFlavTag.cc b/src/LauFlavTag.cc
index 65c1de3..1f25d70 100644
--- a/src/LauFlavTag.cc
+++ b/src/LauFlavTag.cc
@@ -1,467 +1,445 @@
 
 /*
 Copyright 2017 University of Warwick
 
 Licensed under the Apache License, Version 2.0 (the "License");
 you may not use this file except in compliance with the License.
 You may obtain a copy of the License at
 
     http://www.apache.org/licenses/LICENSE-2.0
 
 Unless required by applicable law or agreed to in writing, software
 distributed under the License is distributed on an "AS IS" BASIS,
 WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 See the License for the specific language governing permissions and
 limitations under the License.
 */
 
 /*
 Laura++ package authors:
 John Back
 Paul Harrison
 Thomas Latham
 */
 
 /*! \file LauFlavTag.cc
   \brief File containing implementation of LauFlavTag class.
  */
 
-// TODO - audit these includes, there seem to be a number that are not necessary
-
 #include <iostream>
-#include <iomanip>
-#include <fstream>
 #include <map>
 #include <vector>
 
-#include "TFile.h"
-#include "TMinuit.h"
-#include "TRandom.h"
+#include "TMath.h"
+#include "TString.h"
 #include "TSystem.h"
-#include "TVirtualFitter.h"
-#include "TH1D.h"
 
-#include "LauAbsBkgndDPModel.hh"
-#include "LauAbsCoeffSet.hh"
+#include "Lau1DHistPdf.hh"
 #include "LauAbsPdf.hh"
-#include "LauAsymmCalc.hh"
-#include "LauComplex.hh"
-#include "LauConstants.hh"
-#include "LauDPPartialIntegralInfo.hh"
-#include "LauDaughters.hh"
-#include "LauDecayTimePdf.hh"
-#include "LauFitNtuple.hh"
-#include "LauGenNtuple.hh"
-#include "LauIsobarDynamics.hh"
-#include "LauKinematics.hh"
-#include "LauPrint.hh"
-#include "LauRandom.hh"
-#include "LauScfMap.hh"
 #include "LauFlavTag.hh"
-#include "Lau1DHistPdf.hh"
 
 ClassImp(LauFlavTag)
 
 LauFlavTag::LauFlavTag(const Bool_t useAveDelta, const Bool_t useEtaPrime) : 
 	useAveDelta_(useAveDelta),
 	useEtaPrime_(useEtaPrime)
 {
 }
 
 void LauFlavTag::addTagger(const TString& name, const TString& tagVarName, const TString& mistagVarName, LauAbsPdf* etapdf, const std::pair <Double_t, Double_t> tagEff, const std::pair <Double_t, Double_t> calib_p0, const std::pair <Double_t, Double_t> calib_p1)
 {
 	// Find how many taggers have already been added
 	const ULong_t position { tagVarName_.size() };
 
 	// Update map to relate tagger name and position in the vectors 
 	taggerPosition_[name]=position;
 
 	// Fill vectors
 	tagVarName_.push_back(tagVarName);
 	mistagVarName_.push_back(mistagVarName);
 	if (etapdf){
 		etaPdfs_.push_back(etapdf);
 		Lau1DHistPdf* etahistpdf = dynamic_cast<Lau1DHistPdf*>(etapdf);
 		if (etahistpdf){
 			perEvtAvgMistag_.push_back(etahistpdf->getMean());
 		} else {
 			std::cerr << "WARNING  in LauFlavTag::addTagger : Couldn't determine average eta value from PDF. Setting it to 0.4." << std::endl;
 			perEvtAvgMistag_.push_back(0.4);
 		}
 	} else {
 		std::cerr << "ERROR in LauFlavTag::addTagger : Eta PDF pointer is NULL" << std::endl;
 		gSystem->Exit(EXIT_FAILURE);
 	}
 
 	//Use particle/antiparticle variables
 	if (!useAveDelta_){
 		
 		TString tagEff_b0Name("tagEff_b0_"+name);
 		TString tagEff_b0barName("tagEff_b0bar_"+name);
 		TString calib_p0_b0Name("calib_p0_b0_"+name);
 		TString calib_p0_b0barName("calib_p0_b0bar_"+name);
 		TString calib_p1_b0Name("calib_p1_b0_"+name);
 		TString calib_p1_b0barName("calib_p1_b0bar_"+name);
 
 		LauParameter* tageffb0 = new LauParameter(tagEff_b0Name,tagEff.first,0.0,1.0,kTRUE);
 		tagEff_B0_.push_back(tageffb0);
 		tagEff_B0_[position]->initValue(tagEff.first); tagEff_B0_[position]->genValue(tagEff.first); 
 		tagEff_B0_[position]->fixed(kTRUE); //Update once full code in place
 
 		LauParameter* calibp0b0 = new LauParameter(calib_p0_b0Name,calib_p0.first,-10.0,10.0,kTRUE);
 		calib_p0_B0_.push_back(calibp0b0);
 		calib_p0_B0_[position]->initValue(calib_p0.first); calib_p0_B0_[position]->genValue(calib_p0.first); 
 		calib_p0_B0_[position]->fixed(kTRUE); //Update once full code in place
 
 		LauParameter* calibp1b0 = new LauParameter(calib_p1_b0Name,calib_p1.first,0.0,1.5,kTRUE);
 		calib_p1_B0_.push_back(calibp1b0);
 		calib_p1_B0_[position]->initValue(calib_p1.first); calib_p1_B0_[position]->genValue(calib_p1.first); 
 		calib_p1_B0_[position]->fixed(kTRUE); //Update once full code in place
 		
 		if (tagEff.second==-1.0 && calib_p0.second==-1.0 && calib_p1.second==-1.0){
 			tagEff_B0bar_.push_back(tagEff_B0_[position]->createClone(tagEff_b0barName));
 			calib_p0_B0bar_.push_back(calib_p0_B0_[position]->createClone(calib_p0_b0barName));
 			calib_p1_B0bar_.push_back(calib_p1_B0_[position]->createClone(calib_p1_b0barName));
 		} else {
 			LauParameter* tageffb0bar = new LauParameter(tagEff_b0barName,tagEff.second,0.0,1.0,kTRUE);
 			tagEff_B0bar_.push_back(tageffb0bar);
 			tagEff_B0bar_[position]->initValue(tagEff.second); tagEff_B0bar_[position]->genValue(tagEff.second); 
 			tagEff_B0bar_[position]->fixed(kTRUE); //Update once full code in place
 
 			LauParameter* calibp0b0bar = new LauParameter(calib_p0_b0barName,calib_p0.second,-10.0,10.0,kTRUE);
 			calib_p0_B0bar_.push_back(calibp0b0bar);
 			calib_p0_B0bar_[position]->initValue(calib_p0.second); calib_p0_B0bar_[position]->genValue(calib_p0.second); 
 			calib_p0_B0bar_[position]->fixed(kTRUE); //Update once full code in place
 
 			LauParameter* calibp1b0bar = new LauParameter(calib_p1_b0barName,calib_p1.second,0.0,1.5,kTRUE);
 			calib_p1_B0bar_.push_back(calibp1b0bar);
 			calib_p1_B0bar_[position]->initValue(calib_p1.second); calib_p1_B0bar_[position]->genValue(calib_p1.second); 
 			calib_p1_B0bar_[position]->fixed(kTRUE); //Update once full code in place
 		}
 
 	} else { //Use average and delta variables
 		
 		TString tagEff_aveName("tagEff_ave_"+name);
 		TString tagEff_deltaName("tagEff_delta_"+name);
 		TString calib_p0_aveName("calib_p0_ave_"+name);
 		TString calib_p0_deltaName("calib_p0_delta_"+name);
 		TString calib_p1_aveName("calib_p1_ave_"+name);
 		TString calib_p1_deltaName("calib_p1_delta_"+name);
 
 		LauParameter* tageffave = new LauParameter(tagEff_aveName,tagEff.first,0.0,1.0,kTRUE);
 		tagEff_ave_.push_back(tageffave);
 		tagEff_ave_[position]->initValue(tagEff.first); tagEff_ave_[position]->genValue(tagEff.first); 
 		tagEff_ave_[position]->fixed(kTRUE); //Update once full code in place
 
 		LauParameter* calibp0ave = new LauParameter(calib_p0_aveName,calib_p0.first,-10.0,10.0,kTRUE);
 		calib_p0_ave_.push_back(calibp0ave);
 		calib_p0_ave_[position]->initValue(calib_p0.first); calib_p0_ave_[position]->genValue(calib_p0.first); 
 		calib_p0_ave_[position]->fixed(kTRUE); //Update once full code in place
 
 		LauParameter* calibp1ave = new LauParameter(calib_p1_aveName,calib_p1.first,0.0,1.5,kTRUE);
 		calib_p1_ave_.push_back(calibp1ave);
 		calib_p1_ave_[position]->initValue(calib_p1.first); calib_p1_ave_[position]->genValue(calib_p1.first); 
 		calib_p1_ave_[position]->fixed(kTRUE); //Update once full code in place
 	
 		LauParameter* tageffdelta = new LauParameter(tagEff_deltaName,tagEff.second,-1.0,1.0,kTRUE);
 		tagEff_delta_.push_back(tageffdelta);
 		tagEff_delta_[position]->initValue(tagEff.second); tagEff_delta_[position]->genValue(tagEff.second); 
 		tagEff_delta_[position]->fixed(kTRUE); //Update once full code in place
 
 		LauParameter* calibp0delta = new LauParameter(calib_p0_deltaName,calib_p0.second,-10.0,10.0,kTRUE);
 		calib_p0_delta_.push_back(calibp0delta);
 		calib_p0_delta_[position]->initValue(calib_p0.second); calib_p0_delta_[position]->genValue(calib_p0.second); 
 		calib_p0_delta_[position]->fixed(kTRUE); //Update once full code in place
 
 		LauParameter* calibp1delta = new LauParameter(calib_p1_deltaName,calib_p1.second,-10.0,10.0,kTRUE);
 		calib_p1_delta_.push_back(calibp1delta);
 		calib_p1_delta_[position]->initValue(calib_p1.second); calib_p1_delta_[position]->genValue(calib_p1.second); 
 		calib_p1_delta_[position]->fixed(kTRUE); //Update once full code in place
 
 	}
 
 	std::cout<<"INFO in LauFlavTag::addTagger : Added tagger with name "<< name << std::endl;
 }
 
 void LauFlavTag::cacheInputFitVars(LauFitDataTree* inputFitData)
 {
 	evtTagFlv_.clear();
 	evtMistag_.clear();
 	evtTrueTagFlv_.clear();
 	evtDecayFlv_.clear();
 
 	// Loop over the taggers to check the branches
 	for (ULong_t i=0; i < tagVarName_.size(); ++i){
 		if ( ! inputFitData->haveBranch( tagVarName_[i] ) ) {
 			std::cerr << "ERROR in LauFlavTag::cacheInputFitVars : Input data does not contain branch \"" << tagVarName_[i] << "\"." << std::endl;
 			gSystem->Exit(EXIT_FAILURE);
 		}
 		if ( ! inputFitData->haveBranch( mistagVarName_[i] ) ) {
 			std::cerr << "ERROR in LauFlavTag::cacheInputFitVars : Input data does not contain branch \"" << mistagVarName_[i] << "\"." << std::endl;
 			gSystem->Exit(EXIT_FAILURE);
 		}
 	}
 	if ( ! inputFitData->haveBranch( trueTagVarName_ ) ) {
 		std::cerr << "ERROR in LauFlavTag::cacheInputFitVars : Input data does not contain branch \"" << trueTagVarName_ << "\"." << std::endl;
 		gSystem->Exit(EXIT_FAILURE);
 	}
 
 	const ULong_t nEvents { inputFitData->nEvents() };
 
 	evtTagFlv_.reserve( nEvents );
 	evtMistag_.reserve( nEvents );
 	evtTrueTagFlv_.reserve( nEvents );
 	evtDecayFlv_.reserve( nEvents );
 
 	LauFitData::const_iterator fitdata_iter;
 	for (ULong_t iEvt = 0; iEvt < nEvents; iEvt++) {
 		const LauFitData& dataValues = inputFitData->getData(iEvt);
 
 		//Clear vectors
 		curEvtTagFlv_.clear();
 		curEvtMistag_.clear();
 
 		// For untagged events see if we have a truth tag for normalisation modes
 		curEvtTrueTagFlv_ = ( trueTagVarName_ != "" ) ? static_cast<Int_t>( dataValues.at( trueTagVarName_ ) ) : 0;
 		if ( curEvtTrueTagFlv_ > 1 ) {
 			std::cerr << "WARNING in LauFlavTag::cacheInputFitVars : Invalid true tagging output " << curEvtTrueTagFlv_ << " for event " << iEvt << ", setting it to +1" << std::endl;
 			curEvtTrueTagFlv_ = +1;
 		} else if ( curEvtTrueTagFlv_ <  -1 ){
 			std::cerr << "WARNING in LauFlavTag::cacheInputFitVars : Invalid true tagging output " << curEvtTrueTagFlv_ << " for event " << iEvt << ", setting it to -1" << std::endl;
 			curEvtTrueTagFlv_ = -1;
 		}
 		evtTrueTagFlv_.push_back(curEvtTrueTagFlv_);
 	
 		// Flavour at decay
 		curEvtDecayFlv_ = ( decayFlvVarName_ != "" ) ? static_cast<Int_t>( dataValues.at( decayFlvVarName_ ) ) : 0;
 		if ( curEvtDecayFlv_ > 1 ) {
 			std::cerr << "WARNING in LauFlavTag::cacheInputFitVars : Invalid decay tagging output " << curEvtDecayFlv_ << " for event " << iEvt << ", setting it to +1" << std::endl;
 			curEvtDecayFlv_ = +1;
 		} else if ( curEvtDecayFlv_ <  -1 ){
 			std::cerr << "WARNING in LauFlavTag::cacheInputFitVars : Invalid decay tagging output " << curEvtDecayFlv_ << " for event " << iEvt << ", setting it to -1" << std::endl;
 			curEvtDecayFlv_ = -1;
 		}
 		evtDecayFlv_.push_back(curEvtDecayFlv_);
 		
 		for (ULong_t i=0; i < tagVarName_.size(); ++i){
                 
                 	curEvtTagFlv_.push_back( static_cast<Int_t>( dataValues.at( tagVarName_[i] ) ) );
                 	if ( curEvtTagFlv_[i] > 1 ) {
                 	        std::cerr << "WARNING in LauFlavTag::cacheInputFitVars : Invalid tagging output " << curEvtTagFlv_[i] << " for event " << iEvt << ", setting it to +1" << std::endl;
                 	        curEvtTagFlv_[i] = +1;
                 	} else if ( curEvtTagFlv_[i] < -1 ) {
                 	        std::cerr << "WARNING in LauFlavTag::cacheInputFitVars : Invalid tagging output " << curEvtTagFlv_[i] << " for event " << iEvt << ", setting it to -1" << std::endl;
                 	        curEvtTagFlv_[i] = -1;
                 	}
 
 			curEvtMistag_.push_back( static_cast<Double_t>( dataValues.at( mistagVarName_[i] ) ) );
 			// Calibrated mistag > 0.5 is just a tag flip - handled automatically in getCapitalOmega function
 			if (curEvtMistag_[i] > 0.5){
 				std::cerr<<"WARNING in LauFlavTag::cacheInputFitVars : Mistag value "<<curEvtMistag_[i]<<" is larger than 0.5, setting to 0.5"<<std::endl;
 				curEvtMistag_[i] = 0.5;
 			}
 			if (curEvtMistag_[i] < 0.0){
 				std::cerr<<"WARNING in LauFlavTag::cacheInputFitVars : Mistag value "<<curEvtMistag_[i]<<" is less than 0.0, setting to 0.0"<<std::endl;
 				curEvtMistag_[i] = 0.0;
 			}
 		}
 		evtTagFlv_.push_back(curEvtTagFlv_);
 		evtMistag_.push_back(curEvtMistag_);
 	}
 }
 
 void LauFlavTag::updateEventInfo(const ULong_t iEvt)
 {
 	//Clear vector
 	curEvtTagFlv_.clear();
 	curEvtMistag_.clear();
 
 	//Assign current event variables
         curEvtTagFlv_ = evtTagFlv_[iEvt];
 	curEvtMistag_ = evtMistag_[iEvt];
 	curEvtTrueTagFlv_ = evtTrueTagFlv_[iEvt];
 	curEvtDecayFlv_ = evtDecayFlv_[iEvt];
 }
 
 Double_t LauFlavTag::getLittleOmega(const ULong_t position, const Int_t flag) const
 {
 	if (TMath::Abs(flag) != 1){
 		std::cerr << "ERROR in LauFlavTag::getLittleOmega : Invalid flag, you must request either omega (+1) or omega bar (-1) to be returned" << std::endl;
 		return 0.0;
 	}
 
 	Double_t calibp0(0.), calibp1(0.), calibp0bar(0.), calibp1bar(0.);
 
 	//If we are floating average omega and delta omega we need to use those parameters instead 
 	if (useAveDelta_){
 		calibp0 = calib_p0_ave_[position]->unblindValue() + 0.5*calib_p0_delta_[position]->unblindValue();
 		calibp0bar = calib_p0_ave_[position]->unblindValue() - 0.5*calib_p0_delta_[position]->unblindValue();
 		calibp1 = calib_p1_ave_[position]->unblindValue() + 0.5*calib_p1_delta_[position]->unblindValue();
 		calibp1bar = calib_p1_ave_[position]->unblindValue() - 0.5*calib_p1_delta_[position]->unblindValue();
 	} else {
 		calibp0 = calib_p0_B0_[position]->unblindValue();
 		calibp0bar = calib_p0_B0bar_[position]->unblindValue();
 		calibp1 = calib_p1_B0_[position]->unblindValue();
 		calibp1bar = calib_p1_B0bar_[position]->unblindValue();
 	}
 
 	if (flag == 1){
 		return calibp0 + calibp1 * (curEvtMistag_[position] - perEvtAvgMistag_[position]);
 	}
 	else{
 		return calibp0bar + calibp1bar * (curEvtMistag_[position] - perEvtAvgMistag_[position]);
 	}
 	std::cerr << "ERROR in LauFlavTag::getLittleOmega : Current event flavour tag not defined" << std::endl;
 	return 0.0;
 }
 
 Double_t LauFlavTag::getCapitalOmega(const ULong_t position, const Int_t flag) const
 {
 	if (TMath::Abs(flag) != 1){
 		std::cerr << "ERROR in LauFlavTag::getCapitalOmega : Invalid flag, you must request either Omega (+1) or Omega bar (-1) to be returned" << std::endl;
 		return 0.0;
 	}
 
 	//Delta functions to control which terms contribute
 	Int_t deltap1(0), deltam1(0), delta0(0);
 
 	if (curEvtTagFlv_[position] == -1){
 		deltam1 = 1;
 	} else if(curEvtTagFlv_[position] == 1){
 		deltap1 = 1;
 	} else{
 		delta0 = 1;
 	}
 	
 	//Efficiency
 	Double_t eff=0.0;
 	if (useAveDelta_){
 		if(flag==1){
 			eff = tagEff_ave_[position]->unblindValue() + 0.5*tagEff_delta_[position]->unblindValue();
 		} else {
 			eff = tagEff_ave_[position]->unblindValue() - 0.5*tagEff_delta_[position]->unblindValue();
 		}
 	}else{
 		if(flag==1){
 			eff = tagEff_B0_[position]->unblindValue();
 		}else{
 			eff = tagEff_B0bar_[position]->unblindValue();
 		}
 	}
 	
 	//Little omega
 	Double_t omega = this->getLittleOmega(position, flag);
 	Double_t omegaPrime(0.);
 	//Transform to omega prime - TODO isn't this the inverse, getLittleOmega is actually giving us omegaPrime and on the next line we convert back to omega?
 	if (useEtaPrime_){
 		omegaPrime = (1/(1+TMath::Exp(-1.0*omega)));
 	}else{
 		omegaPrime = omega;
 	}
 
 	//little omega must be between 0 and 1. Force this for now, if the fits keep getting stuck can look more closely at it.
 	if (omegaPrime < 0.0){
 		std::cerr << "WARNING in LauFlavTag::getCapitalOmega the value of little omega is less than 0, shifting to 0" << std::endl;
 		omegaPrime = 0.0;
 	}
 	if (omegaPrime > 1.0){
 		std::cerr << "WARNING in LauFlavTag::getCapitalOmega the value of little omega is greater than 1, shifting to 1" << std::endl;
 		omegaPrime = 1.0;
 	}
 
 	//eta PDF value
 	std::vector<Double_t> abs;
 	abs.push_back(curEvtMistag_[position]);
 	etaPdfs_[position]->calcLikelihoodInfo(abs);
 	const Double_t h { etaPdfs_[position]->getLikelihood() };
 	const Double_t u { 2.0 }; // the PDF value for a uniform PDF between 0.0 and 0.5
 
 	//Put it together
 	if (flag == 1){ //Particle
 		return (deltap1*eff*(1-omegaPrime) + deltam1*eff*omegaPrime)*h + delta0*(1-eff)*u; 
 	} else {
 		return (deltam1*eff*(1-omegaPrime) + deltap1*eff*omegaPrime)*h + delta0*(1-eff)*u;
 	}
 }
 
 Double_t LauFlavTag::getEtaGen(const ULong_t position)
 {
 	//Clear mistag vector for a new event
 	if(position==0){
 		curEvtMistag_.clear();
 	}
         LauFitData data { etaPdfs_[position]->generate(nullptr) }; //TODO Add DP dependence?
         Double_t etagen { data.at(etaPdfs_[position]->varName()) };
 	if (etagen > 0.5){etagen = 0.5;}
 	if (etagen < 0.0){etagen = 0.0;}
 	curEvtMistag_.push_back(etagen);
 	return etagen;
 }
 
 void LauFlavTag::setTrueTagVarName(TString trueTagVarName){
 	trueTagVarName_ = std::move(trueTagVarName);
 }
 
 void LauFlavTag::setDecayFlvVarName(TString decayFlvVarName){
 	decayFlvVarName_ = std::move(decayFlvVarName);
 }
 
 void LauFlavTag::addP0GaussianConstraints(TString name, std::pair <Double_t, Double_t> constraint1, std::pair <Double_t, Double_t> constraint2){
 	//Does key exist?
 	if (taggerPosition_.count(name)==0){
 		std::cerr << "ERROR in LauFlavTag::addP0GaussianConstraints : Tagger name not recognised please check your options" << std::endl;
 		std::cerr << "ERROR in LauFlavTag::addP0GaussianConstraints : Constraints have not been applied" << std::endl;
 		return;
 	}
 	//Find position in the vector from the tagger name
 	Double_t pos  = taggerPosition_.at(name);
 
 	if (!useAveDelta_){
 		calib_p0_B0_[pos]->addGaussianConstraint(constraint1.first,constraint1.second);
 		calib_p0_B0bar_[pos]->addGaussianConstraint(constraint2.first,constraint2.second);
 	}else{
 		calib_p0_ave_[pos]->addGaussianConstraint(constraint1.first,constraint1.second);
 		calib_p0_delta_[pos]->addGaussianConstraint(constraint2.first,constraint2.second);
 	}
 	std::cout << "INFO in LauFlavTag::addP0GaussianConstraints : Added Gaussian constraints for the P0 calibration parameters of tagger " << name << std::endl;
 }
 
 void LauFlavTag::addP1GaussianConstraints(TString name, std::pair <Double_t, Double_t> constraint1, std::pair <Double_t, Double_t> constraint2){
 	//Does key exist?
 	if (taggerPosition_.count(name)==0){
 		std::cerr << "ERROR in LauFlavTag::addP1GaussianConstraints : Tagger name not recognised please check your options" << std::endl;
 		std::cerr << "ERROR in LauFlavTag::addP1GaussianConstraints : Constraints have not been applied" << std::endl;
 		return;
 	}
 	//Find position in the vector from the tagger name
 	Double_t pos  = taggerPosition_.at(name);
 
 	if (!useAveDelta_){
 		calib_p1_B0_[pos]->addGaussianConstraint(constraint1.first,constraint1.second);
 		calib_p1_B0bar_[pos]->addGaussianConstraint(constraint2.first,constraint2.second);
 	}else{
 		calib_p1_ave_[pos]->addGaussianConstraint(constraint1.first,constraint1.second);
 		calib_p1_delta_[pos]->addGaussianConstraint(constraint2.first,constraint2.second);
 	}
 	std::cout << "INFO in LauFlavTag::addP1GaussianConstraints : Added Gaussian constraints for the P1 calibration parameters of tagger " << name << std::endl;
 }
 
 void LauFlavTag::addTagEffGaussianConstraints(TString name, std::pair <Double_t, Double_t> constraint1, std::pair <Double_t, Double_t> constraint2){
 	//Does key exist?
 	if (taggerPosition_.count(name)==0){
 		std::cerr << "ERROR in LauFlavTag::addTagEffGaussianConstraints : Tagger name not recognised please check your options" << std::endl;
 		std::cerr << "ERROR in LauFlavTag::addTagEffGaussianConstraints : Constraints have not been applied" << std::endl;
 		return;
 	}
 	//Find position in the vector from the tagger name
 	Double_t pos  = taggerPosition_.at(name);
 
 	if (!useAveDelta_){
 		tagEff_B0_[pos]->addGaussianConstraint(constraint1.first,constraint1.second);
 		tagEff_B0bar_[pos]->addGaussianConstraint(constraint2.first,constraint2.second);
 	}else{
 		tagEff_ave_[pos]->addGaussianConstraint(constraint1.first,constraint1.second);
 		tagEff_delta_[pos]->addGaussianConstraint(constraint2.first,constraint2.second);
 	}
 	std::cout << "INFO in LauFlavTag::addTagEffGaussianConstraints : Added Gaussian constraints for the tagging efficiency parameters of tagger " << name << std::endl;
 }