#ifndef GUARD_TDSUPERANA
#define GUARD_TDSUPERANA
#include "TDHistSet.hh"
#include "TDSubAna.hh"
#include "TDCanvMan.hh"
#include "TDHistSetSet.hh"
#include "TF1.h"
#include "TPaveText.h"
#include <fstream>
#include <string>
#include "TDVnEEE.hh"
#include "TDMassPeakFit.hh"
#include "TDDilSet.hh"
#include "TDRunInfoSet.hh"
#include "TDhwwlvlvGoodRun.hh"
#include "TDhwwlvlvGoodRun13.hh"
#include "TDRunInfo.hh"
//#include "DilLeptonType.hh"
#include "LeptonCat.hh"
#include "GoodRunv2.h"
#include "TDExcelData.hh"
#include "TDExcelPlot.hh"
#include "TDRunSet.hh"
#include "TMultiLayerPerceptron.h"
#include "TMLPAnalyzer.h"
#include "TList.h"
#include "TEventList.h"
#include "Root_Jetnet.hh"
#include "mclimit_csm.h"
#include "TMVA/src/Factory.h"
#include "TMVA/src/MethodTMlpANN.h"
#include "TMVA/src/Tools.h"
//using namespace TMVA;

class TDSuperAna
{  
  public:
    TDSuperAna(bool debug, char * FileLocation, int diltype, std::ofstream& AnatxtFile, int weighting,  std::ofstream& ExcelTextFile, std::ofstream& ExcelDataFile, std::ofstream& ExcelMassFile, std::ofstream& CutTableFile, std::ofstream& cuttabletexfile, int higgsMass):
      _debug(debug),
      _fileLocation(FileLocation),
      _dilset(diltype),
      //CanvMan1(("results/plots_"+_dilset.GetDilString()+".ps").c_str(),"Ana_Plots",1,2),
      _anaFile(AnatxtFile),
      //_dilWeight(),
      _weighting(weighting),
      _excelFile(ExcelTextFile),
      _excelDataFile(ExcelDataFile),
      _hwwlvlvGoodRun(),
      _runstoignore("MCrunstoignore.txt"),
      _excelMassFile(ExcelMassFile),
      //_cutTableFile(),
      //_cutTableTexFile(),
      _higgsMass(higgsMass),
      _minRun(000000), //_minRun(186599),
      _maxRun(999999),   //_maxRun(203799) (212133)?
      _stn(true),      //TD changes dataset and good run version.
      _mircea(false)   //TD crossections are meaningless with this setting, as they feed off samples with separate channels, and there are no good runs.
     {
        _dilWeight=new DilWeights(DilWeights::JOINTPHYSICSV1,&_hwwlvlvGoodRun13,_minRun,_maxRun);
	/*std::cout<<"++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ lumval dilweight used for scaling: "<<_dilWeight->GetLum()<<std::endl;
	TDRunSet first("1to10");
	TDRunSet second("1to10");
	first.Print();
	second.Print();
	first.SetMinMax(3,7);
	second.SetMin(3);
	second.SetMax(7);
	first.Print();
        second.Print();*///TD TDRunSet Test
	if(_stn)_mircea=false;   //TD if stn do not use mircea settings.
	if(!_stn && _maxRun >186598) _maxRun=186598; //TD in gen 5 the runs only go up to 186598.
	
        char buffer[500];
	char buffer2[500];
	_likelihoodtypes.push_back("$M_H$");  //TD start of the lables column in the likelihood table
	sprintf(buffer,"%d", (int)_higgsMass); 
	_likelihoodvalues.push_back(buffer);//TD start of the values column in the likelihood table 
	
	sprintf(buffer,"results/cutTableTexFile-%d.tex", (int)_higgsMass);

	_cutTableTexFile.open(buffer);
	_cutTableTexFile<<"\\documentclass[12pt]{article}"<<std::endl;
	_cutTableTexFile<<"\\usepackage{epsfig}"<<std::endl;
	_cutTableTexFile<<"\\usepackage{verbatim}"<<std::endl;
	_cutTableTexFile<<"\\usepackage{rotating}"<<std::endl;
	_cutTableTexFile<<"\\def\\rightboxfour#1#2#3#4{\\bigskip\\rightline{\\vbox{\\hbox{#1}"<<std::endl;
	_cutTableTexFile<<"        \\hbox{#2}\\hbox{#3}\\hbox{#4}}}\\bigskip}"<<std::endl;
	_cutTableTexFile<<"        \\parskip12pt plus 1pt minus 1pt"<<std::endl;
	_cutTableTexFile<<"        \\topsep0pt plus 1pt"<<std::endl;
	_cutTableTexFile<<"        \\setlength{\\textheight}{22.0cm}"<<std::endl;
	_cutTableTexFile<<"        \\setlength{\\topmargin}{0.0cm}"<<std::endl;
	_cutTableTexFile<<"        \\setlength{\\textwidth}{15.0cm}"<<std::endl;
	_cutTableTexFile<<"        \\setlength{\\oddsidemargin}{-0.5cm}"<<std::endl;
	_cutTableTexFile<<"        \\setlength{\\topsep}{1mm}"<<std::endl;
	_cutTableTexFile<<"        \\renewcommand{\\topfraction}{1.0}"<<std::endl;
	_cutTableTexFile<<"        \\renewcommand{\\bottomfraction}{1.0}"<<std::endl;
	_cutTableTexFile<<"        \\renewcommand{\\textfraction}{0.2}"<<std::endl;
	_cutTableTexFile<<"        \\renewcommand{\\floatpagefraction}{0.5}"<<std::endl;
	_cutTableTexFile<<"        \\pretolerance=10000"<<std::endl;
	_cutTableTexFile<<""<<std::endl;
	_cutTableTexFile<<""<<std::endl;
	_cutTableTexFile<<""<<std::endl;
	_cutTableTexFile<<"        \\begin{document}"<<std::endl;

	sprintf(buffer,"results/cutTableTexFile-%d.txt", (int)_higgsMass);
	_cutTableFile.open(buffer);
        sprintf(buffer,"results/plots_%s_%d-%d.ps", _dilset.GetDilString().c_str(), _weighting, (int)_higgsMass);
        CanvMan1 = new (TDCanvMan)(buffer, "Ana_Plots", 4, 4);
//std::cout<<"sdfghj"<<std::endl;
        //_anaFile<<"////////////// Start of " + _dilset.GetDilString() + " details /////////////"<<std::endl;
        _anaFile<<std::endl;
        _anaFile<<_dilset.GetDilString()<<std::endl;
	sprintf(buffer,"results/fits_%s_%d-%d.txt",_dilset.GetDilString().c_str(), _weighting, (int)_higgsMass);
        fitinfo.open(buffer);
        //if(_stn)_fileLocation="./stn-21-7-06/";
	//if(_stn)_fileLocation="./stn-21-7-06_gdhww/";
//	if(_stn)_fileLocation="./stn-5-10-06/";
//	if(_stn)_fileLocation="./stn-20-10-06/";
//	if(_stn)_fileLocation="./stn-01-11-06/";
//	if(_stn)_fileLocation="./stn-18-10-06/";
	if(_stn)_fileLocation="./stn-29-11-06_v2_56/";
	if(!_stn)_fileLocation="./top-23-08-06/";
	//if(!_stn)_fileLocation="./top-05-09-06_oldZsamples/"; //TD older Z->mumu,ee samples, because of doubt concerning reproduceability.
	if(_mircea)_fileLocation="./top-05-09-06_Zll/";
		//_fileLocation="./top-24-07-06/";
	//_fileLocation="./top-30-07-06/";
	//_fileLocation="./top-23-08-06/";
	//GoodRuns=(*(_hwwlvlvGoodRun.GetRunSet(_dilType.Int2Type(_dilset.GetRepDilType()))))-_runstoignore;
        //GoodRuns.Intersect(*(_hwwlvlvGoodRun.GetRunSet(_dilType.Int2Type(0))));

       //TDRunSet HwwStn_ZeeDataFile("HwwStn_ZeeDataFile_180_0001_v196dev_runlist");
	/*TDRunSet gen5_1001("goodrunlists/goodrun_mu_nosi.list");
	TDRunSet gen6_1001("v13goodrunquery/goodrun_mu_nosi.list");
	gen6_1001.SetMax(186598);
	TDRunSet gen6_gen5=gen6_1001-gen5_1001;
	gen6_gen5.Print();*/
        //#include "runchoose/runchoose"
        GoodRuns=(*(_hwwlvlvGoodRun.GetRunSet(_dilType.Int2Type(_dilset.GetRepDilType()))));//-_runstoignore;
	//GoodRuns.Intersect(*(_hwwlvlvGoodRun.GetRunSet(_dilType.Int2Type(0))));
        //GoodRuns.SetMinMax()
//	GoodRuns.SetMax(186598);                     //TD gen5 cutoff
//	GoodRuns.SetMax(203799);                    //TD 0h cutoff
//	TDRunSet muonruns("./v1_96_data_m/muonruns.list");   //only runs actually in the data
//	GoodRuns.Intersect(muonruns);                   //only runs actually in the Gen6 data
	//GoodRuns.Intersect(HwwStn_ZeeDataFile);
	//GoodRuns.Intersect(_runstoignore);
	MCGoodRuns=(*(_hwwlvlvGoodRun.GetRunSet(_dilType.Int2Type(0))))-_runstoignore;
	//MCGoodRuns.Intersect(GoodRuns);
	//MCGoodRuns=(*(_hwwlvlvGoodRun.GetRunSet(_dilType.Int2Type(_dilset.GetRepDilType()))))-_runstoignore;      //for stn ver6, there is no runset on the mc.
//	MCGoodRuns.Print();
//	MCGoodRuns.SetMax(203799);                   //TD 0h cutoff
//	MCGoodRuns.Intersect(muonruns);               //only runs actually in the Gen6 data
	//MCGoodRuns.Intersect(*(_hwwlvlvGoodRun.GetRunSet(_dilType.Int2Type(_dilset.GetRepDilType()))));
	//MCGoodRuns.SetMax(186598);                 //TD gen5 cutoff
	//MCGoodRuns.Intersect(HwwStn_ZeeDataFile);
	//for rick 
	//TDRunSet temp=(*(_hwwlvlvGoodRun.GetRunSet(_dilType.Int2Type(_dilset.GetRepDilType()))))-_runstoignore;
	//GoodRuns=temp-(*(_hwwlvlvGoodRun.GetRunSet(_dilType.Int2Type(0))));
	/*TDRunInfoSet _zewk6m_runinfo("../MCruns/zewk6m_runs_1.log");
        _zewk6m_runinfo.AddFile("../MCruns/zewk6m_runs_2.log");
        _zewk6m_runinfo.AddFile("../MCruns/zewk6m_runs_3.log");
        _zewk6m_runinfo.AddFile("../MCruns/zewk6m_runs_4.log");
        _zewk6m_runinfo.AddFile("../MCruns/zewk6m_runs_5.log");
        _zewk6m_runinfo.AddFile("../MCruns/zewk6m_runs_6.log");
        GoodRuns.Intersect(_zewk6m_runinfo.DerivedRunSet());*/
        /*TDRunInfoSet _zewk6d_runinfo("../MCruns/zewk6d_runs_1.log");
        _zewk6d_runinfo.AddFile("../MCruns/zewk6d_runs_2.log");
        _zewk6d_runinfo.AddFile("../MCruns/zewk6d_runs_3.log");
        _zewk6d_runinfo.AddFile("../MCruns/zewk6d_runs_4.log");
        _zewk6d_runinfo.AddFile("../MCruns/zewk6d_runs_5.log");
        _zewk6d_runinfo.AddFile("../MCruns/zewk6d_runs_6.log");
        GoodRuns.Intersect(_zewk6d_runinfo.DerivedRunSet());*/
	//GoodRuns.SetMinMax(150145,156487);
	std::cout<<"creating subanas and hist sets"<<std::endl;
 
        HistSets = new TDHistSetSet(CanvMan1, false);
        //int color=2;
        HistSets->AddHistSet(new TDHistSet("Zee","Zee_plots","Zee plots",CanvMan1,26,false));
        HistSets->AddHistSet(new TDHistSet("Zmumu","Zmumu_plots","Zmumu plots",CanvMan1,38,false));
        //HistSets->AddHistSet(new TDHistSet("wexocw_180","wexocw_180_plots","wexocw 180 plots",CanvMan1,++color,false));
        /*//HistSets->AddHistSet(new TDHistSet("ztoptt_180","ztoptt_180_plots","ztoptt 180 plots",CanvMan1,++color,false));
        //HistSets->AddHistSet(new TDHistSet("b7dp0p_180","b7dp0p_180_plots","b7dp0p 180 plots",CanvMan1,++color,false));*/
        //HistSets->AddHistSet(new TDHistSet("HwwStn-ewk6","HwwStn-ewk6_plots","HwwStn-ewk6 plots",CanvMan1,++color,false));
        //HistSets->AddHistSet(new TDHistSet("HwwStn_zewkae_0001","HwwStn_zewkae_0001_plots","HwwStn_zewkae_0001 plots",CanvMan1,++color,false));
        //HistSets->AddHistSet(new TDHistSet("HwwStn_zewkae_0002","HwwStn_zewkae_0002_plots","HwwStn_zewkae_0002 plots",CanvMan1,++color,false));
	
	HistSets->AddHistSet(new TDHistSet("Ztautau","Ztautau_plots","Ztautau 180 plots",CanvMan1,18,false));
	//HistSets->AddHistSet(new TDHistSet("sexo8t","sexo8t_plots","sexo8t 180 plots",CanvMan1,++color,false));
	HistSets->AddHistSet(new TDHistSet("fake","fake_plots","fake plots",CanvMan1,30,false));
	HistSets->AddHistSet(new TDHistSet("WZ","WZ_plots","WZ plots",CanvMan1,22,false));
        HistSets->AddHistSet(new TDHistSet("ZZ","ZZ_plots","ZZ plots",CanvMan1,36,false));
	HistSets->AddHistSet(new TDHistSet("WW","WW_plots","WW plots",CanvMan1,6,false));
	HistSets->AddHistSet(new TDHistSet("wgamma","wgamma_plots","wgamma plots",CanvMan1,7,false));
        HistSets->AddHistSet(new TDHistSet("ttbar","ttbar_plots","ttbar plots",CanvMan1,5,false));
	HistSets->AddHistSet(new TDHistSet("higgsww","higgsww_plots","higgsww plots",CanvMan1,4,false));

  /*      HistSets->AddHistSet(new TDHistSet("muon_data","muon_plots","muon data plots",CanvMan1,41,true));
        HistSets->AddHistSet(new TDHistSet("electron_data","electron_data","electron data plots",CanvMan1,42,true));
	HistSets->AddHistSet(new TDHistSet("plug_data","plug_data","plug data plots",CanvMan1,43,true));*/
	 HistSets->AddHistSet(new TDHistSet("data","data_plots","data plots",CanvMan1,41,true));
	
        //HistSets->AddHistSet(new TDHistSet("HwwStn_ZeeDataFile","HwwStn_ZeeDataFile","data from HwwStn_ZeeDataFile",CanvMan1,++color,true));
        //HistSets->AddHistSet(new TDHistSet("HwwStn_zewk6ae","HwwStn_zewk6ae","HwwStn_zewk6ae plots",CanvMan1,++color,false));
        //HistSets->AddHistSet(new TDHistSet("HwwStn_zewk6m","HwwStn_zewk6m","HwwStn_zewk6m plots",CanvMan1,++color,false));
	//HistSets->AddHistSet(new TDHistSet("HwwStn_zewk6m_v1_88","HwwStn_zewk6m_v1_88","HwwStn_zewk6m_v1_88 plots",CanvMan1,++color,false));
	//HistSets->AddHistSet(new TDHistSet("HwwStn_zewk6m_v1_90","HwwStn_zewk6m_v1_90","HwwStn_zewk6m_v1_90 plots",CanvMan1,++color,false));

	
       // sprintf(buffer,"%s%s","../MCruns/","Zmumu.*.root"); //TD My Soucrce Files
//	sprintf(buffer2,"%s%s","../MCruns/","zewk6m_zlog"); //TD My Zlog Soucrce Files
        //_fileLocation=
//int DataSourceCounter=0;
	//float events=0;
        float xs=0;

	sprintf(buffer,"%s%s",_fileLocation.c_str(),"data/*.root");
        sprintf(buffer2,"%s%s",_fileLocation.c_str(),"data/data_zlog");
        SubAnas.push_back(new TDSubAna("fake",false,_debug,buffer,CanvMan1, diltype, _weighting, GoodRuns, buffer2, _higgsMass, _cutTableTexFile, _minRun, _maxRun, 0, _stn, _mircea, xs));

//TD---------------------------------------------------------temporary position in order of files
xs=GetHwwXS((int)_higgsMass, _stn);
        sprintf(buffer,"%s%s%d%s",_fileLocation.c_str(),"hww/",(int)_higgsMass,"/*.root");
        sprintf(buffer2,"%s%s%d%s",_fileLocation.c_str(),"hww/",(int)_higgsMass,"/hww_zlog");
         std::cout<<buffer<<buffer2<<std::endl;
        SubAnas.push_back(new TDSubAna("higgsww",false,_debug,buffer,CanvMan1, diltype, _weighting, MCGoodRuns, buffer2, _higgsMass, _cutTableTexFile, _minRun, _maxRun, 1, _stn, _mircea, xs));
        if(!_stn && GetSubAna("higgsww")->GetInputEvents()==0)GetSubAna("higgsww")->SetInputEvents(GetHwwEvents((int)_higgsMass, _stn));
//TD---------------------------------------------------------temporary position in order of files

	sprintf(buffer,"%s%s",_fileLocation.c_str(),"data/*.root");
        sprintf(buffer2,"%s%s",_fileLocation.c_str(),"data/data_zlog");
        SubAnas.push_back(new TDSubAna("data",true,_debug,buffer,CanvMan1, diltype, _weighting, GoodRuns, buffer2, _higgsMass, _cutTableTexFile, _minRun, _maxRun, -1, _stn, _mircea, xs));

        xs=355.2*1.4;
	//if(!_stn)3216951;
	sprintf(buffer,"%s%s",_fileLocation.c_str(),"Zmumu/*.root");
        sprintf(buffer2,"%s%s",_fileLocation.c_str(),"Zmumu/Zmumu_zlog");
        SubAnas.push_back(new TDSubAna("Zmumu",false,_debug,buffer,CanvMan1, diltype, _weighting, MCGoodRuns, buffer2, _higgsMass, _cutTableTexFile, _minRun, _maxRun, 0, _stn, _mircea, xs));
	if(!_stn && GetSubAna("Zmumu")->GetInputEvents()==0)GetSubAna("Zmumu")->SetInputEvents(3216951);

	xs=355.2*1.4;
	if(_mircea)xs=24300*0.029;
	sprintf(buffer,"%s%s",_fileLocation.c_str(),"Zee/*.root"); //TD stn 17-7-06
        sprintf(buffer2,"%s%s",_fileLocation.c_str(),"Zee/Zee_zlog"); //TD stn 17-7-06
        SubAnas.push_back(new TDSubAna("Zee",false,_debug, buffer, CanvMan1, diltype, _weighting, MCGoodRuns, buffer2, _higgsMass, _cutTableTexFile, _minRun, _maxRun, 0, _stn, _mircea, xs));
	if(!_stn && GetSubAna("Zee")->GetInputEvents()==0 && !_mircea)GetSubAna("Zee")->SetInputEvents(3108482);
	if(!_stn && GetSubAna("Zee")->GetInputEvents()==0 && _mircea)GetSubAna("Zee")->SetInputEvents(2965289);

	xs=355.2*1.4;
	sprintf(buffer,"%s%s",_fileLocation.c_str(),"Ztautau/*.root");
	sprintf(buffer2,"%s%s",_fileLocation.c_str(),"Ztautau/Ztautau_zlog");
        SubAnas.push_back(new TDSubAna("Ztautau",false,_debug,buffer,CanvMan1, diltype, _weighting, MCGoodRuns, buffer2, _higgsMass, _cutTableTexFile, _minRun, _maxRun, 0, _stn, _mircea, xs));
        if(!_stn && GetSubAna("Ztautau")->GetInputEvents()==0)GetSubAna("Ztautau")->SetInputEvents(3338836);
	
	xs=12.4;
	if(_stn)xs=1.27;
	sprintf(buffer,"%s%s",_fileLocation.c_str(),"WW/*.root");
        sprintf(buffer2,"%s%s",_fileLocation.c_str(),"WW/WW_zlog");
	SubAnas.push_back(new TDSubAna("WW",false,_debug,buffer,CanvMan1, diltype, _weighting, MCGoodRuns, buffer2, _higgsMass, _cutTableTexFile, _minRun, _maxRun, 0, _stn, _mircea, xs));
        if(!_stn && GetSubAna("WW")->GetInputEvents()==0)GetSubAna("WW")->SetInputEvents(419728);

	xs=3.65;
	if(_stn)xs=0.365*0.76; //TD 0.76 is filter req.
	sprintf(buffer,"%s%s",_fileLocation.c_str(),"WZ/*.root");
        sprintf(buffer2,"%s%s",_fileLocation.c_str(),"WZ/WZ_zlog");
        SubAnas.push_back(new TDSubAna("WZ",false,_debug,buffer,CanvMan1, diltype, _weighting, MCGoodRuns, buffer2, _higgsMass, _cutTableTexFile, _minRun, _maxRun, 0, _stn, _mircea, xs));
        if(!_stn && GetSubAna("WZ")->GetInputEvents()==0)GetSubAna("WZ")->SetInputEvents(407438);
	
	xs=1.39;
	if(_stn)xs=2.01*0.23; //TD 0.23 is filter req.
	sprintf(buffer,"%s%s",_fileLocation.c_str(),"ZZ/*.root");
        sprintf(buffer2,"%s%s",_fileLocation.c_str(),"ZZ/ZZ_zlog");
        SubAnas.push_back(new TDSubAna("ZZ",false,_debug,buffer,CanvMan1, diltype, _weighting, MCGoodRuns, buffer2, _higgsMass, _cutTableTexFile, _minRun, _maxRun, 0, _stn, _mircea, xs));
        if(!_stn && GetSubAna("ZZ")->GetInputEvents()==0)GetSubAna("ZZ")->SetInputEvents(397086);
	
	xs=6.7;
	sprintf(buffer,"%s%s",_fileLocation.c_str(),"ttbar/*.root");
        sprintf(buffer2,"%s%s",_fileLocation.c_str(),"ttbar/ttbar_zlog");
        SubAnas.push_back(new TDSubAna("ttbar",false,_debug,buffer,CanvMan1, diltype, _weighting, MCGoodRuns, buffer2, _higgsMass, _cutTableTexFile, _minRun, _maxRun, 0, _stn, _mircea, xs));
	if(!_stn && GetSubAna("ttbar")->GetInputEvents()==0)GetSubAna("ttbar")->SetInputEvents(1142339);
	
	xs=89.4;
	if(_stn)xs=2*21.48*1.34; //TD each ->lepton sample has 21.48*1.34, and there are two of them => *2.
	sprintf(buffer,"%s%s",_fileLocation.c_str(),"Wg/*.root");
        sprintf(buffer2,"%s%s",_fileLocation.c_str(),"Wg/Wg_zlog");
	SubAnas.push_back(new TDSubAna("wgamma",false,_debug,buffer,CanvMan1, diltype, _weighting, MCGoodRuns, buffer2, _higgsMass, _cutTableTexFile, _minRun, _maxRun, 0, _stn, _mircea, xs));
	if(!_stn && GetSubAna("wgamma")->GetInputEvents()==0)GetSubAna("wgamma")->SetInputEvents(367427);

	/*switch((int)_higgsMass)
	{
	  case(110):sprintf(buffer,"%s%s","./mtuples/","wexobw*.root");break;
	  case(120):sprintf(buffer,"%s%s","./mtuples/","wexocw*.root");break;
	  case(130):sprintf(buffer,"%s%s","./mtuples/","wexo5w*.root");break;
	  case(140):sprintf(buffer,"%s%s","./mtuples/","wexo6w*.root");break;
	  case(150):sprintf(buffer,"%s%s","./mtuples/","wexo7w*.root");break;
	  case(160):sprintf(buffer,"%s%s","./mtuples/","wexo8w*.root");break;
	  case(170):sprintf(buffer,"%s%s","./mtuples/","wexo9w*.root");break;
	  case(180):sprintf(buffer,"%s%s","./mtuples/","wexoaw*.root");break;
	  case(190):sprintf(buffer,"%s%s","./mtuples/","wexodw*.root");break;
	  case(200):sprintf(buffer,"%s%s","./mtuples/","wexoew*.root");break;
          default:sprintf(buffer,"%s%s","./mtuples/","wexobw*.root");break; 
	}*///gen5

/*
	xs=GetHwwXS((int)_higgsMass, _stn);
        //sprintf(buffer,"%s%s%d%s",_fileLocation.c_str(),"hww/",(int)_higgsMass,"*.root"); //line is commented out to stop the unmatched comment
//	std::cout<<buffer<<buffer2<<std::endl;

//	sprintf(buffer,"%s%s","./mtuples/","wexobw_180*.root");
        sprintf(buffer2,"%s%s%d%s",_fileLocation.c_str(),"hww/",(int)_higgsMass,"/hww_zlog");
	 std::cout<<buffer<<buffer2<<std::endl;
	SubAnas.push_back(new TDSubAna("higgsww",false,_debug,buffer,CanvMan1, diltype, _weighting, MCGoodRuns, buffer2, _higgsMass, _cutTableTexFile, _minRun, _maxRun, 1, _stn, _mircea, xs));
	if(!_stn && GetSubAna("higgsww")->GetInputEvents()==0)GetSubAna("higgsww")->SetInputEvents(GetHwwEvents((int)_higgsMass, _stn));
 */                                
        //sprintf(buffer,"%s%s",_fileLocation.c_str(),"data/*.root");
	//SubAnas.push_back(new TDSubAna("data",true,_debug,buffer,CanvMan1, diltype, _weighting, GoodRuns, buffer2, _higgsMass, _cutTableTexFile));

        /*sprintf(buffer,"%s%s",_fileLocation.c_str(),"bhmu0d_180.*.root");               // v1_72 muon data
//	sprintf(buffer,"%s%s","./v1_96_data_m/","*.root");                                //stn  v1_97 muon data
        SubAnas.push_back(new TDSubAna("muon_data",true,_debug,buffer,CanvMan1, diltype, _weighting, GoodRuns));
	
        sprintf(buffer,"%s%s",_fileLocation.c_str(),"bhel0d_180.*.root");              // v1_72 elec data
	//sprintf(buffer,"%s%s","./v1_96_data_e/","*.root");                               //stn  v1_97 elec data
        SubAnas.push_back(new TDSubAna("electron_data",true,_debug,buffer,CanvMan1, diltype, _weighting, GoodRuns));
        sprintf(buffer,"%s%s","bexo4d/","HwwStn_bexo4d_180_*.root");	
	SubAnas.push_back(new TDSubAna("plug_data",true,_debug,buffer,CanvMan1, diltype, _weighting, GoodRuns));
	//sprintf(buffer,"%s%s","/data/nglas08/b/stdenis/","HwwStn_ZeeDataFile_180_0001_v196dev.root");
	//SubAnas.push_back(new TDSubAna("HwwStn_ZeeDataFile",true,_debug,buffer,CanvMan1, diltype, _weighting, GoodRuns));
*/
        /*sprintf(buffer,"%s%s","zewk6ae/","HwwStn_zewkae_180_*.root");
        SubAnas.push_back(new TDSubAna("HwwStn_zewk6ae",false,_debug,buffer,CanvMan1, diltype, _weighting, GoodRuns));
        sprintf(buffer,"%s%s","/data/nglas08/b/stdenis/zewk6m/v1_87/","HwwStn_Zmumu_*.root");
        SubAnas.push_back(new TDSubAna("HwwStn_zewk6m",false,_debug,buffer,CanvMan1, diltype, _weighting, GoodRuns));
	sprintf(buffer,"%s%s","/data/nglas08/b/stdenis/zewk6m/v1_88/","HwwStn_Zmumu_*.root");
	SubAnas.push_back(new TDSubAna("HwwStn_zewk6m_v1_88",false,_debug,buffer,CanvMan1, diltype, _weighting, GoodRuns));
	sprintf(buffer,"%s%s","/data/nglas08/b/stdenis/","HwwStn_HWWTest_180_*.v190dev.root");
	SubAnas.push_back(new TDSubAna("HwwStn_zewk6m_v1_90",false,_debug,buffer,CanvMan1, diltype, _weighting, GoodRuns));*/
        std::cout<<"InitSubAnas()"<<std::endl;
        InitSubAnas();
     }

