diff --git a/examples/GenFitTimeDep_Bd2Dpipi.cc b/examples/GenFitTimeDep_Bd2Dpipi.cc
index 2f84c37..c1c0592 100644
--- a/examples/GenFitTimeDep_Bd2Dpipi.cc
+++ b/examples/GenFitTimeDep_Bd2Dpipi.cc
@@ -1,1021 +1,1023 @@
 
-#include <algorithm>
-#include <filesystem>
-#include <iostream>
-#include <map>
-#include <tuple>
-#include <vector>
-
-#include "TFile.h"
-#include "TH2.h"
-#include "TRandom.h"
-#include "TString.h"
-#include "TSystem.h"
-#include "TF1.h"
-#include "TCanvas.h"
+#include "TimeDep_Dpipi_ProgOpts.hh"
 
 #include "Lau1DHistPdf.hh"
 #include "Lau1DCubicSpline.hh"
 #include "LauBinnedDecayTimeEfficiency.hh"
 #include "LauBkgndDPModel.hh"
 #include "LauCrystalBallPdf.hh"
 #include "LauDaughters.hh"
 #include "LauDecayTime.hh"
 #include "LauDecayTimePdf.hh"
 #include "LauDecayTimePhysicsModel.hh"
 #include "LauDecayTimeResolution.hh"
 #include "LauEffModel.hh"
 #include "LauExponentialPdf.hh"
 #include "LauFlavTag.hh"
 #include "LauGounarisSakuraiRes.hh"
 #include "LauIsobarDynamics.hh"
 #include "LauLHCbNR.hh"
 #include "LauMagPhaseCoeffSet.hh"
 #include "LauRandom.hh"
 #include "LauRealImagCoeffSet.hh"
 #include "LauResonanceMaker.hh"
 #include "LauSplineDecayTimeEfficiency.hh"
 #include "LauSumPdf.hh"
 #include "LauTimeDepFitModel.hh"
 #include "LauVetoes.hh"
 