     void SavePlots()
     {
       char buffer[200];
       if(!(HistSets->InitiateIt()))do{
	 sprintf(buffer,"results/hists_%s-%s_%d-%d.root",(*(HistSets->GetIt()))->GetDataSourcePrefix(),_dilset.GetDilString().c_str(),_weighting,(int)_higgsMass);
         (*(HistSets->GetIt()))->SavePlots(buffer);
       }while(!HistSets->IncrementIt());
     }

     void LoadPlots()
     {
       char buffer[200];
       /*std::cout<<(*(HistSets->GetIt()))->GetDataSourcePrefix()<<std::endl;
       std::cout<<_dilset.GetDilString().c_str()<<std::endl;
       std::cout<<_weighting<<std::endl;
       std::cout<<(int)_higgsMass<<std::endl;*/
       if(!(HistSets->InitiateIt()))do{
         sprintf(buffer,"results/hists_%s-%s_%d-%d.root",(*(HistSets->GetIt()))->GetDataSourcePrefix(),_dilset.GetDilString().c_str(),_weighting,(int)_higgsMass);
	 std::cout<<buffer<<std::endl;
	 std::cout<<"(*(HistSets->GetIt()))->LoadHistSet(buffer);"<<std::endl;
         (*(HistSets->GetIt()))->LoadHistSet(buffer);
	 std::cout<<"end of (*(HistSets->GetIt()))->LoadHistSet(buffer);"<<std::endl;
       }while(!HistSets->IncrementIt());
     }

     ~TDSuperAna()
     {
///       SavePlots();
       if(!(HistSets->InitiateIt()))do{
         delete (*(HistSets->GetIt()));
       }while(!HistSets->IncrementIt());
       
       if(!InitiateIt()) do{
	 (*SubAnaIt)->MyTempDelete();
	 //delete ((mtuple*)(*SubAnaIt)); //do what I can. the leak seems to be caused by using the newer root files.
	 //delete (*SubAnaIt); //memory leak.
       }while(!IncrementIt());
 
       delete HistSets;
       fitinfo.close();
       _anaFile<<std::endl;
       //_anaFile<<"////////////// End of " + _dilset.GetDilString() + " details /////////////"<<std::endl;
       //_anaFile<<std::endl;
       delete CanvMan1;
       _cutTableTexFile<<"        \\end{document}"<<std::endl;
       _cutTableTexFile.close();
       _cutTableFile.close();
       //delete _SecondaryHistSets;
     }
   
     void ScaleHistSets()
     {
/*     double MCMuZpeakHeight = HistSets->GetHistSet("Zmumu")->GetAnaHist("DileptonMass")->GetHist()->GetMaximum();
       double MCElZpeakHeight = HistSets->GetHistSet("Zee")->GetAnaHist("DileptonMass")->GetHist()->GetMaximum();
       double muonZPeakHeight = HistSets->GetHistSet("muon_data")->GetAnaHist("DileptonMass")->GetHist()->GetMaximum();
       double eleZPeakHeight = HistSets->GetHistSet("electron_data")->GetAnaHist("DileptonMass")->GetHist()->GetMaximum(); */

      /* double MCMuZpeakHeight = HistSets->GetHistSet("Zmumu")->GetAnaHist("DileptonMass")->GetHist()->Integral();   //scales zee+zmumu mc to data integrals
       double MCElZpeakHeight = HistSets->GetHistSet("Zee")->GetAnaHist("DileptonMass")->GetHist()->Integral();
       double muonZPeakHeight = HistSets->GetHistSet("muon_data")->GetAnaHist("DileptonMass")->GetHist()->Integral();
       double eleZPeakHeight = HistSets->GetHistSet("electron_data")->GetAnaHist("DileptonMass")->GetHist()->Integral();

       std::cout<<MCMuZpeakHeight<<" "<<muonZPeakHeight<<" "<<MCElZpeakHeight<<" "<<eleZPeakHeight<<std::endl;
 
       double data =eleZPeakHeight+muonZPeakHeight;
       double MC = MCElZpeakHeight+MCMuZpeakHeight;
       if(MC==0.0 || MC!=MC)return;
       if(data==0.0 || data!=data)return;
      
       HistSets->GetHistSet("Zmumu")->Scale(data/MC);
       HistSets->GetHistSet("Zee")->Scale(data/MC);*/ //TD scaling code used before                                 //scales zee+zmumu mc to data integrals
 
/*       double lumcorrd=374.0/(3108482*0.51/(355.2*1.4));
       double lumcorrm=378.0/(3216951*0.51/(355.2*1.4));
       if(lumcorrd!=lumcorrd)lumcorrd=0.0;
       if(lumcorrm!=lumcorrm)lumcorrm=0.0;

       HistSets->GetHistSet("Zmumu")->Scale(lumcorrd);
       HistSets->GetHistSet("Zee")->Scale(lumcorrm);*/
       //GoodRunv2 lumquery;
       //float lum_val=lumquery.Loop(&GoodRuns)*1.019;
       TDRunInfoSet lumins("./v13goodrunquery/allreqruns.list_OffLumi_bhmu.dat");
       TDRunSet dataruns;
       if(_stn)dataruns=(*(_hwwlvlvGoodRun13.GetRunSet(_dilType.Int2Type(0))));//-_runstoignore;
       if(!_stn)dataruns=(*(_hwwlvlvGoodRun.GetRunSet(_dilType.Int2Type(0))));//-_runstoignore;
       dataruns.SetMinMax(_minRun, _maxRun);
	   //186598);
       //dataruns.SetMin(186599);
       TDRunSet mcruns;
       if(_stn)mcruns=(*(_hwwlvlvGoodRun13.GetRunSet(_dilType.Int2Type(0))))-_runstoignore;
       if(!_stn)mcruns=(*(_hwwlvlvGoodRun.GetRunSet(_dilType.Int2Type(0))))-_runstoignore;
       mcruns.SetMinMax(_minRun, _maxRun);
       //mcruns.SetMax(186598);
       //mcruns.SetMin(186599);
       double lum_val=(lumins.GetSumTDRunInfo(dataruns)._goodobsp)*1.019/1000;	  
       std::cout<<"++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ lumval used for scaling: "<<lum_val<<std::endl;
       std::cout<<"++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ lumval dilweight used for scaling: "<<_dilWeight->GetLum()/1000<<std::endl;
       //float lum_val=360;
       if(GetSubAna("Zee")->GetInputEvents()==0)
       {
//	 if(_mircea) HistSets->GetHistSet("Zee")->Scale(lum_val/(3164279/(17354*1.4))); //TD use Zee for the Zll background in the mircea settings.
	 if(_mircea) HistSets->GetHistSet("Zee")->Scale(lum_val/(2965289/(24300*0.029))); //TD peter's values - where do these come from?
	 else HistSets->GetHistSet("Zee")->Scale(lum_val/(3108482/(355.2*1.4)));
	 std::cout<<"old Zee scale:"<<lum_val/(3108482/(355.2*1.4))<<std::endl;
       }else{
	 //std::cout<<"++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"<<GetSubAna("Zee")->GetInputEvents(mcruns)<<std::endl;
	 //HistSets->GetHistSet("Zee")->Scale(lum_val/((GetSubAna("Zee")->GetInputEvents(mcruns))/(355.2*1.4))); //TD using individually weighted events, this work is done event by event.
	 HistSets->GetHistSet("Zee")->Scale(lum_val/((GetSubAna("Zee")->GetInputEvents())/(355.2*1.4)));
	 std::cout<<"old Zee scale:"<<lum_val/((GetSubAna("Zee")->GetInputEvents())/(355.2*1.4))<<std::endl;
       }
       std::cout<<"new scale:"<<GetSubAna("Zee")->GetScaleforNormalisation()<<std::endl;
       if(GetSubAna("Zmumu")->GetInputEvents()==0)
       {
	 HistSets->GetHistSet("Zmumu")->Scale(lum_val/(3216951/(355.2*1.4)));
	 std::cout<<"old Zmumu scale:"<<lum_val/(3216951/(355.2*1.4))<<std::endl;
       }else{
	 //std::cout<<"++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"<<GetSubAna("Zmumu")->GetInputEvents(mcruns)<<std::endl;
	 //HistSets->GetHistSet("Zmumu")->Scale(lum_val/((GetSubAna("Zmumu")->GetInputEvents(mcruns))/(355.2*1.4))); //TD using individually weighted events, this work is done event by event.
	 HistSets->GetHistSet("Zmumu")->Scale(lum_val/((GetSubAna("Zmumu")->GetInputEvents())/(355.2*1.4)));
	 std::cout<<"old Zmumu scale:"<<lum_val/((GetSubAna("Zmumu")->GetInputEvents())/(355.2*1.4))<<std::endl;
       }
       std::cout<<"new scale:"<<GetSubAna("Zmumu")->GetScaleforNormalisation()<<std::endl;
       if(GetSubAna("Ztautau")->GetInputEvents()==0)
       {
	 //std::cout<<"++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"<<"Ztautau is normal"<<std::endl; 
         HistSets->GetHistSet("Ztautau")->Scale(lum_val/(3338836/(355.2*1.4)));
	 std::cout<<"old Ztautau scale:"<<lum_val/(3338836/(355.2*1.4))<<std::endl;
       }else{
	 std::cout<<"++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"<<GetSubAna("Ztautau")->GetInputEvents(mcruns)<<std::endl;
	 //HistSets->GetHistSet("Ztautau")->Scale(lum_val/((GetSubAna("Ztautau")->GetInputEvents(mcruns))/(355.2*1.4))); //TD using individually weighted events, this work is done event by event.
	 HistSets->GetHistSet("Ztautau")->Scale(lum_val/((GetSubAna("Ztautau")->GetInputEvents())/(355.2*1.4)));
	 std::cout<<"old Ztautau scale:"<<lum_val/((GetSubAna("Ztautau")->GetInputEvents())/(355.2*1.4))<<std::endl;
       }
       std::cout<<"new scale:"<<GetSubAna("Ztautau")->GetScaleforNormalisation()<<std::endl;
       //       HistSets->GetHistSet("sexo8t")->Scale(lum_val/(3164279/(17354.0*1.4)));
       if(GetSubAna("WW")->GetInputEvents()==0 || !_stn)
       {
	 HistSets->GetHistSet("WW")->Scale(lum_val/(419728/12.4));
	 std::cout<<"old WW scale:"<<(float)(lum_val)/((float)419728/(float)12.4)<<std::endl;
	 std::cout<<(float)(lum_val)<<" "<<(float)419728<<" "<<(float)12.4<<" "<<((float)419728/(float)12.4)<<" "<<(float)(lum_val)/((float)419728/(float)12.4)<<std::endl;
       }else{
	 //HistSets->GetHistSet("WW")->Scale(lum_val/((GetSubAna("WW")->GetInputEvents(mcruns))/(1.27))); //TD using individually weighted events, this work is done event by event.
	 HistSets->GetHistSet("WW")->Scale(lum_val/((GetSubAna("WW")->GetInputEvents())/(1.27)));
	 std::cout<<"old WW scale input events:"<<lum_val/((GetSubAna("WW")->GetInputEvents())/(1.27))<<std::endl;
       }
       std::cout<<"new scale:"<<GetSubAna("WW")->GetScaleforNormalisation()<<" xs: "<<GetSubAna("WW")->GetXS()<<" input events "<<GetSubAna("WW")->GetInputEvents()<<std::endl;
       if(GetSubAna("WZ")->GetInputEvents()==0 || !_stn)
       {
	 HistSets->GetHistSet("WZ")->Scale(lum_val/(407438/3.65));
	 std::cout<<"old WZ scale:"<<lum_val/(407438/3.65)<<std::endl;
       }else{
	 //HistSets->GetHistSet("WZ")->Scale((lum_val*0.76)/((GetSubAna("WW")->GetInputEvents(mcruns))/0.365)); //TD using individually weighted events, this work is done event by event.
	 HistSets->GetHistSet("WZ")->Scale((lum_val*0.76)/((GetSubAna("WZ")->GetInputEvents())/0.365));
	 std::cout<<"old WZ scale:"<<(lum_val*0.76)/((GetSubAna("WZ")->GetInputEvents())/0.365)<<std::endl;
       }
       std::cout<<"new scale:"<<GetSubAna("WZ")->GetScaleforNormalisation()<<std::endl;
       if(GetSubAna("ZZ")->GetInputEvents()==0 || !_stn)
       {
	 HistSets->GetHistSet("ZZ")->Scale(lum_val/(397086/1.39));
	 std::cout<<"old ZZ scale:"<<lum_val/(397086/1.39)<<std::endl;
       }else{
	 //HistSets->GetHistSet("ZZ")->Scale((lum_val*0.23)/((GetSubAna("ZZ")->GetInputEvents(mcruns))/2.01)); //TD using individually weighted events, this work is done event by event.
	 HistSets->GetHistSet("ZZ")->Scale((lum_val*0.23)/((GetSubAna("ZZ")->GetInputEvents())/2.01));
	 std::cout<<"old ZZ scale:"<<(lum_val*0.23)/((GetSubAna("ZZ")->GetInputEvents())/2.01)<<std::endl;
       }
       std::cout<<"new scale:"<<GetSubAna("ZZ")->GetScaleforNormalisation()<<std::endl;
       //HistSets->GetHistSet("ktope")->Scale(lum_val/(170746/22.35));
       //HistSets->GetHistSet("ktopm")->Scale(lum_val/(196681/22.35));
       if(GetSubAna("ttbar")->GetInputEvents()==0 || !_stn)
       {
	 HistSets->GetHistSet("ttbar")->Scale(lum_val/(1142339/6.7));
	 std::cout<<"old ttbar scale:"<<lum_val/(1142339/6.7)<<std::endl;
       }else{
	 //HistSets->GetHistSet("ttbar")->Scale(lum_val/(((GetSubAna("ttbar")->GetInputEvents(mcruns))/6.7))); //TD using individually weighted events, this work is done event by event.
	 HistSets->GetHistSet("ttbar")->Scale(lum_val/(((GetSubAna("ttbar")->GetInputEvents())/6.7)));
	 std::cout<<"old ttbar scale:"<<lum_val/(((GetSubAna("ttbar")->GetInputEvents())/6.7))<<std::endl;
       }
       std::cout<<"new scale:"<<GetSubAna("ttbar")->GetScaleforNormalisation()<<std::endl;
       if(GetSubAna("wgamma")->GetInputEvents()==0 || !_stn)
       {
	 HistSets->GetHistSet("wgamma")->Scale(lum_val/(367427/89.4));
	 std::cout<<"old wgamma scale:"<<lum_val/(367427/89.4)<<std::endl;
       }else{
	 //HistSets->GetHistSet("wgamma")->Scale(lum_val/(((GetSubAna("wgamma")->GetInputEvents(mcruns))/(21.48*1.34)))); //TD using individually weighted events, this work is done event by event.
	 HistSets->GetHistSet("wgamma")->Scale(lum_val/(((GetSubAna("wgamma")->GetInputEvents())/(2*21.48*1.34)))); //TD each smaples has 21.48*1.34, and there are two of them => *2.
	 std::cout<<"old wgamma scale:"<<lum_val/(((GetSubAna("wgamma")->GetInputEvents())/(2*21.48*1.34)))<<std::endl;
       }
       std::cout<<"new scale:"<<GetSubAna("wgamma")->GetScaleforNormalisation()<<std::endl;
       //fakes:
       float fakes[11]={ 0., 1.595, 1.672, 1.780, 1.526, 1.507, 1.334, 1.242, 0.884,
	                     0.800, 0.722};
         float wgamT[11]={ 0,  2.837, 2.895, 2.217, 1.967, 1.477, 1.142, 0.969, 0.816,
	                       0.654, 0.477};
	 int arraynum=(int)(((double)((int(_higgsMass))-100))/10.0);
       //HistSets->GetHistSet("wgamma")->Scale(1+fakes[arraynum]/wgamT[arraynum]);  //mass 160==[6];//no more scaling of the Wgamma for fakes!
       //wexo8w

       float events=0; 
       float xs=0;
       if(!_stn)
       switch((int)_higgsMass)//gen 5 top mc if using hww gen6 mc use below!!
       {
	 case(110):events=1143420;xs=0.004;break;
	 case(120):events=1137034;xs=0.009;break;
	 case(130):events=1199597;xs=0.017;break;
	 case(140):events=1136447;xs=0.022;break;
	 case(150):events=1162684;xs=0.025;break;
	 case(160):events=1118660;xs=0.028;break;
	 case(170):events=1106451;xs=0.025;break;
	 case(180):events=1153345;xs=0.020;break;
	 case(190):events=1093662;xs=0.014;break;
	 case(200):events=1092527;xs=0.011;break;
	 default:events=1118660;xs=0.028;break;
       }
       
       if(_stn)
       switch((int)_higgsMass)//Gen 6
       {
         case(110):events=1339174;xs=0.006;break;
         case(120):events=1138940;xs=0.014;break;
         case(130):events=1305478;xs=0.024;break;
         case(140):events=1263099;xs=0.033;break;
         case(150):events=476447;xs=0.038;break;
         case(160):events=1339351;xs=0.041;break;
         case(170):events=881884;xs=0.036;break;
         case(180):events=1091842;xs=0.029;break;
         case(190):events=1031206;xs=0.020;break;
         case(200):events=947534;xs=0.016;break;
         default:events=1339351;xs=0.041;break;
       }

       if(xs!=GetHwwXS((int)_higgsMass, _stn) || events!=GetHwwEvents((int)_higgsMass, _stn))
	 std::cout<<"getHwwEvents or xs is not working"<<std::endl;
       
       if(GetSubAna("higgsww")->GetInputEvents()==0 || !_stn)
       {
         HistSets->GetHistSet("higgsww")->Scale(lum_val/(events/xs));
	 std::cout<<"old higgsww scale:"<<lum_val/(events/xs)<<std::endl;
       }else{
	 //HistSets->GetHistSet("higgsww")->Scale(lum_val/(((GetSubAna("higgsww")->GetInputEvents(mcruns))/xs))); //TD using individually weighted events, this work is done event by event.
	 HistSets->GetHistSet("higgsww")->Scale(lum_val/(((GetSubAna("higgsww")->GetInputEvents())/xs)));
	 std::cout<<"old higgsww scale:"<<lum_val/(((GetSubAna("higgsww")->GetInputEvents())/xs))<<std::endl;
       }
       std::cout<<"new scale:"<<GetSubAna("higgsww")->GetScaleforNormalisation()<<std::endl;
     }