-#include "TimeDep_Dpipi_ProgOpts.hh"
+#include "TFile.h"
+#include "TH2.h"
+#include "TRandom.h"
+#include "TString.h"
+#include "TSystem.h"
+#include "TF1.h"
+#include "TCanvas.h"
+
+#include "nlohmann/json.hpp"
+
+#include <algorithm>
+#include <filesystem>
+#include <iostream>
+#include <map>
+#include <tuple>
+#include <vector>
 
 // Helper function for creating Double-Crystal-Ball PDFs
 LauAbsPdf* makeDCBPDF( const TString& varName, const Double_t varMin, const Double_t varMax, const TString& componentName, const std::map<TString,Double_t> paramVals );
 
 // Helper function for creating sets of DTA splines
 // TODO - not yet used
 using LauCubicSplineDecayTimeEfficiencyPtr = std::unique_ptr<LauSplineDecayTimeEfficiency<LauSplineOrder::Cubic>>;
 using LauSixthSplineDecayTimeEfficiencyPtr = std::unique_ptr<LauSplineDecayTimeEfficiency<LauSplineOrder::Sixth>>;
 auto createSplineDTASet( const TString& jsonFileName, const TString& controlSplineName, const TString& signalSplineName, const std::vector<TString>& bkgndSplineNames )
 {
 	auto controlSpline = Lau1DCubicSpline::readFromJson( jsonFileName, controlSplineName );
 	LauCubicSplineDecayTimeEfficiencyPtr controlEff { std::make_unique<LauSplineDecayTimeEfficiency<LauSplineOrder::Cubic>>( controlSplineName, std::move( controlSpline ) ) };
 
 	LauSixthSplineDecayTimeEfficiencyPtr signalEff;
 	if ( signalSplineName != "" ) {
 		auto signalSpline = Lau1DCubicSpline::readFromJson( jsonFileName, signalSplineName );
 		signalEff = std::make_unique<LauSplineDecayTimeEfficiency<LauSplineOrder::Sixth>>( signalSplineName, *controlEff, std::move( signalSpline ) );
 	}
 
 	std::vector<LauSixthSplineDecayTimeEfficiencyPtr> bkgndEffs;
 	bkgndEffs.reserve( bkgndSplineNames.size() );
 	for ( auto& bkgndSplineName : bkgndSplineNames ) {
 		auto bkgndSpline = Lau1DCubicSpline::readFromJson( jsonFileName, bkgndSplineName );
 		bkgndEffs.push_back( std::make_unique<LauSplineDecayTimeEfficiency<LauSplineOrder::Sixth>>( bkgndSplineName, *controlEff, std::move( bkgndSpline ) ) );
 	}
 
 	return std::make_tuple( std::move(controlEff), std::move(signalEff), std::move(bkgndEffs) );
 }
 
 class BgTagDecayFlvAsym {
 	private:
 		const Double_t eff_decayB0_taggedB0       {0.};
 		const Double_t eff_decayB0_taggedB0bar    {0.};
 		const Double_t eff_decayB0bar_taggedB0    {0.};
 		const Double_t eff_decayB0bar_taggedB0bar {0.};
 		const Double_t eff_decayB0_ave            {0.};
 		const Double_t eff_decayB0bar_ave         {0.};
 		const Double_t eff_decayB0_delta          {0.};
 		const Double_t eff_decayB0bar_delta       {0.};
 
 	public:
 		BgTagDecayFlvAsym() = default;
 
 		BgTagDecayFlvAsym( //This constructor assumes a tagging efficincy of 100%
 			const Double_t n_decayB0_taggedB0,
 			const Double_t n_decayB0_taggedB0bar,
 			const Double_t n_decayB0bar_taggedB0,
 			const Double_t n_decayB0bar_taggedB0bar
 		):
 			eff_decayB0_taggedB0      {n_decayB0_taggedB0   /(n_decayB0_taggedB0+n_decayB0_taggedB0bar)},
 			eff_decayB0_taggedB0bar   {n_decayB0_taggedB0bar/(n_decayB0_taggedB0+n_decayB0_taggedB0bar)},
 			eff_decayB0bar_taggedB0   {n_decayB0bar_taggedB0   /(n_decayB0bar_taggedB0+n_decayB0bar_taggedB0bar)},
 			eff_decayB0bar_taggedB0bar{n_decayB0bar_taggedB0bar/(n_decayB0bar_taggedB0+n_decayB0bar_taggedB0bar)},
 
 			eff_decayB0_ave           {0.5*(eff_decayB0_taggedB0    + eff_decayB0_taggedB0bar)},
 			eff_decayB0bar_ave        {0.5*(eff_decayB0bar_taggedB0 + eff_decayB0bar_taggedB0bar)},
 			eff_decayB0_delta         {eff_decayB0_taggedB0    - eff_decayB0_taggedB0bar},
 			eff_decayB0bar_delta      {eff_decayB0bar_taggedB0 - eff_decayB0bar_taggedB0bar}
 		{}
 
+		BgTagDecayFlvAsym( const nlohmann::json& j ) : BgTagDecayFlvAsym( j[0], j[1], j[2], j[3] )
+		{}
+
 		std::pair<Double_t,Double_t> decayB0()       const {return std::make_pair(eff_decayB0_taggedB0   ,eff_decayB0_taggedB0bar   );}
 		std::pair<Double_t,Double_t> decayB0bar()    const {return std::make_pair(eff_decayB0bar_taggedB0,eff_decayB0bar_taggedB0bar);}
 		std::pair<Double_t,Double_t> aveDeltaB0()    const {return std::make_pair(eff_decayB0_ave   , eff_decayB0_delta)   ;}
 		std::pair<Double_t,Double_t> aveDeltaB0bar() const {return std::make_pair(eff_decayB0bar_ave, eff_decayB0bar_delta);}
 
 	//TODO add constructor for if tagging eff isn't 100%?
 };
 
 
 int main(const int argc, const  char ** argv)
 {
 	const TestDpipi_ProgramSettings settings{argc,argv};
 	if ( ! settings.parsedOK ) {
 		return EXIT_FAILURE;
 	}
 
 	if ( settings.helpRequested ) {
 		return EXIT_SUCCESS;
 	}
 
 	if(settings.dataFit and settings.command == Command::Generate)
 	{
 		std::cerr << "I can't generate data! Wait for Run 3" << std::endl;
 		return EXIT_FAILURE;
 	}
 	if(settings.dataFit and not settings.blindFit and not settings.fixPhiMix)
 	{
 		std::cerr << "We don't have permission to do unblinded fits yet!" << std::endl;
 		return EXIT_FAILURE;
 	}
 	if(settings.run == 0 or settings.run > 2)
 	{
 		std::cerr << "The Run number must be either 1 or 2" << std::endl;
 		return EXIT_FAILURE;
 	}
 
 	const UInt_t runNo = settings.run;
+	const std::string runStr = Form("Run%u", runNo);
+
 	LauRandom::setSeed(settings.RNGseed);
 
 	LauDaughters* daughtersB0bar = new LauDaughters("B0_bar", "pi+", "pi-", "D0",    kTRUE);
 	LauDaughters* daughtersB0    = new LauDaughters("B0",     "pi+", "pi-", "D0_bar",kTRUE);
 
 	// efficiency
 	auto vetoes { LauVetoes::readFromJson( "Bd2D0pipi_DP_Vetoes.json", "Bd2D0pipi_Vetoes" ) };
 	LauEffModel* effModelB0bar = new LauEffModel(daughtersB0bar, vetoes.get());
 	LauEffModel* effModelB0 = new LauEffModel(daughtersB0, vetoes.get());
 
 	TFile* effFile(nullptr);
 	TH2* effHist(nullptr);
 
 	if (settings.directory == "B2Dhh_Kpi"){
 		effFile = TFile::Open(Form("%s/%s/eff_map_mprime_thprime_model_d2kpi_run%u.root",settings.eosRoot.c_str(),settings.eosConfigDir.c_str(),runNo));
 		effHist = dynamic_cast<TH2*>(effFile->Get(Form("eff_map_mprime_thprime_model_d2kpi_run%u",runNo)));
 	} else if (settings.directory == "B2Dhh_KK"){
 		effFile = TFile::Open(Form("%s/%s/eff_map_mprime_thprime_model_d2kk_run%u.root",settings.eosRoot.c_str(),settings.eosConfigDir.c_str(),runNo));
 		effHist = dynamic_cast<TH2*>(effFile->Get(Form("eff_map_mprime_thprime_model_d2kk_run%u",runNo)));
 	} else if (settings.directory == "B2Dhh_pipi"){
 		effFile = TFile::Open(Form("%s/%s/eff_map_mprime_thprime_model_d2pipi_run%u.root",settings.eosRoot.c_str(),settings.eosConfigDir.c_str(),runNo));
 		effHist = dynamic_cast<TH2*>(effFile->Get(Form("eff_map_mprime_thprime_model_d2pipi_run%u",runNo)));
 	} else {std::cerr << "bad dir given!" << std::endl; return EXIT_FAILURE;}
 
 	effModelB0->setEffHisto(effHist, kTRUE, kFALSE, 0.0, 0.0, kFALSE, kTRUE);
 	effModelB0bar->setEffHisto(effHist, kTRUE, kFALSE, 0.0, 0.0, kFALSE, kTRUE);
 
 	// background types
 	std::map<TString,LauFlavTag::BkgndType> BkgndTypes {
 		{"comb",    LauFlavTag::BkgndType::Combinatorial},
 		{"Bs2DKpi", LauFlavTag::BkgndType::FlavourSpecific},
 		{"Bd2DKpi", LauFlavTag::BkgndType::FlavourSpecific},
 		{"Lb2Dppi", LauFlavTag::BkgndType::FlavourSpecific}
 	};
 	auto notInBkgndList = [&settings](const std::map<TString,LauFlavTag::BkgndType>::value_type& e) {
 		return ( settings.bkgndList.find( e.first.Data() ) == settings.bkgndList.end() );
 	};
 	std::vector<TString> toRemove;
 	for ( const auto& elem : BkgndTypes ) {
 		if ( notInBkgndList(elem) ) {
 			toRemove.push_back( elem.first );
 		}
 	}
 	for ( const auto& name : toRemove ) {
 		BkgndTypes.erase( name );
 	}
 
 	// setup flavour tagging
 	const Bool_t useAveDeltaCalibVals { settings.useAveDeltaCalibVals };
 	const Bool_t useEtaPrime { kFALSE };
 	//LauFlavTag* flavTag = new LauFlavTag(useAveDeltaCalibVals,useEtaPrime);
 	LauFlavTag* flavTag = new LauFlavTag(useAveDeltaCalibVals,useEtaPrime,BkgndTypes);
 	if (settings.dType == LauTimeDepFitModel::CPEigenvalue::QFS) {
 		flavTag->setDecayFlvVarName("decayFlv");
 	}
 	if(not settings.dataFit)
 	{
 		flavTag->setTrueTagVarName("trueTag");
 	}
 
 	TFile* etaFile = TFile::Open(Form("%s/%s/ETAhists.root",settings.eosRoot.c_str(),settings.eosConfigDir.c_str()));
 	TH1* etaHist{nullptr};
 	TH1* etaHistB0{nullptr};
 	TH1* etaHistB0bar{nullptr};
 
 	if (settings.directory == "B2Dhh_Kpi"){
 		etaHistB0 = dynamic_cast<TH1*>(etaFile->Get(Form("signalMC_Kpi_Run%u_etaHist_dcyB",runNo)));
 		etaHistB0bar = dynamic_cast<TH1*>(etaFile->Get(Form("signalMC_Kpi_Run%u_etaHist_dcyBbar",runNo)));
 		etaHist = dynamic_cast<TH1*>(etaHistB0->Clone(Form("signalMC_Kpi_Run%u_etaHist",runNo)));
 		etaHist->Add(etaHistB0bar);
 	} else if (settings.directory == "B2Dhh_KK"){
 		etaHist = dynamic_cast<TH1*>(etaFile->Get(Form("signalMC_KK_Run%u_etaHist",runNo)));
 	} else if (settings.directory == "B2Dhh_pipi"){
 		etaHist = dynamic_cast<TH1*>(etaFile->Get(Form("signalMC_pipi_Run%u_etaHist",runNo)));
 	} else {std::cerr << "bad dir given!" << std::endl; return EXIT_FAILURE;}
 
 	// Crude check as to whether we're doing perfect vs realistic mis-tag
 	// - in the former case all entries should be in the first bin
 	// If the tagging is perfect then don't interpolate the eta histogram
 	// and also make it perfectly efficient, otherwise do interpolate and
 	// make it 50% efficient
 	Lau1DHistPdf* etaHistPdf = new Lau1DHistPdf( "eta", etaHist, 0.0, 0.5, kTRUE, kFALSE );
 	const Double_t meanEta { etaHistPdf->getMean() };
 
 	const Double_t tagEffVal{1.0};
 	std::pair<Double_t,Double_t> tagEff {tagEffVal, useAveDeltaCalibVals ? 0.0 : tagEffVal};
 
 	// use a null calibration for the time being, so p0 = <eta> and p1 = 1
 	std::pair<Double_t,Double_t> calib0 {meanEta, useAveDeltaCalibVals ? 0.0 : meanEta};
 	std::pair<Double_t,Double_t> calib1 {1.0, useAveDeltaCalibVals ? 0.0 : 1.0};
 
 	if(settings.addTaggers) {
 		flavTag->addTagger(Form("IFT_Run%u",runNo), "IFT_TAGDEC", "IFT_TAGETA", etaHistPdf, tagEff, calib0, calib1);
 		if (settings.floatCalibPars) {
 			flavTag->floatAllCalibPars();
 		}
 	}
 
 	// flavour tagging setup for backgrounds
 	std::map< TString, std::pair<Lau1DHistPdf*,Lau1DHistPdf*> > BkgndEtas;
 	TH1* bkgndEtaHist1{nullptr};
 	Lau1DHistPdf* bkgndEtaHistPDF1{nullptr};
 	TH1* bkgndEtaHistMinus1{nullptr};
 	Lau1DHistPdf* bkgndEtaHistPDFMinus1{nullptr};
 
 	for(auto& [name, _] : BkgndTypes)
 	{
 		const char* decayName = name == "comb" ? "sideband" : name.Data();
 		if (settings.directory == "B2Dhh_Kpi"){
 			bkgndEtaHist1      = dynamic_cast<TH1*>(etaFile->Get( Form("%s_Kpi_Run%u_etaHist_dcyB" ,decayName,runNo) ));
 			bkgndEtaHistPDF1      = new Lau1DHistPdf("eta",bkgndEtaHist1     ,0.0,0.5,kTRUE,kFALSE);
 			bkgndEtaHistMinus1 = dynamic_cast<TH1*>(etaFile->Get( Form("%s_Kpi_Run%u_etaHist_dcyBbar",decayName,runNo) ));
 			bkgndEtaHistPDFMinus1 = new Lau1DHistPdf("eta",bkgndEtaHistMinus1,0.0,0.5,kTRUE,kFALSE);
 		} else if (settings.directory == "B2Dhh_KK"){
 			bkgndEtaHist1      = dynamic_cast<TH1*>(etaFile->Get( Form("%s_KK_Run%u_etaHist" ,decayName,runNo) ));
 			bkgndEtaHistPDF1      = new Lau1DHistPdf("eta",bkgndEtaHist1     ,0.0,0.5,kTRUE,kFALSE);
 		} else if (settings.directory == "B2Dhh_pipi"){
 			bkgndEtaHist1      = dynamic_cast<TH1*>(etaFile->Get( Form("%s_pipi_Run%u_etaHist" ,decayName,runNo) ));
 			bkgndEtaHistPDF1      = new Lau1DHistPdf("eta",bkgndEtaHist1     ,0.0,0.5,kTRUE,kFALSE);
 		} else {std::cerr << "bad dir given!" << std::endl; return EXIT_FAILURE;}
 
 		BkgndEtas.emplace( std::make_pair(name, std::make_pair(bkgndEtaHistPDF1,bkgndEtaHistPDFMinus1) ) );
 	}
 
+	// NB deliberately not using uniform initialisation here because of this issue:
+	// https://json.nlohmann.me/home/faq/#brace-initialization-yields-arrays
+	const TString combTagAsymJsonFileName { "Bd2D0pipi_CombTagAsym.json" };
+	if ( not std::filesystem::exists(combTagAsymJsonFileName.Data()) ) {
+		std::cerr << "Warning : couldn't find the json file of combinatorial tagging asymmetry parameters; attempting to fetch from eos ..." << std::endl;
+		int r = system( Form("xrdcp %s/%s/%s .", settings.eosRoot.c_str(), settings.eosConfigDir.c_str(), combTagAsymJsonFileName.Data()) );
+		if (r != 0) { std::cerr << "FATAL : couldn't get the json file from remote!" << std::endl; return 1;}
+	}
+	const nlohmann::json combTagAsymJson = LauJsonTools::readJsonFile( combTagAsymJsonFileName.Data(), "", LauJsonTools::JsonType::Object );
+	if ( combTagAsymJson.is_null() ) {
+		std::cerr << "ERROR : unable to retrieve JSON object from root element of file \"" << combTagAsymJsonFileName << "\"" << std::endl;
+		return 1;
+	}
+
 	for(auto& [name,hists] : BkgndEtas)
 	{
-		if (settings.dType == LauTimeDepFitModel::CPEigenvalue::QFS){
-			if(name == "comb")
-			{
-				const BgTagDecayFlvAsym bgt =
-					runNo == 1 ? BgTagDecayFlvAsym(535.,576.,456.,588.) : BgTagDecayFlvAsym(1514.,1497.,1220.,1575.);
-				if(useAveDeltaCalibVals)
-				{
-					flavTag->setBkgndParams(name, Form("IFT_Run%u",runNo),
-							hists.first  , bgt.aveDeltaB0(),
-							hists.second , bgt.aveDeltaB0bar() );
-				}
-				else
-				{
-					flavTag->setBkgndParams(name, Form("IFT_Run%u",runNo),
-							hists.first  , bgt.decayB0(),
-							hists.second , bgt.decayB0bar());
-				}
+		if ( settings.dType == LauTimeDepFitModel::CPEigenvalue::QFS ) {
+			if ( name == "comb" ) {
+				const BgTagDecayFlvAsym bgt { combTagAsymJson.at(runStr).at(settings.directory) };
+				const std::pair<Double_t,Double_t> effsB0    { useAveDeltaCalibVals ? bgt.aveDeltaB0()    : bgt.decayB0()    };
+				const std::pair<Double_t,Double_t> effsB0bar { useAveDeltaCalibVals ? bgt.aveDeltaB0bar() : bgt.decayB0bar() };
+				flavTag->setBkgndParams( name, Form("IFT_Run%u",runNo),
+							 hists.first,  effsB0,
+							 hists.second, effsB0bar );
+			} else {
+				flavTag->setBkgndParams( name, Form("IFT_Run%u",runNo),
+							 hists.first,  tagEff,
+							 hists.second, tagEff );
 			}
-			else{flavTag->setBkgndParams(name, Form("IFT_Run%u",runNo), hists.first, tagEff, hists.second, tagEff );}
-			//flavTag->setBkgndParams(name, Form("IFT_Run%u",runNo), hist, name == "comb" ? std::make_pair(0.5,0.038) : tagEff );
 		} else if (settings.dType == LauTimeDepFitModel::CPEigenvalue::CPEven){
-			if(name == "comb")
-			{
-				BgTagDecayFlvAsym* bgt = nullptr;
-				if(settings.directory == "B2Dhh_KK")
-				{
-					bgt = runNo == 1 ? new BgTagDecayFlvAsym(344.,442.,344.,442.) : new BgTagDecayFlvAsym(710.,814.,710.,814.);
-				}
-				else
-				{
-					bgt = runNo == 1 ? new BgTagDecayFlvAsym(55.,68.,55.,68.)     : new BgTagDecayFlvAsym(210.,241.,210.,241.);
-				}
-
-				std::pair<Double_t,Double_t> effs{bgt->decayB0()};
-				if(useAveDeltaCalibVals){effs = bgt->aveDeltaB0();}
-				flavTag->setBkgndParams(name, Form("IFT_Run%u",runNo), hists.first, effs);
-				delete bgt;
+			if ( name == "comb" ) {
+				const BgTagDecayFlvAsym bgt { combTagAsymJson.at(runStr).at(settings.directory) };
+				const std::pair<Double_t,Double_t> effs { useAveDeltaCalibVals ? bgt.aveDeltaB0() : bgt.decayB0() };
+				flavTag->setBkgndParams( name, Form("IFT_Run%u",runNo), hists.first, effs);
+			} else {
+				flavTag->setBkgndParams( name, Form("IFT_Run%u",runNo), hists.first, tagEff );
 			}
-			else{flavTag->setBkgndParams(name, Form("IFT_Run%u",runNo), hists.first, tagEff);}
-//			else{flavTag->setBkgndParams(name, Form("IFT_Run%u",runNo), hists.first, tagEff, hists.second, tagEff );} FIXME once we have hists split by decay flv, uncomment this line and remove the above... actually don't FIXME, this is fine as above since the data isn't flavour specific in the CP modes. Keep it combined
 		}
 	}
 
 	// signal dynamics
 	if ( not std::filesystem::exists("Bd2D0pipi_DP_Model_Coeffs.json") ) {
 		std::cerr << "Warning : couldn't find the json file of DP model coeffs; attempting to fetch from eos ..." << std::endl;
 		int r = system( Form("xrdcp %s/%s/Bd2D0pipi_DP_Model_Coeffs.json .", settings.eosRoot.c_str(), settings.eosConfigDir.c_str()) );
 		if (r != 0) { std::cerr << "FATAL : couldn't get the json file from remote!" << std::endl; return 1;}
 	}
 
 	LauIsobarDynamics* sigModelB0bar = new LauIsobarDynamics(daughtersB0bar, effModelB0bar);
 	if (settings.directory == "B2Dhh_Kpi"){
 		sigModelB0bar->setIntFileName(Form("integ_Kpi_Run%u_B0bar.dat",runNo));
 	} else if (settings.directory == "B2Dhh_KK"){
 		sigModelB0bar->setIntFileName(Form("integ_KK_Run%u_B0bar.dat",runNo));
 	} else if (settings.directory == "B2Dhh_pipi"){
 		sigModelB0bar->setIntFileName(Form("integ_pipi_Run%u_B0bar.dat",runNo));
 	}
 	sigModelB0bar->constructModelFromJson( "Bd2D0pipi_DP_Model_Coeffs.json", "Bd2D0pipi_B0bar_Model" );
 
 	LauIsobarDynamics* sigModelB0 = new LauIsobarDynamics(daughtersB0, effModelB0);
 	if (settings.directory == "B2Dhh_Kpi"){
 		sigModelB0->setIntFileName(Form("integ_Kpi_Run%u_B0.dat",runNo));
 	} else if (settings.directory == "B2Dhh_KK"){
 		sigModelB0->setIntFileName(Form("integ_KK_Run%u_B0.dat",runNo));
 	} else if (settings.directory == "B2Dhh_pipi"){
 		sigModelB0->setIntFileName(Form("integ_pipi_Run%u_B0.dat",runNo));
 	}
 	sigModelB0->constructModelFromJson( "Bd2D0pipi_DP_Model_Coeffs.json", "Bd2D0pipi_B0_Model" );
 
 	// fit model
 	LauTimeDepFitModel* fitModel = new LauTimeDepFitModel(sigModelB0bar,sigModelB0,flavTag);
 	fitModel->breakOnBadLikelihood(kFALSE);
 	fitModel->sendFixedYieldsToFitter(kTRUE);
 
 	std::vector<std::unique_ptr<LauAbsCoeffSet>> coeffset { LauAbsCoeffSet::readFromJson( "Bd2D0pipi_DP_Model_Coeffs.json", "Bd2D0pipi_Coeffs" ) };
 
 	for ( auto& coeff : coeffset ) {
 		fitModel->setAmpCoeffSet( std::move(coeff) );
 	}
 
 	// background DP models
 	LauBkgndDPModel* BkgndModel{nullptr};
 	LauBkgndDPModel* BkgndModelbar{nullptr};
 	TFile* BkgndHistFile = TFile::Open(Form("%s/%s/BG_SDPs.root",settings.eosRoot.c_str(),settings.eosConfigDir.c_str()));
 	TFile* BkgndHistFile2 = TFile::Open(Form("%s/%s/sidebandSDP.root",settings.eosRoot.c_str(),settings.eosConfigDir.c_str()));
 	TH2* BkgndHistB0{nullptr};
 	TH2* BkgndHistB0bar{nullptr};
 	for(auto& [name, _] : BkgndTypes)
 	{
 		//Only have Kpi histos for now for misid
 		if( name=="comb"){
 			if( settings.directory=="B2Dhh_Kpi"){
 				BkgndHistB0 = dynamic_cast<TH2*>(BkgndHistFile2->   Get(Form("sideband_Run%u_Kpi_SDP_decayFlv1",runNo)));
 				BkgndHistB0bar = dynamic_cast<TH2*>(BkgndHistFile2->Get(Form("sideband_Run%u_Kpi_SDP_decayFlvMinus1",runNo)));
 			} else if( settings.directory=="B2Dhh_KK"){
 				BkgndHistB0 = dynamic_cast<TH2*>(BkgndHistFile2->Get(Form("sideband_Run%u_KK_SDP",runNo)));
 			} else if( settings.directory=="B2Dhh_pipi"){
 				BkgndHistB0 = dynamic_cast<TH2*>(BkgndHistFile2->Get(Form("sideband_Run%u_pipi_SDP",runNo)));
 			}
 		} else {
 			TString mode = settings.directory;
 			mode.ReplaceAll("B2Dhh_","");
 			TString histName = Form("%s_%s_Run%u_SDP_decFlv",name.Data(),mode.Data(),runNo);
 			BkgndHistB0 = dynamic_cast<TH2*>(BkgndHistFile->Get(histName+"1"));
 			BkgndHistB0bar = dynamic_cast<TH2*>(BkgndHistFile->Get(histName+"Minus1"));
 		}
 
 		BkgndModel = new LauBkgndDPModel(daughtersB0, vetoes.get());
 		BkgndHistB0->Smooth();
 		BkgndModel->setBkgndSpline(BkgndHistB0, kFALSE, kFALSE, kTRUE );
 
 		//For CP comb we don't split based on decay flavour since it's impossible
 		if(name == "comb" and settings.directory!="B2Dhh_Kpi")
 		{BkgndModelbar=nullptr;}
 		else
 		{
 			BkgndModelbar = new LauBkgndDPModel(daughtersB0bar, vetoes.get());
 			BkgndHistB0bar->Smooth();
 			BkgndModelbar->setBkgndSpline(BkgndHistB0bar, kFALSE, kFALSE, kTRUE );
 		}
 		fitModel->setBkgndDPModels( name, BkgndModel, BkgndModelbar );
 	}
 
 	// decay type and mixing parameter
 	const Bool_t fixPhiMix{ settings.fixPhiMix || settings.dType == LauTimeDepFitModel::CPEigenvalue::QFS };
 	const Bool_t useSinCos{ settings.useSinCos };
 
 	fitModel->setCPEigenvalue( settings.dType );
 	fitModel->setPhiMix( 2.0*LauConstants::beta, fixPhiMix, useSinCos );
 
 	if ( useSinCos ) {
 		if ( settings.blindFit ) {
 			//Blinding strings generated from: https://www.thewordfinder.com/random-sentence-generator/
 			const TString sinBlindString{settings.dataFit ? "Tom realized he could be making a big mistake." : "Don't Look a Gift Horse In The Mouth"};
 			const TString cosBlindString{settings.dataFit ? "He admitted defeat." : "Long In The Tooth"};
 			fitModel->blindPhiMix( sinBlindString, cosBlindString );
 		}
 		std::cout << "Is sin(phiMix) blinded? " << std::boolalpha << fitModel->sinPhiMixBlinded() << std::endl;
 		std::cout << "Is cos(phiMix) blinded? " << std::boolalpha << fitModel->cosPhiMixBlinded() << std::endl;
 	} else {
 		if ( settings.blindFit ) {
 			const TString blindString{"A Cut Above"};
 			fitModel->blindPhiMix( blindString );
 		}
 		std::cout << "Is phiMix blinded? " << std::boolalpha << fitModel->phiMixBlinded() << std::endl;
 	}
 
 	// production asymmetries
 	fitModel->setAsymmetries( 0.0, kTRUE );
 	for(auto& [name, _] : BkgndTypes) {
 		fitModel->setBkgndAsymmetries( name, 0.0, kTRUE );
 	}
 
 	// decay time PDFs
 	const Double_t minDt{0.2};
 	const Double_t maxDt{15.0};
 	const Double_t minDtErr{0.0};
 	const Double_t maxDtErr{0.15};
 
 
 	// physics models
 	auto dtPhysicsModels { LauDecayTimePhysicsModel::readFromJson( "Bd2D0pipi_DecayTime.json", "PhysicsModels" ) };
 
 	auto& dtPhysModel { dtPhysicsModels.at("signal") };
 	auto lifetimePar = dynamic_cast<LauParameter*>( dtPhysModel->findParameter("dt_sig_tau_Bd") );
 	if ( lifetimePar ) {
 		lifetimePar->fixed( settings.fixLifetime );
 	}
 	auto deltaMPar = dynamic_cast<LauParameter*>( dtPhysModel->findParameter("dt_sig_deltaM_Bd") );
 	if ( deltaMPar ) {
 		deltaMPar->fixed( settings.fixDeltaM );
 	}
 
 	// resolution models
 	auto dtResolutionModels { LauDecayTimeResolution::readFromJson( "Bd2D0pipi_DecayTime.json", Form("ResolutionModels_Run%d", runNo) ) };
 
 	auto& dtResoModel { dtResolutionModels.at("signal") };
 
 	// Decay time error histogram
 	// (always set this so that it gets generated properly, whether we're using it in the PDF or not)
 	TFile* dtErrFile = TFile::Open(Form("%s/%s/dte-hist.root",settings.eosRoot.c_str(),settings.eosConfigDir.c_str()));
 	TH1* dtErrHist = dynamic_cast<TH1*>(dtErrFile->Get("dte_hist"));
 
 	LauDecayTimePdf * dtPdf{nullptr};
 	if ( settings.timeResolution ) {
 		if ( settings.perEventTimeErr ) {
 			std::cerr << "Warning : the resolution model has not been constructed with per-event error scaling, but you have requested it" << std::endl;
 			dtPdf = new LauDecayTimePdf( "decayTime", minDt, maxDt, "decayTimeErr", minDtErr, maxDtErr, std::move(dtPhysModel), std::move(dtResoModel), dtErrHist );
 		} else {
 			dtPdf = new LauDecayTimePdf( "decayTime", minDt, maxDt, std::move(dtPhysModel), std::move(dtResoModel) );
 		}
 	} else {
 		dtPdf = new LauDecayTimePdf( "decayTime", minDt, maxDt, std::move(dtPhysModel) );
 	}
 
 	// Decay time acceptance histogram
 	TFile* dtaFile = TFile::Open(Form("%s/%s/DTAhists.root",settings.eosRoot.c_str(),settings.eosConfigDir.c_str()));
 	TH1* dtaHist{nullptr};
 
 	if (settings.directory == "B2Dhh_Kpi"){
 		dtaHist = dynamic_cast<TH1*>(dtaFile->Get(Form("signalMC_Kpi_Run%u_DTAhist",runNo)));
 	} else if (settings.directory == "B2Dhh_KK"){
 		dtaHist = dynamic_cast<TH1*>(dtaFile->Get(Form("signalMC_KK_Run%u_DTAhist",runNo)));
 	} else if (settings.directory == "B2Dhh_pipi"){
 		dtaHist = dynamic_cast<TH1*>(dtaFile->Get(Form("signalMC_pipi_Run%u_DTAhist",runNo)));
 	} else {std::cerr << "bad dir given!" << std::endl; return EXIT_FAILURE;}
 
 	if( not std::filesystem::exists("run1DTAsplines.json") and runNo == 1 )
 	{
 		std::cerr << "Warning : couldn't find the json file of splines in Run 1; attempting to fetch from eos ..." << std::endl;
 		int r = system( Form("xrdcp %s/%s/run1DTAsplines.json .", settings.eosRoot.c_str(), settings.eosConfigDir.c_str()) );
 		if(r != 0){ std::cerr << "FATAL : couldn't get the json file from remote!" << std::endl; return 1;}
 	}
 	if( not std::filesystem::exists("run2DTAsplines.json") and runNo == 2 )
 	{
 		std::cerr << "Warning : couldn't find the json file of splines in Run 2; attempting to fetch from eos ..." << std::endl;
 		int r = system( Form("xrdcp %s/%s/run2DTAsplines.json .", settings.eosRoot.c_str(), settings.eosConfigDir.c_str()) );
 		if(r != 0){ std::cerr << "FATAL : couldn't get the json file from remote!" << std::endl; return 1;}
 	}
 
 	//The below is replaced with the new Lau1DCubicSpline::readFromJson
 	// Create the spline knot positions and
 	// starting Y values, to be fit to dtaHist
 	//const std::vector<Double_t> dtvals  {  0.0,   0.5,   1.0,   2.0,   3.0,   4.0,   7.0,  10.0,  15.0};
 	//const std::vector<Double_t> effvals {0.000, 0.010, 0.022, 0.035, 0.042, 0.050, 0.051, 0.052, 0.055};
 	//const std::vector<Double_t> corvals { 1.10,  1.08,  1.06,  1.04,  1.02,  1.00,  0.98,  0.96,  0.94};
 
 	std::unique_ptr<Lau1DCubicSpline> dtEffSpline = nullptr;
 	Lau1DCubicSpline* dtEffSplinePtr = nullptr;
 	std::unique_ptr<Lau1DCubicSpline> dtCorrSpline = nullptr;
 	Lau1DCubicSpline* dtCorrSplinePtr = nullptr;
 
 	LauSplineDecayTimeEfficiency<LauSplineOrder::Cubic>* sigDTA = nullptr;
 
 	std::unique_ptr<LauSplineDecayTimeEfficiency<LauSplineOrder::Cubic>> dtaModel;
 
 	switch(settings.timeEffModel)
 	{
 		case LauDecayTime::EfficiencyMethod::Spline:
 			{
 				fitModel->setASqMaxValue(0.01);
 				if ( settings.dType == LauTimeDepFitModel::CPEigenvalue::QFS ) {
 					fitModel->setASqMaxValue(0.009);
 				}
 
 				//auto dtEffSpline = std::make_unique<Lau1DCubicSpline>( dtvals, effvals, Lau1DCubicSpline::SplineType::AkimaSpline );
 
 				TString jsonFilename = Form("run%uDTAsplines.json",runNo);
 				dtEffSpline = Lau1DCubicSpline::readFromJson(jsonFilename, Form("signalMC_Kpi_Run%u_DTAhist",runNo));
 				dtEffSplinePtr = dtEffSpline.get();
 
 				dtaModel = std::make_unique<LauSplineDecayTimeEfficiency<LauSplineOrder::Cubic>>( Form("dteff_QFS_Run%u",runNo), std::move(dtEffSpline) );
 				dynamic_cast<LauParameter*>(dtaModel->getParameters().front())->maxValue(1e-1);
 				//dtaModel->fitToTH1(dtaHist);
 
 				// Set which knots to float and which to fix (at least 1 knot must be fixed, not the first one)
 				// Knots should only be floating if requested AND the B lifetime is fixed!
 				if ( settings.fixSplineKnots or not settings.fixLifetime ) {
 					dtaModel->fixKnots();
 				} else {
 					dtaModel->floatKnots();
 					dtaModel->fixKnot( 0, true );
 					dtaModel->fixKnot( 3, true );
 				}
 				sigDTA = dtaModel.get();
 				if ( settings.dType == LauTimeDepFitModel::CPEigenvalue::QFS ) {
 					// For the QFS mode we just use the cubic model as it is
 					sigDTA = dtaModel.get();
 					dtPdf->setSplineEfficiency( std::move(dtaModel) );
 				} else {
 					// For the CP modes we modify it using a corrective spline
 					// Json available?
 					//auto dtCorrSpline = std::make_unique<Lau1DCubicSpline>( dtvals, corvals, Lau1DCubicSpline::SplineType::AkimaSpline );
 					std::unique_ptr<LauSplineDecayTimeEfficiency<LauSplineOrder::Sixth>> dtaCPModel;
 					if (settings.directory == "B2Dhh_KK"){
 						dtCorrSpline = Lau1DCubicSpline::readFromJson(jsonFilename, Form("signalMC_KK_Run%u_DTAhist",runNo));
 						dtCorrSplinePtr = dtCorrSpline.get();
 						dtaCPModel = std::make_unique<LauSplineDecayTimeEfficiency<LauSplineOrder::Sixth>>( Form("dteff_CPKK_Run%u",runNo), *dtaModel, std::move( dtCorrSpline ) );
 					} else if (settings.directory == "B2Dhh_pipi"){
 						dtCorrSpline = Lau1DCubicSpline::readFromJson(jsonFilename, Form("signalMC_pipi_Run%u_DTAhist",runNo));
 						dtCorrSplinePtr = dtCorrSpline.get();
 						dtaCPModel = std::make_unique<LauSplineDecayTimeEfficiency<LauSplineOrder::Sixth>>( Form("dteff_CPpipi_Run%u",runNo), *dtaModel, std::move( dtCorrSpline ) );
 					}
 					dtPdf->setSplineEfficiency( std::move(dtaCPModel) );
 				}
 
 				break;
 			}
 
 		case LauDecayTime::EfficiencyMethod::Binned:
 			{
 				fitModel->setASqMaxValue(0.06);
 
 				auto dtaBinnedModel = std::make_unique<LauBinnedDecayTimeEfficiency>( *dtaHist );
 				dtPdf->setBinnedEfficiency( std::move(dtaBinnedModel) );
 
 				break;
 			}
 
 		case LauDecayTime::EfficiencyMethod::Uniform:
 			{
 				fitModel->setASqMaxValue(4.45);
 				break;
 			}
 	}
 
 	fitModel->setSignalDtPdf( dtPdf );
 
 	std::map<TString, LauDecayTime::FuncType> bkgndDecayTimeTypes = {
 		{"comb",    LauDecayTime::FuncType::Hist},
 		{"Bd2DKpi", LauDecayTime::FuncType::ExpTrig},
 		{"Bs2DKpi", LauDecayTime::FuncType::ExpHypTrig},
 		{"Lb2Dppi", LauDecayTime::FuncType::Exp}
 	};
 
 	TH1* bgDTA{nullptr};
 	std::map<TString, std::pair<Double_t, Double_t > > decayFlvYields = {
 		{"comb"   , runNo == 1 ? std::make_pair( 1111, 1044) : std::make_pair( 3011, 2795)},
 		{"Bd2DKpi", runNo == 1 ? std::make_pair( 9552, 8933) : std::make_pair(13842,12036)},
 		{"Bs2DKpi", runNo == 1 ? std::make_pair(16989,15634) : std::make_pair(18294,17005)},
 		{"Lb2Dppi", runNo == 1 ? std::make_pair( 3695, 3431) : std::make_pair( 2389, 2501)}
 	};
 	TString DTAnameSuffix{""};
 
 	if (settings.directory == "B2Dhh_Kpi"){
 		DTAnameSuffix = Form("_Kpi_Run%u_DTAhist",runNo);
 	} else if (settings.directory == "B2Dhh_KK"){
 		DTAnameSuffix = Form("_KK_Run%u_DTAhist",runNo);
 	} else if (settings.directory == "B2Dhh_pipi"){
 		DTAnameSuffix = Form("_pipi_Run%u_DTAhist",runNo);
 	}
 
 	for (auto& [bg, _] : BkgndTypes) {
 
 		if ( bkgndDecayTimeTypes[bg] == LauDecayTime::FuncType::Hist ) {
 
 			TString histname = bg + DTAnameSuffix;
 			if(histname.Contains("comb"))
 			{
 				histname = "Sideband" + DTAnameSuffix;
 				histname.ReplaceAll("DTAhist","DThist");
 			}
 			TH1* dt_hist = dynamic_cast<TH1*>(dtaFile->Get(histname.Data()));
 			TH1* dt_err_hist = nullptr; // TODO get these!
 
 			if ( not dt_hist ) { std::cerr << "PROBLEM FINDING THE LIFETIME HIST:" << histname << std::endl; return EXIT_FAILURE; }
 			if ( settings.perEventTimeErr and not dt_err_hist ) { std::cerr << "PROBLEM FINDING THE LIFETIME ERROR HIST:" << bg << std::endl; return EXIT_FAILURE; }
 
 			LauDecayTimePdf* backgroundPdf {nullptr};
 			if ( settings.perEventTimeErr ) {
 				backgroundPdf = new LauDecayTimePdf("decayTime", minDt, maxDt, "decayTimeErr", minDtErr, maxDtErr, dt_hist, dt_err_hist );
 			} else {
 				backgroundPdf = new LauDecayTimePdf("decayTime", minDt, maxDt, dt_hist );
 			}
 
 			fitModel->setBkgndDtPdf(bg,backgroundPdf);
 
 		} else {
 
 			auto& bgPhysModel { dtPhysicsModels.at(bg) };
 			auto& bgResoModel { dtResolutionModels.at(bg) };
 
 			TH1* dt_err_hist = nullptr; // TODO get these!
 			if ( settings.perEventTimeErr and not dt_err_hist ) { std::cerr << "PROBLEM FINDING THE LIFETIME ERROR HIST:" << bg << std::endl; return EXIT_FAILURE; }
 
 			LauDecayTimePdf* backgroundPdf {nullptr};
 			if ( settings.timeResolution ) {
 				if ( settings.perEventTimeErr ) {
 					std::cerr << "Warning : the resolution model has not been constructed with per-event error scaling, but you have requested it" << std::endl;
 					backgroundPdf = new LauDecayTimePdf( "decayTime", minDt, maxDt, "decayTimeErr", minDtErr, maxDtErr, std::move(bgPhysModel), std::move(bgResoModel), dt_err_hist );
 				} else {
 					backgroundPdf = new LauDecayTimePdf( "decayTime", minDt, maxDt, std::move(bgPhysModel), std::move(bgResoModel) );
 				}
 			} else {
 				backgroundPdf = new LauDecayTimePdf( "decayTime", minDt, maxDt, std::move(bgPhysModel) );
 			}
 
 			// background decay time acceptance histogram
 			TString histname = bg + DTAnameSuffix;
 			bgDTA = dynamic_cast<TH1*>(dtaFile->Get(histname));
 			if ( not bgDTA and bg!="comb" ) {
 				std::cerr << "Couldn't find DTA histogram for background: " << bg << std::endl;
 				return EXIT_FAILURE;
 			}
 
 			//XXX if we decide to rebin, do it here
 			bgDTA->Divide(dtaHist);
 			for (Int_t bin = 1; bin <= bgDTA->GetNbinsX(); ++bin)
 			{
 				if( bgDTA->GetBinContent(bin)==0. ){bgDTA->SetBinContent(bin,1.);}
 			}
 
 			TString jsonFilename = Form("run%iDTAsplines.json",runNo);
 			auto dtBkgndCorrSpline = Lau1DCubicSpline::readFromJson(jsonFilename,histname);
 
 			std::unique_ptr<LauSplineDecayTimeEfficiency<LauSplineOrder::Sixth>> bgDTAModel;
 			if (sigDTA==nullptr){
 				std::cout << "Null pointer" << std::endl;
 			}
 			if (settings.directory == "B2Dhh_Kpi"){
 				bgDTAModel = std::make_unique<LauSplineDecayTimeEfficiency<LauSplineOrder::Sixth>>( Form("%s_DTA_ratio_Run%u_Kpi",bg.Data(), runNo), *sigDTA, std::move( dtBkgndCorrSpline ) );
 			} else if (settings.directory == "B2Dhh_KK"){
 				bgDTAModel = std::make_unique<LauSplineDecayTimeEfficiency<LauSplineOrder::Sixth>>( Form("%s_DTA_ratio_Run%u_KK",bg.Data(), runNo), *sigDTA, std::move( dtBkgndCorrSpline ) );
 			} else if (settings.directory == "B2Dhh_pipi"){
 				bgDTAModel = std::make_unique<LauSplineDecayTimeEfficiency<LauSplineOrder::Sixth>>( Form("%s_DTA_ratio_Run%u_pipi",bg.Data(),runNo), *sigDTA, std::move( dtBkgndCorrSpline ) );
 			}
 
 			//bgDTAModel->fitToTH1(bgDTA, LauOutputLevel::Quiet, true);
 			backgroundPdf->setSplineEfficiency( std::move(bgDTAModel) );
 
 			fitModel->setBkgndDtPdf(bg,backgroundPdf);
 		}
 	}
 
 	// Signal and bkgnd yields
 	Double_t nSigEvents{0};
 	Double_t nCombEvents{0};
 	std::map<TString, Double_t> bkgndYields;
 	if (settings.directory == "B2Dhh_Kpi") {
 		nSigEvents = (runNo == 1) ? 28362 : 115032;
 		nCombEvents = (runNo == 1) ? 2415 : 7637;
 		if ( ! settings.scaleYields ) {
 			nCombEvents += (runNo == 1) ? (303 + 111 + 123 + 143 + 2) : (974 + 390 + 256 + 389 + 4);
 		}
 		bkgndYields = {{"comb", nCombEvents },
 			  {"Bs2DKpi", runNo == 1 ? 796 : 1713  },
 			  {"Bd2DKpi", runNo == 1 ? 254 : 639   },
 			  {"Lb2Dppi", runNo == 1 ? 245 : 304  }};
 	} else if (settings.directory == "B2Dhh_KK"){
 		if (settings.setNegToZero and runNo == 1) {
 			// Mass fit with mis-ID yields fixed to 0
 			nSigEvents = 3098;
 			nCombEvents = 800;
 			if( ! settings.scaleYields){
 				nCombEvents += 20 + 22 + 53 + 39 + 0;
 			}
 			bkgndYields = {{"comb", nCombEvents },
 				{"Bs2DKpi", 0  },
 				{"Bd2DKpi", 0  },
 				{"Lb2Dppi", 0  }};
 		}
 		else if (settings.constrainedYields and runNo == 1)
 		{
 
 			// Mass fit with mis-ID yields constrained based on
 			// Run1/Run2 in Kpi and pipi and Run2 yield of KK
 			nSigEvents = 3093;
 			nCombEvents = 777;
 			if( ! settings.scaleYields){
 				nCombEvents += 19 + 21 + 52 + 41 + 0;
 			}
 			bkgndYields = {{"comb", nCombEvents },
 				{"Bs2DKpi", 19  },
 				{"Bd2DKpi", 6  },
 				{"Lb2Dppi", 7  }};
 		} else {
 			// Original mass fit
 			nSigEvents = runNo == 1 ? 3140 : 12099;
 			nCombEvents = (runNo == 1) ? 1016 : 1501;
 			if ( ! settings.scaleYields ) {
 				nCombEvents += (runNo == 1) ? (63 + 23 + 26 + 13 + 1) : (107 + 43 + 28 + 5 + 1);
 			}
 			bkgndYields = {{"comb", nCombEvents },
 				{"Bs2DKpi", runNo == 1 ? -172 : 167  },
 				{"Bd2DKpi", runNo == 1 ? -64  : 73   },
 				{"Lb2Dppi", runNo == 1 ? -61  : 33   }};
 		}
 
 	} else if (settings.directory == "B2Dhh_pipi"){
 		nSigEvents = runNo == 1 ? 1175 : 4566;
 		nCombEvents = (runNo == 1) ? 115 : 301;
 		if ( ! settings.scaleYields ) {
 			nCombEvents += (runNo == 1) ? (7 + 3 + 3 + 79 + 0) : (36 + 15 + 10 + 138 + 1);
 		}
 		bkgndYields = {{"comb", nCombEvents },
 			  {"Bs2DKpi", runNo == 1 ? 57 : 99  },
 			  {"Bd2DKpi", runNo == 1 ? 18 : 41   },
 			  {"Lb2Dppi", runNo == 1 ? 21 : 19   }};
 	} else {std::cerr << "bad dir given!" << std::endl; return EXIT_FAILURE;}
 
 	if(settings.command == Command::Generate)
 	{
 		for(auto& [mode, yield] : bkgndYields)
 		{
 			if(yield < 0.)
 			{
 				std::cerr << "WARNING : setting negative yield for mode " << mode << " to 0 for generation" << std::endl;
 				yield = 0.;
 			}
 		}
 	}
 
 	const Bool_t fixYields { settings.bkgndList.empty() or not settings.floatYields };
 
 	TString yieldName;
 	for (auto& [bg, _] : BkgndTypes) {
 		const Double_t nBkg = bkgndYields[bg];
 		std::cout<<bg<<" Events = "<<nBkg<<std::endl;
 		yieldName = Form("%sEvents_%s_Run%u", bg.Data(), settings.directory.c_str(), runNo);
 		LauParameter* nBkgndEvents = new LauParameter(yieldName,nBkg,-2.0*TMath::Abs(nBkg),2.0*TMath::Abs(nBkg),fixYields);
 		fitModel->setNBkgndEvents( nBkgndEvents );
 	}
 
 	std::cout<<"nSigEvents = "<<nSigEvents<<std::endl;
 	yieldName = Form("signalEvents_%s_Run%u", settings.directory.c_str(), runNo);
 	LauParameter* nSigPar = new LauParameter(yieldName, nSigEvents, -2.0*nSigEvents, 2.0*nSigEvents, fixYields);
 	fitModel->setNSigEvents(nSigPar);
 
 	if ( settings.command == Command::Generate or not fixYields ) {
 		// mB PDFs
 		// PDFs for B+ and B- are set to be identical here but could be different in general
 		const TString mBName { "mB" };
 		const Double_t mBMin { 5.200 };
 		const Double_t mBMax { 5.700 };
 
 		// Signal PDF is a double Crystal Ball function
 		const std::map<TString,Double_t> sigParVals {
 			{"mB_sig_mean",    5.2802},
 				{"mB_sig_sigma",   0.01125},
 				{"mB_sig_alphaL",  1.226},
 				{"mB_sig_alphaR", -1.722},
 				{"mB_sig_orderL",  1.817},
 				{"mB_sig_orderR",  2.861},
 				{"mB_sig_frac",    0.480}
 		};
 		fitModel->setSignalPdfs( makeDCBPDF( mBName, mBMin, mBMax, "sig", sigParVals ) );
 
 		// Combinatorial PDF is an exponential
 		if ( BkgndTypes.find("comb") != BkgndTypes.end() ) {
 			LauParameter* mB_comb_slope = new LauParameter("mB_comb_slope", -0.00419);
 			LauAbsPdf* mB_comb_pdf = new LauExponentialPdf(mBName, {mB_comb_slope}, mBMin, mBMax);
 			fitModel->setBkgndPdf("comb", mB_comb_pdf);
 		}
 
 		// Bd2DKpi PDF is double Crystal Ball function
 		if ( BkgndTypes.find("Bd2DKpi") != BkgndTypes.end() ) {
 			const std::map<TString,Double_t> Bd2DKpiParVals {
 				{"mB_Bd2DKpi_mean",    5.2401},
 					{"mB_Bd2DKpi_sigma",   0.01580},
 					{"mB_Bd2DKpi_alphaL",  0.560},
 					{"mB_Bd2DKpi_alphaR", -2.212},
 					{"mB_Bd2DKpi_orderL",  1.076},
 					{"mB_Bd2DKpi_orderR",  2.453},
 					{"mB_Bd2DKpi_frac",    0.544}
 			};
 			fitModel->setBkgndPdf("Bd2DKpi", makeDCBPDF( mBName, mBMin, mBMax, "Bd2DKpi", Bd2DKpiParVals ) );
 		}
 
 		// Bs2DKpi PDF is double Crystal Ball function
 		if ( BkgndTypes.find("Bs2DKpi") != BkgndTypes.end() ) {
 			const std::map<TString,Double_t> Bs2DKpiParVals {
 				{"mB_Bs2DKpi_mean",    5.3244},
 					{"mB_Bs2DKpi_sigma",   0.01827},
 					{"mB_Bs2DKpi_alphaL",  0.865},
 					{"mB_Bs2DKpi_alphaR", -0.001},
 					{"mB_Bs2DKpi_orderL", 10.000},
 					{"mB_Bs2DKpi_orderR", 10.000},
 					{"mB_Bs2DKpi_frac",    0.954}
 			};
 			fitModel->setBkgndPdf("Bs2DKpi", makeDCBPDF( mBName, mBMin, mBMax, "Bs2DKpi", Bs2DKpiParVals ) );
 		}
 
 		// Lb2Dppi PDF is double Crystal Ball function
 		if ( BkgndTypes.find("Lb2Dppi") != BkgndTypes.end() ) {
 			const std::map<TString,Double_t> Lb2DppiParVals {
 				{"mB_Lb2Dppi_mean",     5.4884},
 					{"mB_Lb2Dppi_sigma",    0.02561},
 					{"mB_Lb2Dppi_alphaL",   0.091},
 					{"mB_Lb2Dppi_alphaR", -18.053},
 					{"mB_Lb2Dppi_orderL",   2.906},
 					{"mB_Lb2Dppi_orderR",   6.087},
 					{"mB_Lb2Dppi_frac",     0.968}
 			};
 			fitModel->setBkgndPdf("Lb2Dppi", makeDCBPDF( mBName, mBMin, mBMax, "Lb2Dppi", Lb2DppiParVals ) );
 		}
 	}
 
 	// set the number of experiments
 	if (settings.command == Command::Generate) {
 		fitModel->setNExpts(settings.nExptGen, settings.firstExptGen);
 	} else {
 		fitModel->setNExpts(settings.nExptFit, settings.firstExptFit);
 	}
 
 	fitModel->useAsymmFitErrors(kFALSE);
 	fitModel->useRandomInitFitPars(kFALSE);
 	fitModel->writeLatexTable(kFALSE);
 	const Bool_t haveBkgnds = ( fitModel->nBkgndClasses() > 0 );
 	fitModel->doPoissonSmearing(haveBkgnds);
 	fitModel->doEMLFit(haveBkgnds);
 
 	TString dTypeStr;
 	switch (settings.dType) {
 		case LauTimeDepFitModel::CPEigenvalue::CPEven :
 			dTypeStr = "CPEven";
 			break;
 		case LauTimeDepFitModel::CPEigenvalue::CPOdd :
 			dTypeStr = "CPOdd";
 			break;
 		case LauTimeDepFitModel::CPEigenvalue::QFS :
 			dTypeStr = "QFS";
 			break;
 	}
 
 	TString dataFile("");
 	TString treeName("fitTree");
 	TString rootFileName("");
 	TString tableFileName("");
 	TString fitToyFileName("");
 	TString splotFileName("");
 
 	dataFile = "TEST-Dpipi";
 
 	dataFile += "_"+dTypeStr+"_"+settings.directory;
 
 	switch(settings.timeEffModel)
 	{
 		case LauDecayTime::EfficiencyMethod::Spline:
 			dataFile += "_Spline";
 			break;
 		case LauDecayTime::EfficiencyMethod::Binned:
 			dataFile += "_Binned";
 			break;
 		case LauDecayTime::EfficiencyMethod::Uniform:
 			dataFile += "_Uniform";
 			break;
 	}
 
 	if (settings.timeResolution) {
 		if (settings.perEventTimeErr) {
 			dataFile += "_DTRperevt";
 		} else {
 			dataFile += "_DTRavg";
 		}
 	} else {
 		dataFile += "_DTRoff";
 	}
 
 	dataFile += "_expts";
 	dataFile += settings.firstExptGen;
 	dataFile += "-";
 	dataFile += settings.firstExptGen+settings.nExptGen-1;
 	dataFile += Form("-Run%u",runNo);
 
 	dataFile += ".root";
 
 	if (settings.dataFit)
 	{
 		if(runNo == 2)
 		{
 			dataFile = Form("%s/%s/Run2/B2Dhh-CollisionCombined-MagCombined-Stripping24r2-withDNN-withPIDcorr-withMVA-withPIDMVA-WithDpiMatching-withSwappedMassHypotheses-withKKvetoes-withSelection-LauraPrepped.root",settings.eosRoot.c_str(),settings.eosDataDir.c_str());
 		} else {
 			dataFile = Form("%s/%s/Run1/B2Dhh-CollisionCombined-MagCombined-Stripping21r1p2-withDNN-withPIDcorr-withMVA-withPIDMVA-WithDpiMatching-withSwappedMassHypotheses-withKKvetoes-withSelection-LauraPrepped.root",settings.eosRoot.c_str(),settings.eosDataDir.c_str());
 		}
 		treeName = Form("%s/B2DhhReco",settings.directory.data());
 	}
 
 	if (settings.command == Command::Generate) {
 		rootFileName = "dummy.root";
 		tableFileName = "genResults";
 	} else {
 		rootFileName = "fit"; rootFileName += settings.iFit;
 		rootFileName += "_Results_"; rootFileName += dTypeStr;
 		rootFileName += "_"; rootFileName += settings.directory;
 		rootFileName += "_expts"; rootFileName += settings.firstExptFit; rootFileName += "-"; rootFileName += settings.firstExptFit+settings.nExptFit-1;
 		rootFileName += Form("_Run%u",runNo);
 		rootFileName += ".root";
 
 		tableFileName = "fit"; tableFileName += settings.iFit;
 		tableFileName += "_Results_"; tableFileName += dTypeStr;
 		tableFileName += "_"; tableFileName += settings.directory;
 		tableFileName += "_expts"; tableFileName += settings.firstExptFit; tableFileName += "-"; tableFileName += settings.firstExptFit+settings.nExptFit-1;
 		tableFileName += Form("_Run%u",runNo);
 
 		fitToyFileName = "fit"; fitToyFileName += settings.iFit;
 		fitToyFileName += "_ToyMC_"; fitToyFileName += dTypeStr;
 		fitToyFileName += "_"; fitToyFileName += settings.directory;
 		fitToyFileName += "_expts"; fitToyFileName += settings.firstExptFit; fitToyFileName += "-"; fitToyFileName += settings.firstExptFit+settings.nExptFit-1;
 		fitToyFileName += Form("_Run%u",runNo);
 		fitToyFileName += ".root";
 
 		splotFileName = "fit"; splotFileName += settings.iFit;
 		splotFileName += "_sPlot_"; splotFileName += dTypeStr;
 		splotFileName += "_"; splotFileName += settings.directory;
 		splotFileName += "_expts"; splotFileName += settings.firstExptFit; splotFileName += "-"; splotFileName += settings.firstExptFit+settings.nExptFit-1;
 		splotFileName += Form("_Run%u",runNo);
 		splotFileName += ".root";
 	}
 
 	if( settings.genToy and not settings.blindFit )
 	{
 		fitModel->compareFitData(50, fitToyFileName);
 	}
 
 	// Generate toy from the fitted parameters
 	//fitModel->compareFitData(1, fitToyFileName);
 
 	// Write out per-event likelihoods and sWeights
 	//fitModel->writeSPlotData(splotFileName, "splot", kFALSE);
 
 	// Execute the generation/fit
 	switch (settings.command) {
 		case Command::Generate :
 			fitModel->run( "gen", dataFile, treeName, rootFileName, tableFileName );
 			break;
 		case Command::Fit :
 			fitModel->run( "fit", dataFile, treeName, rootFileName, tableFileName );
 			break;
 		case Command::SimFit :
 			fitModel->runTask( dataFile, treeName, rootFileName, tableFileName, "localhost", settings.port );
 			break;
 	}
 
 	return EXIT_SUCCESS;
 }
 
 LauAbsPdf* makeDCBPDF( const TString& varName, const Double_t varMin, const Double_t varMax, const TString& componentName, const std::map<TString,Double_t> paramVals )
 {
 	const TString meanName   { Form( "%s_%s_mean",   varName.Data(), componentName.Data() ) };
 	const TString sigmaName  { Form( "%s_%s_sigma",  varName.Data(), componentName.Data() ) };
 	const TString alphaLName { Form( "%s_%s_alphaL", varName.Data(), componentName.Data() ) };
 	const TString alphaRName { Form( "%s_%s_alphaR", varName.Data(), componentName.Data() ) };
 	const TString orderLName { Form( "%s_%s_orderL", varName.Data(), componentName.Data() ) };
 	const TString orderRName { Form( "%s_%s_orderR", varName.Data(), componentName.Data() ) };
 	const TString fracName   { Form( "%s_%s_frac",   varName.Data(), componentName.Data() ) };
 
 	LauParameter* mean   = new LauParameter( meanName,   paramVals.at( meanName )   );
         LauParameter* sigma  = new LauParameter( sigmaName,  paramVals.at( sigmaName )  );
 	LauParameter* alphaL = new LauParameter( alphaLName, paramVals.at( alphaLName ) );
         LauParameter* alphaR = new LauParameter( alphaRName, paramVals.at( alphaRName ) );
 	LauParameter* orderL = new LauParameter( orderLName, paramVals.at( orderLName ) );
         LauParameter* orderR = new LauParameter( orderRName, paramVals.at( orderRName ) );
 	LauParameter* frac   = new LauParameter( fracName,   paramVals.at( fracName )   );
 
 	std::vector<LauAbsRValue*> pars;
 	pars.reserve(4);
 
 	pars.clear();
         pars.push_back(mean);
         pars.push_back(sigma);
         pars.push_back(alphaL);
         pars.push_back(orderL);
 	LauAbsPdf* pdf_L = new LauCrystalBallPdf(varName, pars, varMin, varMax);
 
 	pars.clear();
         pars.push_back(mean);
         pars.push_back(sigma);
         pars.push_back(alphaR);
         pars.push_back(orderR);
 	LauAbsPdf* pdf_R = new LauCrystalBallPdf(varName, pars, varMin, varMax);
 
 	LauAbsPdf* pdf = new LauSumPdf(pdf_L, pdf_R, frac);
 
 	return pdf;
 }
diff --git a/examples/runTimeDepTest.sh b/examples/runTimeDepTest.sh
index 6749f0c..fb2971c 100755
--- a/examples/runTimeDepTest.sh
+++ b/examples/runTimeDepTest.sh
@@ -1,254 +1,260 @@
 #!/bin/bash
 #SBATCH --partition=epp
 #SBATCH --ntasks=1
 #SBATCH --cpus-per-task=6
 #SBATCH --mem-per-cpu=3997
 #SBATCH --time=02:00:00
 #Run with sbatch runTimeDepTest.sh
 
 # Copyright 2020 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
 
 processes=()
 
 # trap ^C and call ctrl_c()
 trap ctrl_c INT
 
 function ctrl_c()
 {
 	echo " ** Caught ^C event!"
 	for process in ${processes[@]}
 	do
 		if ps -u $USER | grep -q "$process"
 		then
 			echo "Killing : $process"
 			kill $process
 		fi
 	done
 	exit 1
 }
 
 # Check that we have a kerberos token, since nothing will work without one!
 klist -s > /dev/null 2>&1
 if [ $? -ne 0 ]
 then
     echo "You don't have a valid Kerberos ticket!"
     exit 1
 fi
 
 iFit=0
 dta_model="spline"
 dtr=1
 dtr_perevent=0
 fixTau=1
 fixDeltaM=1
 seed=140279 # $(pwd | xargs basename)
 dataFit=1
 runNo=12
 
 # Modify these if you want to run more experiments
 nExpt=1
 firstExpt=0
 
 logName_genQFS_run1=Logs/gen-QFS-run1-dtamodel_$dta_model-dtr_$dtr-dtrperevent_$dtr_perevent.log
 logName_genQFS_run2=Logs/gen-QFS-run2-dtamodel_$dta_model-dtr_$dtr-dtrperevent_$dtr_perevent.log
 logName_genCPE_KK_run1=Logs/gen-CPEvenKK-run1-dtamodel_$dta_model-dtr_$dtr-dtrperevent_$dtr_perevent.log
 logName_genCPE_KK_run2=Logs/gen-CPEvenKK-run2-dtamodel_$dta_model-dtr_$dtr-dtrperevent_$dtr_perevent.log
 logName_genCPE_pipi_run1=Logs/gen-CPEvenpipi-run1-dtamodel_$dta_model-dtr_$dtr-dtrperevent_$dtr_perevent.log
 logName_genCPE_pipi_run2=Logs/gen-CPEvenpipi-run2-dtamodel_$dta_model-dtr_$dtr-dtrperevent_$dtr_perevent.log
 
 logName_coord=Logs/coordinator-dtamodel_$dta_model-dtr_$dtr-dtrperevent_$dtr_perevent-$iFit.log
 logName_taskQFS_run1=Logs/task-QFS-run1-dtamodel_$dta_model-dtr_$dtr-dtrperevent_$dtr_perevent-$iFit.log
 logName_taskQFS_run2=Logs/task-QFS-run2-dtamodel_$dta_model-dtr_$dtr-dtrperevent_$dtr_perevent-$iFit.log
 logName_taskCPE_KK_run1=Logs/task-CPEvenKK-run1-dtamodel_$dta_model-dtr_$dtr-dtrperevent_$dtr_perevent-$iFit.log
 logName_taskCPE_KK_run2=Logs/task-CPEvenKK-run2-dtamodel_$dta_model-dtr_$dtr-dtrperevent_$dtr_perevent-$iFit.log
 logName_taskCPE_pipi_run1=Logs/task-CPEvenpipi-run1-dtamodel_$dta_model-dtr_$dtr-dtrperevent_$dtr_perevent-$iFit.log
 logName_taskCPE_pipi_run2=Logs/task-CPEvenpipi-run2-dtamodel_$dta_model-dtr_$dtr-dtrperevent_$dtr_perevent-$iFit.log
 
 mkdir -p Logs
 
 # make sure we have the most up-to-date json files
 eosArea="root://eoslhcb.cern.ch//eos/lhcb/wg/b2oc/B2D0pipi/timedep_fit_histos"
 rm -f run1DTAsplines.json run2DTAsplines.json Bd2D0pipi_DP_Model_Coeffs.json Bd2D0pipi_DP_Vetoes.json Bd2D0pipi_DecayTime.json
 xrdcp $eosArea/run1DTAsplines.json .
 if [ $? -ne 0 ]
 then
     echo "Problem downloading json files from EOS, exiting"
     exit 1
 fi
 xrdcp $eosArea/run2DTAsplines.json .
 if [ $? -ne 0 ]
 then
     echo "Problem downloading json files from EOS, exiting"
     exit 1
 fi
 xrdcp $eosArea/Bd2D0pipi_DP_Model_Coeffs.json .
 if [ $? -ne 0 ]
 then
     echo "Problem downloading json files from EOS, exiting"
     exit 1
 fi
 xrdcp $eosArea/Bd2D0pipi_DP_Vetoes.json .
 if [ $? -ne 0 ]
 then
     echo "Problem downloading json files from EOS, exiting"
     exit 1
 fi
 xrdcp $eosArea/Bd2D0pipi_DecayTime.json .
 if [ $? -ne 0 ]
 then
     echo "Problem downloading json files from EOS, exiting"
     exit 1
 fi
+xrdcp $eosArea/Bd2D0pipi_CombTagAsym.json .
+if [ $? -ne 0 ]
+then
+    echo "Problem downloading json files from EOS, exiting"
+    exit 1
+fi
 
 if [ $iFit == 0 -a $dataFit == 0 ]
 then
     #run generation of the modes in parallel (put the tasks in the bg)
     echo "Generating samples..."
     if [ $runNo -eq 1 -o $runNo -eq 12 ]
     then
         GenFitTimeDep_Bd2Dpipi gen --dtype QFS    --sigOnly $sigOnly --fixTau $fixTau --seed $seed --dta-model $dta_model --dtr $dtr --dtr-perevent $dtr_perevent --firstExptGen $firstExpt --nExptGen $nExpt --run 1 --dir B2Dhh_Kpi > $logName_genQFS_run1 2>&1 &
         task=$!
         processes+=( $task )
 
         GenFitTimeDep_Bd2Dpipi gen --dtype CPEven --sigOnly $sigOnly --fixTau $fixTau --seed $seed --dta-model $dta_model --dtr $dtr --dtr-perevent $dtr_perevent --firstExptGen $firstExpt --nExptGen $nExpt --run 1 --dir B2Dhh_KK > $logName_genCPE_KK_run1 2>&1 &
         task=$!
         processes+=( $task )
 
         GenFitTimeDep_Bd2Dpipi gen --dtype CPEven --sigOnly $sigOnly --fixTau $fixTau --seed $seed --dta-model $dta_model --dtr $dtr --dtr-perevent $dtr_perevent --firstExptGen $firstExpt --nExptGen $nExpt --run 1 --dir B2Dhh_pipi > $logName_genCPE_pipi_run1 2>&1 &
         task=$!
         processes+=( $task )
     fi
 
     if [ $runNo -eq 2 -o $runNo -eq 12 ]
     then
         GenFitTimeDep_Bd2Dpipi gen --dtype QFS    --sigOnly $sigOnly --fixTau $fixTau --seed $seed --dta-model $dta_model --dtr $dtr --dtr-perevent $dtr_perevent --firstExptGen $firstExpt --nExptGen $nExpt --run 2 --dir B2Dhh_Kpi > $logName_genQFS_run2 2>&1 &
         task=$!
         processes+=( $task )
 
         GenFitTimeDep_Bd2Dpipi gen --dtype CPEven --sigOnly $sigOnly --fixTau $fixTau --seed $seed --dta-model $dta_model --dtr $dtr --dtr-perevent $dtr_perevent --firstExptGen $firstExpt --nExptGen $nExpt --run 2 --dir B2Dhh_KK > $logName_genCPE_KK_run2 2>&1 &
         task=$!
         processes+=( $task )
 
         GenFitTimeDep_Bd2Dpipi gen --dtype CPEven --sigOnly $sigOnly --fixTau $fixTau --seed $seed --dta-model $dta_model --dtr $dtr --dtr-perevent $dtr_perevent --firstExptGen $firstExpt --nExptGen $nExpt --run 2 --dir B2Dhh_pipi > $logName_genCPE_pipi_run2 2>&1 &
         task=$!
         processes+=( $task )
     fi
 
     wait #wait for the generation to complete
 fi
 
 echo "Running fit $iFit"
 if [ $runNo -eq 12 ]
 then
     nTasks=6
 else
     nTasks=3
 fi
 
 SimFitCoordinator $iFit $nExpt $firstExpt $nTasks > $logName_coord 2>&1 &
 task=$!
 echo "SimFitCoordinator process : $task"
 processes+=( $task )
 
 echo "Initialised coordinator process for fit $iFit"
 
 sleep 5
 
 port=""
 while [ "x$port" == "x" ]
 do
     sleep 2
     port=`tail $logName_coord | grep "Waiting for connection" | awk '{print $NF}'`
 done
 
 echo "Coordinator is listening on port $port"
 echo "Initialising tasks..."
 
 if [ $runNo -eq 1 -o $runNo -eq 12 ]
 then
     GenFitTimeDep_Bd2Dpipi simfit --blindFit $dataFit --dataFit $dataFit --dtype QFS --seed $seed --port $port --iFit $iFit --fixTau $fixTau --dta-model $dta_model --dtr $dtr --dtr-perevent $dtr_perevent --firstExpt $firstExpt --nExpt $nExpt --firstExptGen $firstExpt --nExptGen $nExpt --run 1  --dir B2Dhh_Kpi --fixDm $fixDeltaM > $logName_taskQFS_run1 2>&1 &
     task=$!
     echo "QFS Run 1 process : $task"
     processes+=( $task )
     sleep 5
     GenFitTimeDep_Bd2Dpipi simfit --blindFit $dataFit --dataFit $dataFit --dtype CPEven --seed $seed --port $port --iFit $iFit --fixTau $fixTau --dta-model $dta_model --dtr $dtr --dtr-perevent $dtr_perevent --firstExpt $firstExpt --nExpt $nExpt --firstExptGen $firstExpt --nExptGen $nExpt --run 1  --dir B2Dhh_KK --fixDm $fixDeltaM > $logName_taskCPE_KK_run1 2>&1 &
     task=$!
     echo "CP KK Run 1 process : $task"
     processes+=( $task )
     sleep 5
     GenFitTimeDep_Bd2Dpipi simfit --blindFit $dataFit --dataFit $dataFit --dtype CPEven --seed $seed --port $port --iFit $iFit --fixTau $fixTau --dta-model $dta_model --dtr $dtr --dtr-perevent $dtr_perevent --firstExpt $firstExpt --nExpt $nExpt --firstExptGen $firstExpt --nExptGen $nExpt --run 1  --dir B2Dhh_pipi --fixDm $fixDeltaM > $logName_taskCPE_pipi_run1 2>&1 &
     task=$!
     echo "CP pipi Run 1 process : $task"
     processes+=( $task )
 fi
 if [ $runNo -eq 2 -o $runNo -eq 12 ]
 then
     GenFitTimeDep_Bd2Dpipi simfit --blindFit $dataFit --dataFit $dataFit --dtype QFS --seed $seed --port $port --iFit $iFit --fixTau $fixTau --dta-model $dta_model --dtr $dtr --dtr-perevent $dtr_perevent --firstExpt $firstExpt --nExpt $nExpt --firstExptGen $firstExpt --nExptGen $nExpt --run 2  --dir B2Dhh_Kpi --fixDm $fixDeltaM > $logName_taskQFS_run2 2>&1 &
     task=$!
     echo "QFS Run 2 process : $task"
     processes+=( $task )
     sleep 5
     GenFitTimeDep_Bd2Dpipi simfit --blindFit $dataFit --dataFit $dataFit --dtype CPEven --seed $seed --port $port --iFit $iFit --fixTau $fixTau --dta-model $dta_model --dtr $dtr --dtr-perevent $dtr_perevent --firstExpt $firstExpt --nExpt $nExpt --firstExptGen $firstExpt --nExptGen $nExpt --run 2  --dir B2Dhh_KK --fixDm $fixDeltaM > $logName_taskCPE_KK_run2 2>&1 &
     task=$!
     echo "CP KK Run 2 process : $task"
     processes+=( $task )
     sleep 5
     GenFitTimeDep_Bd2Dpipi simfit --blindFit $dataFit --dataFit $dataFit --dtype CPEven --seed $seed --port $port --iFit $iFit --fixTau $fixTau --dta-model $dta_model --dtr $dtr --dtr-perevent $dtr_perevent --firstExpt $firstExpt --nExpt $nExpt --firstExptGen $firstExpt --nExptGen $nExpt --run 2  --dir B2Dhh_pipi --fixDm $fixDeltaM > $logName_taskCPE_pipi_run2 2>&1 &
     task=$!
     echo "CP pipi Run 2 process : $task"
     processes+=( $task )
 fi
 
 wait #wait for the fits and coordinator to terminate
 
 echo "Fit completed"
 
 echo "Warnings/Errors from coordinator:"
 grep -e ERROR -e WARNING -e Error -e Warning $logName_coord
 
 echo "Warnings/Errors from CPeven KK tasks:"
 if [ $runNo -eq 1 -o $runNo -eq 12 ]
 then
     grep -e ERROR -e WARNING -e Error -e Warning $logName_taskCPE_KK_run1
 fi
 if [ $runNo -eq 2 -o $runNo -eq 12 ]
 then
     grep -e ERROR -e WARNING -e Error -e Warning $logName_taskCPE_KK_run2
 fi
 
 echo "Warnings/Errors from CPeven pipi tasks:"
 if [ $runNo -eq 1 -o $runNo -eq 12 ]
 then
     grep -e ERROR -e WARNING -e Error -e Warning $logName_taskCPE_pipi_run1
 fi
 if [ $runNo -eq 2 -o $runNo -eq 12 ]
 then
     grep -e ERROR -e WARNING -e Error -e Warning $logName_taskCPE_pipi_run2
 fi
 
 echo "Warnings/Errors from QFS tasks:"
 if [ $runNo -eq 1 -o $runNo -eq 12 ]
 then
     grep -e ERROR -e WARNING -e Error -e Warning $logName_taskQFS_run1
 fi
 if [ $runNo -eq 2 -o $runNo -eq 12 ]
 then
     grep -e ERROR -e WARNING -e Error -e Warning $logName_taskQFS_run2
 fi