     void InitSubAnas()
     {
       if(!(HistSets->InitiateIt()))do{
         (*(HistSets->GetIt()))->SetDebug(false);
       }while(!HistSets->IncrementIt());

       bool matchedhistset;
       if(!InitiateIt())do{
         matchedhistset=false;
         DebugOut("here");
         if(!(HistSets->InitiateIt()))do{
           if(_debug)std::cout<<(*SubAnaIt)->GetDataSource()<<(*(HistSets->GetIt()))->GetDataSource()<<std::endl;
           if((std::string)((*SubAnaIt)->GetDataSource())==(std::string)((*(HistSets->GetIt()))->GetDataSource()))
           {
             (*SubAnaIt)->SetHistSet(*(HistSets->GetIt()));
             (*SubAnaIt)->InitHists();
             matchedhistset=true;
             DebugOut("corresponding histset found");
           }
         }while(!(HistSets->IncrementIt()) || !matchedhistset);
       }while(!IncrementIt());
       //DebugOut("end of histset init");
       float lowestlum=9999999.0;
       if(!InitiateIt())do{
	// std::cout<<"in loop: "<<(*SubAnaIt)->GetDataSource()<<std::endl;
         if(!((*SubAnaIt)->isData()))
         {
	   //std::cout<<"in loop: "<<(*SubAnaIt)->GetDataSource()<<std::endl;
	   //std::cout<<(*SubAnaIt)->GetInputEvents()<<" "<<(*SubAnaIt)->GetXS()<<std::endl;
           if((*SubAnaIt)->GetInputEvents()/(*SubAnaIt)->GetXS() < lowestlum)
           {
             lowestlum=(*SubAnaIt)->GetInputEvents()/(*SubAnaIt)->GetXS();
             //std::cout<<"lowest luminosity is: "<<lowestlum<<" from "<<(*SubAnaIt)->GetDataSource()<<std::endl;
           }
           //sampletreelist.Add((*SubAnaIt)->GetDescVar());
         }
       }while(!IncrementIt());
       //std::cout<<"lowest luminosity is: "<<lowestlum<<std::endl;
       if(!InitiateIt())do{
	 (*SubAnaIt)->Set_lowestlumback(lowestlum);
       }while(!IncrementIt());
//       Set_lowestlumback
       DebugOut("end of histset init");
     }
 
     void RunSubAnas()
     {
       //HistSets->InitiateIt();do{
       //  (*(HistSets->GetIt()))->PrintValues();
       //}while(!HistSets->IncrementIt());

       DebugOut("RunSubAnas");

       if(!InitiateIt())do{
         (*SubAnaIt)->Loop();
       }while(!IncrementIt());

        DebugOut("before EndAna");

       if(!InitiateIt()) do{
         (*SubAnaIt)->EndAna();
       }while(!IncrementIt());
       
       DebugOut("after EndAna");
     }
   
     void DrawPlots(char * options="")
     {
       DrawUnscaledPlots(options);
       ScaleHistSets();
       DrawScaledPlots(options);
     }

     void DrawUnscaledPlots(char * options="")
     {
       char buffer[500];
       //HistSets->SetInnateColors();
       //HistSets->DrawPlots(options);
       //HistSets->StackPlot(options);
       //HistSets->SetSimpleColors();

       /*HistSets->InitiateIt();do{DebugOut("after draw plots");
	 (*(HistSets->GetIt()))->PrintValues();
	 }while(!HistSets->IncrementIt());*/
       //XSroutine();

       TDMassPeakFit testmasspeak(CanvMan1, fitinfo, _debug);
       TDVnEEE elecmean(0);

       if(!(HistSets->InitiateIt()))do{
	 //testmasspeak.Massfitter( *((*(HistSets->GetIt()))->GetAnaHist("DileptonMass")->GetHist()) );
	 //testmasspeak.DoubleMassfitter( *((*(HistSets->GetIt()))->GetAnaHist("DileptonMass")->GetHist()) );
	 //testmasspeak.RangeMassfitter( *((*(HistSets->GetIt()))->GetAnaHist("DileptonMass")->GetHist()), 87.0, 95.0);
	 testmasspeak.RangeMassfitter(*((TH1F*)((*(HistSets->GetIt()))->GetAnaHist("DileptonMass")->GetHist())), 87.0, 95.0);
	 sprintf(buffer,"%s_%d_%s ", _dilset.GetDilString().c_str(), _weighting, (*(HistSets->GetIt()))->GetName());
	 _excelMassFile<<buffer;
	 elecmean=testmasspeak.GetVnEEEMean();
	 elecmean.PrintValues();
	 ExcelDataOut(elecmean, _excelMassFile);
	 _excelMassFile<<std::endl;
       }while(!(HistSets->IncrementIt()));       

       /*testmasspeak.RangeMassfitterSame(*(HistSets->GetHistSet("Zee")->GetAnaHist("DileptonMass")->GetHist()), //Ricks plots
	*(HistSets->GetHistSet("HwwStn_zewk6ae")->GetAnaHist("DileptonMass")->GetHist()),
	87.0, 95.0);
	testmasspeak.RangeMassfitterSame(*(HistSets->GetHistSet("HwwStn_zewk6ae")->GetAnaHist("DileptonMass")->GetHist()), //Ricks plots
	*(HistSets->GetHistSet("Zee")->GetAnaHist("DileptonMass")->GetHist()),
	87.0, 95.0);
	testmasspeak.RangeMassfitterSame(*(HistSets->GetHistSet("Zmumu")->GetAnaHist("DileptonMass")->GetHist()), //Ricks plots
	*(HistSets->GetHistSet("HwwStn_zewk6m_v1_88")->GetAnaHist("DileptonMass")->GetHist()),
	87.0, 95.0);
	testmasspeak.RangeMassfitterSame(*(HistSets->GetHistSet("HwwStn_zewk6m_v1_88")->GetAnaHist("DileptonMass")->GetHist()), //Ricks plots
	*(HistSets->GetHistSet("Zmumu")->GetAnaHist("DileptonMass")->GetHist()),
	87.0, 95.0);
	testmasspeak.RangeMassfitterSame(*(HistSets->GetHistSet("Zmumu")->GetAnaHist("DileptonMass")->GetHist()), //Ricks plots
	*(HistSets->GetHistSet("HwwStn_zewk6m_v1_90")->GetAnaHist("DileptonMass")->GetHist()),
	87.0, 95.0);
	testmasspeak.RangeMassfitterSame(*(HistSets->GetHistSet("HwwStn_zewk6m_v1_90")->GetAnaHist("DileptonMass")->GetHist()), //Ricks plots
	*(HistSets->GetHistSet("Zmumu")->GetAnaHist("DileptonMass")->GetHist()),
	87.0, 95.0);*/

       //ScaleHistSets();
       //HistSets->SynchronousPlot(options);
     }
     void DrawScaledPlots(char * options="")
     {
       HistSets->SynchronousPlot(options);
       //HistSets->CumulativeStackPlot(options);
       //HistSets->StackPlot(options);
       //CanvMan1->NextCanvasSlot();
       //gPad->SetLogy();
       HistSets->CumuAnaStackPlot(options);
       HistSets->CumuAnaStackPlot(options, true);
      

       //CanvMan1->NextCanvasSlot();
       
       //HistSets->GetHistSet("Zmumu")->GetAnaHist("DileptonMass")->GetHist()->Draw("HIST");
       //HistSets->GetHistSet("muon_data")->GetAnaHist("DileptonMass")->GetHist()->Draw("ESAME");
       //HistSets->GetHistSet("Zmumu")->GetAnaHist("DileptonMass")->GetHist()->Draw("SAMEAXIS");

       //CanvMan1->NextCanvasSlot();

       //HistSets->GetHistSet("Zee")->GetAnaHist("DileptonMass")->GetHist()->Draw("HIST");
       //HistSets->GetHistSet("electron_data")->GetAnaHist("DileptonMass")->GetHist()->Draw("ESAME");
       //HistSets->GetHistSet("Zee")->GetAnaHist("DileptonMass")->GetHist()->Draw("SAMEAXIS");


       //testmasspeak.RangeMassfitterSame(*(HistSets->GetHistSet("electron_data")->GetAnaHist("DileptonMass")->GetHist()),
        //                               *(HistSets->GetHistSet("Zee")->GetAnaHist("DileptonMass")->GetHist()),
         //                              87.0, 95.0);
       //testmasspeak.RangeMassfitterSame(*(HistSets->GetHistSet("muon_data")->GetAnaHist("DileptonMass")->GetHist()),
       //                                *(HistSets->GetHistSet("Zmumu")->GetAnaHist("DileptonMass")->GetHist()),
       //                                87.0, 95.0);
       if(!InitiateIt()) do{
         if((*SubAnaIt)->GetAnaCutSet()->GetCutAnaCut()->Gettried()!=0)
         {
           _cutTableFile<<_dilset.GetDilString()<<" "<<(*SubAnaIt)->GetDataSource()<<" ";
           (*SubAnaIt)->GetAnaCutSet()->GetCutAnaCut()->PrintCumulativeTreeinBlockForm(_cutTableFile);
           _cutTableFile<<std::endl;
         }
       }while(!IncrementIt());
       
       CanvMan1->NextCanvasSlot();
       //Dphiplots2();
       Dphiplots();
       MuonEtaPhiPlots();
       //std::cout<<"++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ lumval dilweight used for scaling: "<<_dilWeight->GetLum()<<std::endl;
     }

     void DrawOneCumuAnaPlot(TDHistSetSet * hisetset, const char * plotname, bool Log=false, bool title=true, bool plotdata=true, bool legend=true)
     {
       char buffer[200];
       sprintf(buffer,"results/individual_%s_%s_%d_%d_plot.eps",plotname,_dilset.GetDilString().c_str(),_weighting,(int)_higgsMass);
       TDCanvMan DrawRickPlotscanvman(buffer, "plot", 1, 1, 400, 400);
       hisetset->SetCanvMan(&DrawRickPlotscanvman);
       hisetset->Individual_CumuAnaStackPlot(plotname, "HIST", Log, title, plotdata, legend);
     }

     void DrawRickPlots()
     {
       /*if(!(HistSets->InitiateIt()))do{
         if(!((*(HistSets->GetIt()))->InitiateIt()))do{
           (*((*(HistSets->GetIt()))->GetAnaIt()))->GetHist()->SetTitle("");
         }while(!((*(HistSets->GetIt()))->IncrementIt()));
       }while(!(HistSets->IncrementIt()));*/

       DrawOneCumuAnaPlot(HistSets, "LeadJetEt_anyN", true, false);
       DrawOneCumuAnaPlot(HistSets, "LeadJetEt_1Jet", true, false);
       DrawOneCumuAnaPlot(HistSets, "LeadJetEt_2Jet", true, false);
       DrawOneCumuAnaPlot(HistSets, "SecondJetEt_2Jet", true, false);
       DrawOneCumuAnaPlot(HistSets, "JetNumBeforeCut", true, false);
       DrawOneCumuAnaPlot(HistSets, "MissingETMag", true, false);
       DrawOneCumuAnaPlot(HistSets, "dilmassmax", true, false);
       DrawOneCumuAnaPlot(HistSets, "TotalEt", true, false);
/*       DrawOneCumuAnaPlot(HistSets, "WWmetsig", false);
       DrawOneCumuAnaPlot(HistSets, "WWdilmass", false);
       DrawOneCumuAnaPlot(HistSets, "WWmetmag", false);
       DrawOneCumuAnaPlot(HistSets, "WWsumEt",false);
       DrawOneCumuAnaPlot(HistSets, "WWlepPt",false);*/
       HistSets->SetCanvMan(CanvMan1);
     }

     void MuonEtaPhiPlots()
     {
       TLine * line;
       std::vector<std::string> etaphiplotnames;
       etaphiplotnames.push_back("CMUEtaPhi");
       etaphiplotnames.push_back("CMUPEtaPhi");
       etaphiplotnames.push_back("CMPEtaPhi");
       for(int name=0; name!=(int)etaphiplotnames.size(); name++)
       {
         CanvMan1->NextCanvasSlot();
         //(*(HistSets->GetIt()))->GetAnaHist("CMUEtaPhi")->GetHist()->Draw();
	 std::cout<<(etaphiplotnames[name]).c_str()<<std::endl;
	 HistSets->GetHistSet("Zmumu")->GetAnaHist((etaphiplotnames[name]).c_str())->Draw();
	 for(int nline=1; nline!=12; nline++)
	 {
           line = new TLine(-1.5,(3.14*2)*nline/12,1.5,(3.14*2)*nline/12);
           line->SetLineColor(kRed);
	   line->Draw();
	   /*line = new TLine(0.26,0,0.26,3.14*2);
	   line->SetLineColor(kRed);
	   line->Draw();
	   line = new TLine(-0.26,0,-0.26,3.14*2);
	   line->SetLineColor(kRed);
	   line->Draw();*/
	   line = new TLine(0.6,0,0.6,3.14*2);
	   line->SetLineColor(kRed);
	   line->Draw();
	   line = new TLine(-0.6,0,-0.6,3.14*2);
	   line->SetLineColor(kRed);
	   line->Draw();
	   //delete line;   //TD memory leaks! correct this!
	 }
       }
       char buffer[200];
       for(int name=0; name!=(int)etaphiplotnames.size(); name++)
       {
	 sprintf(buffer,"results/individual_%s_%s_%d_%d_plot.eps",(etaphiplotnames[name]).c_str(),_dilset.GetDilString().c_str(),_weighting,(int)_higgsMass);
	 TDCanvMan DrawRickPlotscanvman(buffer, "plot", 1, 1, 400, 400);
	 HistSets->SetCanvMan(&DrawRickPlotscanvman);
         DrawRickPlotscanvman.NextCanvasSlot();
	 //HistSets->GetHistSet("Zmumu")->SetDebug(true);
	 HistSets->GetHistSet("Zmumu")->GetAnaHist((etaphiplotnames[name]).c_str())->Draw("HIST BOX");
         for(int nline=1; nline!=12; nline++)
         {
           line = new TLine(-1.5,(3.14*2)*nline/12,1.5,(3.14*2)*nline/12);
           line->SetLineColor(kRed);
           line->Draw();
           line = new TLine(0.6,0,0.6,3.14*2);
           line->SetLineColor(kRed);
           line->Draw();
           line = new TLine(-0.6,0,-0.6,3.14*2);
           line->SetLineColor(kRed);
           line->Draw();
           //delete line;   //TD memory leaks! correct this!
         }
       }
       HistSets->SetCanvMan(CanvMan1);
       //HistSets->GetHistSet("Zmumu")->SetDebug(false);
     }

     void Dphiplots2()
     {
       char buffer[500];
       if(!(HistSets->GetHistSet("Zee")->InitiateIt()))do{
	 sprintf(buffer,"%s",(*(HistSets->GetHistSet("Zee")->GetAnaIt()))->GetHistName());
	 std::cout<<"Dphiplots2 hist: "<<buffer<<std::endl;
	 TDHistSet * dphiplots=new (TDHistSet)(buffer,buffer,buffer, CanvMan1);
	 TDAnaHist * temp= new TDAnaHist(**(HistSets->GetHistSet("Zee")->GetAnaIt()));
	 temp->GetHist()->SetFillColor(17);
	 temp->GetHist()->Add(HistSets->GetHistSet("Zmumu")->GetAnaHist(buffer)->GetHist());
	 temp->GetHist()->Add(HistSets->GetHistSet("Ztautau")->GetAnaHist(buffer)->GetHist());
	 dphiplots->CumulativeAddAnaHist(temp);
	 temp =  new TDAnaHist(*(HistSets->GetHistSet("wgamma")->GetAnaHist(buffer)));
	 temp->GetHist()->Add(HistSets->GetHistSet("fake")->GetAnaHist(buffer)->GetHist());
	 temp->GetHist()->SetFillColor(30);
	 dphiplots->CumulativeAddAnaHist(temp);
	 temp =  new TDAnaHist(*(HistSets->GetHistSet("WW")->GetAnaHist(buffer)));
	 temp->GetHist()->SetFillColor(38);
	 dphiplots->CumulativeAddAnaHist(temp);
	 temp =  new TDAnaHist(*(HistSets->GetHistSet("WZ")->GetAnaHist(buffer)));
	 temp->GetHist()->Add(HistSets->GetHistSet("ZZ")->GetAnaHist(buffer)->GetHist());
	 temp->GetHist()->Add(HistSets->GetHistSet("ttbar")->GetAnaHist(buffer)->GetHist());
	 temp->GetHist()->SetFillColor(4);
	 dphiplots->CumulativeAddAnaHist(temp);
	 temp =  new TDAnaHist(*(HistSets->GetHistSet("higgsww")->GetAnaHist(buffer)));
	 temp->GetHist()->SetFillColor(2);
	 dphiplots->CumulativeAddAnaHist(temp);
	 dphiplots->DrawPlotsSame("HIST");
//	 CanvMan1->NextCanvasSlot()
	 dphiplots->DrawPlotsSame("HIST", true);
	 delete dphiplots;
       }while(!(HistSets->GetHistSet("Zee")->IncrementIt()));
      CanvMan1->NextCanvasSlot();
     }

     void Dphiplots()
     {
       std::cout<<"Dphiplots1"<<std::endl;
       TDHistSetSet Dphiplotshists(CanvMan1, false);
       std::cout<<"Dphiplots2"<<std::endl; 
       TDHistSet DY(*(HistSets->GetHistSet("Zee")));
       //std::cout<<"Dphiplots2 printvalues"<<std::endl;
       //DY.PrintValues();
       std::cout<<"Dphiplots3"<<std::endl;
       DY+=(*(HistSets->GetHistSet("Zmumu")));
       std::cout<<"Dphiplots4"<<std::endl;
       DY+=(*(HistSets->GetHistSet("Ztautau")));
       std::cout<<"Dphiplots5"<<std::endl;
       DY.SetColor(17);
       //DY.SetDebug(true);
       //DY.DrawPlots("HIST");
       std::cout<<"Dphiplots6"<<std::endl;
       Dphiplotshists.AddHistSet(&DY);
       std::cout<<"Dphiplots7"<<std::endl;
       DY.SetDataSourcePrefix("Drell Yan");

       TDHistSet Wgamma(*(HistSets->GetHistSet("wgamma")));
       std::cout<<"Dphiplots8"<<std::endl;
       Wgamma+=(*(HistSets->GetHistSet("fake")));
       Wgamma.SetColor(30);
       //Wgamma.SetDebug(true);
       std::cout<<"Dphiplots9"<<std::endl;
       Dphiplotshists.AddHistSet(&Wgamma);
       std::cout<<"Dphiplots10"<<std::endl;
       //Wgamma.DrawPlots("HIST");
       
       std::cout<<"Dphiplots11"<<std::endl;
       TDHistSet WW(*(HistSets->GetHistSet("WW")));
       std::cout<<"Dphiplots12"<<std::endl;
       WW.SetColor(38);
       //WW.SetDebug(true);
       std::cout<<"Dphiplots13"<<std::endl;
       Dphiplotshists.AddHistSet(&WW);
       std::cout<<"Dphiplots14"<<std::endl;

       std::cout<<"Dphiplots15"<<std::endl;
//       (*(HistSets->GetHistSet("WZ"))).PrintValues();
       TDHistSet WZZZTT(*(HistSets->GetHistSet("WZ")));
       std::cout<<"Dphiplots16"<<std::endl;
       //WZZZTT+=(*(HistSets->GetHistSet("ZZ")));
       std::cout<<"Dphiplots17"<<std::endl;
       //WZZZTT+=(*(HistSets->GetHistSet("ttbar")));
       std::cout<<"Dphiplots18"<<std::endl;
       WZZZTT.SetColor(4);
       //WZZZTT.SetDebug(true);
       std::cout<<"Dphiplots19"<<std::endl;
       Dphiplotshists.AddHistSet(&WZZZTT);
       std::cout<<"Dphiplots20"<<std::endl;
       WZZZTT.SetDataSourcePrefix("WZ, ZZ, ttbar");

       std::cout<<"Dphiplots21"<<std::endl;
       TDHistSet HWW(*(HistSets->GetHistSet("higgsww")));
       std::cout<<"Dphiplots22"<<std::endl;
       HWW.SetColor(2);
       //HWW.SetDebug(true);
       std::cout<<"Dphiplots23"<<std::endl;
       Dphiplotshists.AddHistSet(&HWW);
       std::cout<<"Dphiplots24"<<std::endl;

       TDHistSet DATA(*(HistSets->GetHistSet("data")));
       std::cout<<"Dphiplots24.1"<<std::endl;
       //DATA+=(*(HistSets->GetHistSet("electron_data")));
       //DATA+=(*(HistSets->GetHistSet("plug_data")));
       //DATA.SetDebug(true);
       DATA.SetColor(6);
       std::cout<<"Dphiplots24.2"<<std::endl;
       Dphiplotshists.AddHistSet(&DATA);
       
       std::cout<<"Dphiplots25"<<std::endl;
       Dphiplotshists.CumuAnaStackPlot("HIST");
       std::cout<<"Dphiplots26"<<std::endl;
       CanvMan1->NextCanvasSlot();
       Dphiplotshists.CumuAnaStackPlot("HIST",true);
       CanvMan1->NextCanvasSlot();
       DrawOneCumuAnaPlot(&Dphiplotshists, "dilDphi", false, true, false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "dilDphiee", false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "dilDphiemu", false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "dilDphimumu", false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "WWmetsig", false, false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "WWdilmass", false, false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "WWmetmag", false, false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "WWsumEt",false, false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "WWlepPt",false, false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "MVA_Likelihood", true, false, true, false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "MVA2_Likelihood", true, false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "MVA_LikelihoodD", true, false, true, false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "MVA2_LikelihoodD", true, false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "MVA_PDERS", true, false, true, false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "MVA2_PDERS", true, false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "MVA_HMatrix", true, false, true, false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "MVA2_HMatrix", true, false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "MVA_Fisher", true, false, true, false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "MVA2_Fisher", true, false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "MVA_BDT", true, false, true, false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "MVA2_BDT", true, false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "MVA_CFMlpANN", true, false, true, false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "MVA2_CFMlpANN", true, false);
       DrawOneCumuAnaPlot(&Dphiplotshists, "MVA2_CFMlpANN_cut2_7", true, false);
       
     }

     TTree * MakeAnaTree()
     {
       TList sampletreelist;
       TTree * CumuTree;
       TTree * TempTree;
       TTree * TempTree2;
       float lowestlum=99999999;
       std::cout<<"starting lowest lumcalc"<<std::endl;
       if(!InitiateIt())do{
	 if(!((*SubAnaIt)->isData()))
	 {
           if((*SubAnaIt)->GetInputEvents()/(*SubAnaIt)->GetXS() < lowestlum)
	   {
	     lowestlum=(*SubAnaIt)->GetInputEvents()/(*SubAnaIt)->GetXS();
             std::cout<<"lowest luminosity is: "<<lowestlum<<" from "<<(*SubAnaIt)->GetDataSource()<<std::endl;
	   }
	   //sampletreelist.Add((*SubAnaIt)->GetDescVar());
	 }
       }while(!IncrementIt());
       std::cout<<"finished lowest lumcalc"<<std::endl;
       bool first=true;
       if(!InitiateIt())do{
         if(!((*SubAnaIt)->isData()))
         {
           float thislum=(*SubAnaIt)->GetInputEvents()/(*SubAnaIt)->GetXS();
	   std::cout<<" luminosity is: "<<thislum<<" from "<<(*SubAnaIt)->GetDataSource()<<std::endl;
           char ratio[200];
           sprintf(ratio,"rand<%1.6f", lowestlum/thislum);
	   std::cout<<"ratio: "<<ratio<<std::endl;
	   if(!first)
	   {
	     TempTree2=CumuTree->CloneTree();
	     std::cout<<"tree2 cloning ok"<<std::endl;
	     if(strcmp((*SubAnaIt)->GetDataSource(),"higgsww")==0)TempTree=(*SubAnaIt)->GetDescVar()->CloneTree();
	     else TempTree=(*SubAnaIt)->GetDescVar()->CopyTree(ratio);
	     std::cout<<"tree copying ok"<<std::endl;
	     std::cout<<"TempTree->GetEntries()"<<TempTree->GetEntries()<<std::endl;
	     std::cout<<"TempTree2->GetEntries()"<<TempTree2->GetEntries()<<std::endl;
	     delete CumuTree;
	     std::cout<<"deleted cumu tree ok"<<std::endl;
	     sampletreelist.Add(TempTree);
	     sampletreelist.Add(TempTree2);
	     std::cout<<"sampletreelist.GetEntries()";
	     std::cout<<sampletreelist.GetEntries()<<std::endl;
	     CumuTree=(TTree::MergeTrees(&sampletreelist));//TD memory leak, I think.. (maybe not)
	   }
	   if(first)CumuTree=(*SubAnaIt)->GetDescVar()->CopyTree(ratio);
	   std::cout<<"CumuTree->GetEntries()";
	   std::cout<<CumuTree->GetEntries()<<std::endl;
	   std::cout<<"remade cumu tree ok"<<std::endl;
           sampletreelist.Clear();
           if(!first)
	   {
	     delete TempTree;
             delete TempTree2;
	   }
	   first=false;
         }
       }while(!IncrementIt());
       return CumuTree;
     }

     TTree * MakeAnaTree2()
     {
       TList sampletreelist;
       TTree * CumuTree;
       TTree * TempTree;
       TTree * TempTree2;
       float lowestlum=99999999;
       std::cout<<"starting lowest lumcalc"<<std::endl;
       if(!InitiateIt())do{
	 if(!((*SubAnaIt)->isData()))
	 {
           if((*SubAnaIt)->GetInputEvents()/(*SubAnaIt)->GetXS() < lowestlum)
	   {
	     lowestlum=(*SubAnaIt)->GetInputEvents()/(*SubAnaIt)->GetXS();
             std::cout<<"lowest luminosity is: "<<lowestlum<<" from "<<(*SubAnaIt)->GetDataSource()<<std::endl;
	   }
	   //sampletreelist.Add((*SubAnaIt)->GetDescVar2());
	 }
       }while(!IncrementIt());
       std::cout<<"finished lowest lumcalc"<<std::endl;
       bool first=true;
       if(!InitiateIt())do{
         if(!((*SubAnaIt)->isData()))
         {
           float thislum=(*SubAnaIt)->GetInputEvents()/(*SubAnaIt)->GetXS();
	   std::cout<<" luminosity is: "<<thislum<<" from "<<(*SubAnaIt)->GetDataSource()<<std::endl;
           char ratio[200];
           sprintf(ratio,"rand<%1.6f", lowestlum/thislum);
	   std::cout<<"ratio: "<<ratio<<std::endl;
	   if(!first)
	   {
	     TempTree2=CumuTree->CloneTree();
	     std::cout<<"tree2 cloning ok"<<std::endl;
	     if(strcmp((*SubAnaIt)->GetDataSource(),"higgsww")==0)TempTree=(*SubAnaIt)->GetDescVar2()->CloneTree();
	     else TempTree=(*SubAnaIt)->GetDescVar2()->CopyTree(ratio);
	     std::cout<<"tree copying ok"<<std::endl;
	     std::cout<<"TempTree->GetEntries()"<<TempTree->GetEntries()<<std::endl;
	     std::cout<<"TempTree2->GetEntries()"<<TempTree2->GetEntries()<<std::endl;
	     delete CumuTree;
	     std::cout<<"deleted cumu tree ok"<<std::endl;
	     sampletreelist.Add(TempTree);
	     sampletreelist.Add(TempTree2);
	     std::cout<<"sampletreelist.GetEntries()";
	     std::cout<<sampletreelist.GetEntries()<<std::endl;
	     CumuTree=(TTree::MergeTrees(&sampletreelist));//TD memory leak, I think.. (maybe not)
	   }
	   if(first)CumuTree=(*SubAnaIt)->GetDescVar2()->CopyTree(ratio);
	   std::cout<<"CumuTree->GetEntries()";
	   std::cout<<CumuTree->GetEntries()<<std::endl;
	   std::cout<<"remade cumu tree ok"<<std::endl;
           sampletreelist.Clear();
           if(!first)
	   {
	     delete TempTree;
             delete TempTree2;
	   }
	   first=false;
         }
       }while(!IncrementIt());
       return CumuTree;
     }

     void NN()
     {
       char buffer[200];
       std::cout<<"starting NN now"<<std::endl;

       TList sampletreelist2;
       sampletreelist2.Add( GetSubAna("Zee")->GetDescVar());
       sampletreelist2.Add( GetSubAna("Zmumu")->GetDescVar());
       sampletreelist2.Add( GetSubAna("Ztautau")->GetDescVar());
       sampletreelist2.Add( GetSubAna("higgsww")->GetDescVar());
       TTree * totaltree=TTree::MergeTrees(&sampletreelist2);
       std::cout<<"made tree"<<std::endl;
       
      /* TEventList train("train");
       TEventList test("test");
       totaltree->Draw(">>train","rand<0.05 && datasource==0");
       totaltree->Draw(">>+train","rand<0.002 && rand>0.001 && datasource==0");
       totaltree->Draw(">>+train","rand<0.003 && rand>0.002 && datasource==0");  
       totaltree->Draw(">>+train","rand<0.004 && rand>0.003 && datasource==0");
       totaltree->Draw(">>+train","rand<0.005 && rand>0.004 && datasource==0");
       std::cout<<"ntraining events: "<<train.GetN()<<std::endl;
       totaltree->Draw(">>+train","rand<0.05 && datasource==1");
       
       totaltree->Draw(">>test","rand>0.95&& datasource==0");
       totaltree->Draw(">>+test","rand>0.998 && rand<0.999 && datasource==0");
       totaltree->Draw(">>+test","rand>0.997 && rand<0.998 && datasource==0");
       totaltree->Draw(">>+test","rand>0.996 && rand<0.997 && datasource==0");
       totaltree->Draw(">>+test","rand>0.995 && rand<0.996 && datasource==0");
       std::cout<<"ntest events: "<<test.GetN()<<std::endl;
       totaltree->Draw(">>+test","rand>0.95 && datasource==1");*/
       sprintf(buffer,"results/NN1_%d.root", (int)_higgsMass);
       //TFile* f= new TFile("results/NNtest.root","recreate");
       TFile* f= new TFile(buffer,"recreate");
       //totaltree->SetEventList(&train);
       totaltree->Write("",TObject::kOverwrite);
       f->Close();
       delete totaltree;

       /*//totaltree->Scan("*", "","",10);
       //totaltree->Scan("*", "");
       //totaltree->Scan("*", "rand<0.01 && datasource==1","",10);
       //totaltree->Scan("*", "rand>0.999 && datasource==0","",10);
       //totaltree->Scan("*", "rand>0.99 && datasource==1","",10);*/
       TList sampletreelist;
       //if(!InitiateIt())do{
//       if(!((*SubAnaIt)->isData()))
//       {
 //          sampletreelist.Add((*SubAnaIt)->GetDescVar());
//       }
  //     }while(!IncrementIt());
       sampletreelist.Add( GetSubAna("WW")->GetDescVar2());
       sampletreelist.Add( GetSubAna("higgsww")->GetDescVar2());
       TTree * totaltree2=TTree::MergeTrees(&sampletreelist);

       std::cout<<"made tree2"<<std::endl;
       sprintf(buffer,"results/NN2_%d.root", (int)_higgsMass);
       //TFile* f= new TFile("results/NNtest.root","recreate");
       TFile* f2= new TFile(buffer,"recreate");
       //totaltree->SetEventList(&train);
       totaltree2->Write("",TObject::kOverwrite);
       f2->Close();
       //delete totaltree2;
      
       double sigent=GetSubAna("higgsww")->GetDescVar2()->GetEntries();
       double backent=GetSubAna("WW")->GetDescVar2()->GetEntries();
       sprintf(buffer,"weight*((type==0)*%f)", sigent/backent);
      /* std::cout<<"made evtlists"<<std::endl;*/
       //TMultiLayerPerceptron myNN("LepAEn,LepBEn,dilmass,metmag,ChargeProduct,addEt,dPhiMetLJ,ntightjets,jet1_Et,jet2_Et:22:10:5:datasource",totaltree,&train,&test);
  //     char buffer[200];
//       sprintf(buffer,"\"((Entry$+1)%2) && Entry<%d \"",(GetSubAna("WW")->GetDescVar()->GetEntries())*2);
       //TMultiLayerPerceptron myNN("LepAEn,LepBEn,dilmass,metmag,ChargeProduct,addEt,dPhiMetLJ,ntightjets,jet1_Et,jet2_Et,dildphi,metsig,lepR,diltype:15:14:type","weight+((type==0)*20)",totaltree2,"Entry$%2","((Entry$+1)%2)");
       //myNN.SetEventWeight("weight");
       //std::cout<<"made NN; training..."<<std::endl;
       //std::cout<<"ntraining events: "<<train.GetN()<<std::endl;
       //std::cout<<"ntest events: "<<test.GetN()<<std::endl;
       //myNN.Train(200,"text, graph, update=1");
       //sprintf(buffer,"./results/TMLP_%d.wgt", (int)_higgsMass);
       //myNN.DumpWeights(buffer);
       //CanvMan1->NextCanvasSlot();
       //TMLPAnalyzer ana(myNN);
       //ana.GatherInformations();
       //ana.DrawNetwork(0,"type==1","type==0");
       //myNN.Draw();
       //myNN.DumpWeights("NNout.txt");
       std::cout<<"done training..."<<std::endl;
       delete totaltree2;
     }

     void TMVAnalysis()
     {
       std::cout<<"TMV Analysis"<<std::endl;//not using at the moment - getit right sometime, though..
       /*// ---------------------------------------------------------------
       // choose MVA methods to be trained + tested
       bool Use_Cuts            = 0;
       bool Use_Likelihood      = 1;
       bool Use_LikelihoodD     = 1;
       bool Use_PDERS           = 0;
       bool Use_HMatrix         = 1;
       bool Use_Fisher          = 1;
       bool Use_CFMlpANN        = 1;
       bool Use_TMlpANN         = 0;
       bool Use_BDT_GiniIndex   = 0; // default BDT method
       bool Use_BDT_CrossEntro  = 0;
       bool Use_BDT_SdivStSpB   = 0;
       bool Use_BDT_MisClass    = 0;
       bool Use_BDT_Bagging_Gini= 0;
       // ---------------------------------------------------------------
       bool EvaluateVariables  = 0; // perform evaluation for each input variable
       // ---------------------------------------------------------------

       TFile* target = TFile::Open( "./results/TMVA.root", "RECREATE" );
       std::cout<<"making factory.."<<std::endl;
       TMVA::Factory *factory = new TMVA::Factory( "TMVAnalysis", target, "" ) ;
       std::cout<<"factory made.."<<std::endl;
       vector<TString>* inputVars = new vector<TString>;
       inputVars->push_back("LepAEn");
       inputVars->push_back("LepBEn");
       inputVars->push_back("dilmass");
       inputVars->push_back("metmag");
       //inputVars->push_back("ChargeProduct");
       inputVars->push_back("addEt");
       inputVars->push_back("dPhiMetLJ");
       inputVars->push_back("ntightjets");
       inputVars->push_back("jet1_Et");
       inputVars->push_back("jet2_Et");
       inputVars->push_back("dildphi");
       factory->SetInputTrees(GetSubAna("WW")->GetDescVar(), GetSubAna("higgsww")->GetDescVar());
       std::cout<<"Input Trees set"<<std::endl;
       factory->SetInputVariables( inputVars );
       TCut mycut = "";
       std::cout<<"factory PrepareTrainingAndTestTree..."<<std::endl;
       factory->PrepareTrainingAndTestTree( mycut, 2000, 4000 );
       std::cout<<"factory PrepareTrainingAndTestTree done. Booking types"<<std::endl;
       if (Use_Cuts)factory->BookMethod( "MethodCuts",  "V:GA:EffSel:30:3:10:5:0.95" );
       if (Use_Likelihood)factory->BookMethod( TMVA::Types::Likelihood, "Spline2:3"           );
       if (Use_LikelihoodD)factory->BookMethod( TMVA::Types::Likelihood, "Spline2:10:25:D");
       if (Use_Fisher)factory->BookMethod( TMVA::Types::Fisher, "Fisher" ); // Fisher method ("Fi" or "Ma")
       if (Use_CFMlpANN)factory->BookMethod( TMVA::Types::CFMlpANN, "5000:N:N"  ); // n_cycles:#nodes:#nodes:...  
       if (Use_TMlpANN)factory->BookMethod( TMVA::Types::TMlpANN, "200:N+1:N"  ); // n_cycles:#nodes:#nodes:...
       if (Use_HMatrix)factory->BookMethod( TMVA::Types::HMatrix ); // H-Matrix (chi2-squared) method
       if (Use_PDERS)factory->BookMethod( TMVA::Types::PDERS, "Adaptive:50:100:50:0.99" );
       if (Use_BDT_GiniIndex)factory->BookMethod( TMVA::Types::BDT, "200:AdaBoost:GiniIndex:10:0.:20" );
       if (Use_BDT_CrossEntro)factory->BookMethod( TMVA::Types::BDT, "200:AdaBoost:CrossEntropy:10:0.:20" );
       if (Use_BDT_SdivStSpB)factory->BookMethod( TMVA::Types::BDT, "200:AdaBoost:SdivSqrtSplusB:10:0.:20" );
       if (Use_BDT_MisClass)factory->BookMethod( TMVA::Types::BDT, "200:AdaBoost:MisClassificationError:10:0.:20" );
       if (Use_BDT_Bagging_Gini)factory->BookMethod( TMVA::Types::BDT, "200:Bagging:GiniIndex:10:0.:20","bagging" );
       std::cout<<"methods set; training..."<<std::endl;
       factory->TrainAllMethods();
       std::cout<<"testing..."<<std::endl;
       factory->TestAllMethods();
       std::cout<<"testing done."<<std::endl;
       if (EvaluateVariables) factory->EvaluateAllVariables();
       factory->EvaluateAllMethods();
       target->Close();
       delete factory;
       delete inputVars;*/
     }

     void JetNet()
     {
       Root_Jetnet* net = new Root_Jetnet;
       net->setNode(11, 22, 1);
       net->setFileNum(2);     
       net->setFile(0, "./results/rjn-higgsww-ALL-2.txt");
       net->setFileDescription(0,"hww");
       net->setFile(1, "./results/rjn-WW-ALL-2.txt");
       net->setFileDescription(1,"WWbkg");
       net->setTarget(0,0,1.0);
       net->setTarget(1,0,0.0);
       net->setTrain(0, 3400, 2000); 
       net->setTrain(1, 3400, 2000); 
       net->setExeDir("/home/tdavies/rootjetnet/Root_Jetnet/src/");
       net->setOutDir("./results/");
       net->setNet(0, 0, 1000, 0); // (min. type, init weights, epoch #, print spew) 
       net->setUpdateFreq(0);
       net->setLearnRate(0.0);
       net->setMomentum(999.9);
       net->setInverseTemp(0.0);
       net->setPruning(1);
       net->setNetFile("./netspecs.dat");
       net->getNet(); 
       net->jetnet("11111111111");
       net->plotInput1D();
       net->plotPerform();
       net->plotNN();
     }

     void likelihood(bool);
     void TomLikelihood7708();
     void TomLikelihoodDetailed(int nbin);
     void TomLikelihoodTop(const char * thehist="TomdilDphi");
     void mirclikelihood(bool sim);
     void TomLikehisdata();
     void TomLikelihoodemu();
     void WWXS();
     //void MassPlot();

     std::ostream& PrintAcceptances(std::ostream& os=std::cout);
     std::vector<std::string> getDatasets();
     std::vector<std::string> getAcceptanceVector();


     void DebugOut(char * statement, bool condition=true)
     {
       if(_debug&&condition)std::cout<<statement<<std::endl;
     }

     void XSroutine();
     void XSroutine2();
     void XSroutine3();
     void XSroutine4();
     double GetHwwXS(int HMass, bool Stn);
     int GetHwwEvents(int HMass, bool Stn);
     
     //void XSroutine(char * dataname, char * mc_name, double dataintegral, TDRunInfoSet mc_runinfo);
     TDVnEEE crossection(const char * channel, TDVnEEE dataint, TDVnEEE mcint, TDVnEEE backgrnd,  TDVnEEE Lum, TDVnEEE mcEventsGenerated, bool print=false, ostream& os=std::cout);
    
     void ExcelXS(const char * channel, TDVnEEE dataint, TDVnEEE mcint, TDVnEEE backgrnd,  TDVnEEE Lum, TDVnEEE mcEventsGenerated, double meanweight, ostream& os=std::cout);
     void ExcelXSErr(const char * channel, TDVnEEE dataint, TDVnEEE mcint, TDVnEEE backgrnd,  TDVnEEE Lum, TDVnEEE mcEventsGenerated, double meanweight, ostream& os);
     void ExcelDataOut(TDVnEEE xs, ostream& os=std::cout);
     
    std::vector<TDSubAna*>::iterator GetSubAnaIt () {return SubAnaIt;}
    //void InitiateIt  (){SubAnaIt=SubAnas.begin();}
    bool InitiateIt  ()
    {
      SubAnaIt=SubAnas.begin();
      if(SubAnaIt==SubAnas.end())return true;
      return false;
    }
    bool IncrementIt ()
    {
      if(SubAnaIt==SubAnas.end())return true;
      ++SubAnaIt;
      if(SubAnaIt==SubAnas.end())
      {
        InitiateIt();
        return true;
      }
      return false;
    }

    TDSubAna* GetSubAna(char * subananame)
    {
      InitiateIt();do{
        if(std::string((*SubAnaIt)->GetDataSource())==std::string(subananame))
        return(*SubAnaIt);
      }while(!IncrementIt());
      return NULL;
    }

    TDHistSetSet * GetHistSets(){return HistSets;}
    std::vector<std::string> GetLikelihoodTypes(){std::cout<<"GetLikelihoodTypes; size:"<<_likelihoodtypes.size()<<std::endl; return _likelihoodtypes;}
    std::vector<std::string> GetLikelihoodValues(){std::cout<<"GetLikelihoodValues; size:"<<_likelihoodvalues.size()<<std::endl; return _likelihoodvalues;}
    
   protected:
     std::vector<TDSubAna*> SubAnas;
     std::vector<TDSubAna*>::iterator SubAnaIt;
     //TDCanvMan CanvMan1;
     //TDHistSet * Hists1;
     //TDHistSet * Hists2;
     bool _debug;
     TDHistSetSet * HistSets;
     std::string _fileLocation;
     ofstream fitinfo;//("fitinfo.txt");
     TDDilSet _dilset;
     TDCanvMan * CanvMan1;
     std::ofstream& _anaFile;
     DilWeights * _dilWeight;
     int _weighting;
     std::ofstream& _excelFile;
     std::ofstream& _excelDataFile;
     TDhwwlvlvGoodRun13 _hwwlvlvGoodRun13;
     TDhwwlvlvGoodRun _hwwlvlvGoodRun;	  
     DilType _dilType;
     TDRunSet _runstoignore;
     TDRunSet GoodRuns;
     TDRunSet MCGoodRuns;
     std::ofstream& _excelMassFile;
     std::ofstream _cutTableFile;
     std::ofstream _cutTableTexFile;
     float _higgsMass;
//     TDHistSetSet * _SecondaryHistSets;
     int _minRun;
     int _maxRun;
     bool _stn;
     bool _mircea;
     std::vector<std::string> _likelihoodtypes;
     std::vector<std::string> _likelihoodvalues;
};
#endif