C++ TLorentzVector::Px方法代码示例

inline float kinematics::betaiSystem( TLorentzVector* pa, TLorentzVector* pb, int i)
{
	TLorentzVector pTmp = (*pa)+(*pb);
	float E = pTmp.E();
	if(E<=0.) return -99999999.;

	if     (i==1) return pTmp.Px()/E;
	else if(i==2) return pTmp.Py()/E;
	else if(i==3) return pTmp.Pz()/E;
	else _WARNING("i needs to be 1,2,3 (x,y,z), returning -99999999.");
	return -99999999.;
}

void rochcor2012jan22::momcor_mc( TLorentzVector& mu, float charge, int runopt, float& qter){
  
  //sysdev == num : deviation = num

  float ptmu = mu.Pt();
  float muphi = mu.Phi();
  float mueta = mu.Eta(); // same with mu.Eta() in Root

  float px = mu.Px();
  float py = mu.Py();
  float pz = mu.Pz();
  float e = mu.E();

  int mu_phibin = phibin(muphi);
  int mu_etabin = etabin(mueta);
  
  if(mu_phibin>=0 && mu_etabin>=0){
    
    float Mf = (mcor_bf[mu_phibin][mu_etabin] + mptsys_mc_dm[mu_phibin][mu_etabin]*mcor_bfer[mu_phibin][mu_etabin])/(mpavg[mu_phibin][mu_etabin]+mmavg[mu_phibin][mu_etabin]);
    float Af = ((mcor_ma[mu_phibin][mu_etabin]+mptsys_mc_da[mu_phibin][mu_etabin]*mcor_maer[mu_phibin][mu_etabin]) - Mf*(mpavg[mu_phibin][mu_etabin]-mmavg[mu_phibin][mu_etabin]));     
    
    float cor = 1.0/(1.0 + 2.0*Mf + charge*Af*ptmu);
    
    //for the momentum tuning - eta,phi,Q correction
    px *= cor;
    py *= cor;
    pz *= cor;
    e  *= cor;
    
    float gscler = mgscl_stat;    
    float gscl = (genm_smr/mrecm);
    
    px *= (gscl + gscler_mc_dev*gscler);
    py *= (gscl + gscler_mc_dev*gscler);
    pz *= (gscl + gscler_mc_dev*gscler);
    e  *= (gscl + gscler_mc_dev*gscler);
    
    float momscl = sqrt(px*px + py*py)/ptmu;
    
    float tune = gsf[mu_etabin]*eran.Gaus(1.0,sf[mu_etabin]);
    
    px *= (tune); 
    py *= (tune);  
    pz *= (tune);  
    e  *= (tune);   
    
    qter *= sqrt(momscl*momscl + (1.0-tune)*(1.0-tune));
  }
  
  mu.SetPxPyPzE(px,py,pz,e);
  
}

void Boost_To_Stop_Rest_Frame(TLorentzVector& stop4, TLorentzVector& chargino4, TLorentzVector& b4, TLorentzVector& neutralino4, TLorentzVector& W4, TLorentzVector& up4, TLorentzVector& down4, TLorentzVector& s4)
{
    TVector3 betaV(-stop4.Px()/stop4.Energy(),-stop4.Py()/stop4.Energy(),-stop4.Pz()/stop4.Energy());
    stop4.Boost(betaV);
    chargino4.Boost(betaV);
    b4.Boost(betaV);
    neutralino4.Boost(betaV);
    W4.Boost(betaV);
    up4.Boost(betaV);
    down4.Boost(betaV);
    s4.SetE(chargino4.P()/chargino4.M());
    s4.SetVect(chargino4.Vect().Unit()*chargino4.Gamma());
}

///------------------------------------------------------------------------
/// MEKD::computeKD - compute KD and MEs for the input processes A and B
///------------------------------------------------------------------------
int MEKD::computeKD( TString processA, TString processB,
                     TLorentzVector lept1P, int lept1Id, TLorentzVector lept2P, int lept2Id,
                     TLorentzVector lept3P, int lept3Id, TLorentzVector lept4P, int lept4Id,
                     double& kd, double& me2processA, double& me2processB )
{
	/// Prepare 4-momenta in the required format
	lept1P_i[0] = lept1P.E();
	lept1P_i[1] = lept1P.Px();
	lept1P_i[2] = lept1P.Py();
	lept1P_i[3] = lept1P.Pz();
	
	lept2P_i[0] = lept2P.E();
	lept2P_i[1] = lept2P.Px();
	lept2P_i[2] = lept2P.Py();
	lept2P_i[3] = lept2P.Pz();
	
	lept3P_i[0] = lept3P.E();
	lept3P_i[1] = lept3P.Px();
	lept3P_i[2] = lept3P.Py();
	lept3P_i[3] = lept3P.Pz();
	
	lept4P_i[0] = lept4P.E();
	lept4P_i[1] = lept4P.Px();
	lept4P_i[2] = lept4P.Py();
	lept4P_i[3] = lept4P.Pz();
	
	/// Load internal containers
	four_particle_Ps_i[0] = lept1P_i;
	four_particle_Ps_i[1] = lept2P_i;
	four_particle_Ps_i[2] = lept3P_i;
	four_particle_Ps_i[3] = lept4P_i;
	
	four_particle_IDs_i[0] = lept1Id;
	four_particle_IDs_i[1] = lept2Id;
	four_particle_IDs_i[2] = lept3Id;
	four_particle_IDs_i[3] = lept4Id;
	
	return computeKD( (string) processA.Data(), (string) processB.Data(), four_particle_Ps_i, four_particle_IDs_i, kd, me2processA, me2processB );
}

float GetWeightWjetsPolarizationF0(TLorentzVector _p4W, TLorentzVector _p4lepton,float PercentVariation, bool isWplus){

  LorentzVector p4W, p4lepton;
  p4W.SetPx(_p4W.Px());
  p4W.SetPy(_p4W.Py());
  p4W.SetPz(_p4W.Pz());
  p4W.SetE(_p4W.E());
  p4lepton.SetPx(_p4lepton.Px());
  p4lepton.SetPy(_p4lepton.Py());
  p4lepton.SetPz(_p4lepton.Pz());
  p4lepton.SetE(_p4lepton.E());

  float final_weight=1;

  float cos_theta = WjetPolarizationAngle(p4W,p4lepton);
  //  final_weight  = GetWeight( cos_theta,PercentVariation );
  //final_weight  = GetWeightFLminusFR( cos_theta,PercentVariation,p4W, isWplus );
  final_weight  = GetWeightF0( cos_theta,PercentVariation,p4W,isWplus );
  
  return final_weight;
  
  
}//end of function

///------------------------------------------------------------------------
/// MEKD::computeMEs - compute MEs for a multiple reuse
///------------------------------------------------------------------------
int MEKD::computeMEs( TLorentzVector lept1P, int lept1Id, TLorentzVector lept2P, int lept2Id,
					TLorentzVector lept3P, int lept3Id, TLorentzVector lept4P, int lept4Id )
{
	/// Prepare 4-momenta in the required format
	lept1P_i[0] = lept1P.E();
	lept1P_i[1] = lept1P.Px();
	lept1P_i[2] = lept1P.Py();
	lept1P_i[3] = lept1P.Pz();
	
	lept2P_i[0] = lept2P.E();
	lept2P_i[1] = lept2P.Px();
	lept2P_i[2] = lept2P.Py();
	lept2P_i[3] = lept2P.Pz();
	
	lept3P_i[0] = lept3P.E();
	lept3P_i[1] = lept3P.Px();
	lept3P_i[2] = lept3P.Py();
	lept3P_i[3] = lept3P.Pz();
	
	lept4P_i[0] = lept4P.E();
	lept4P_i[1] = lept4P.Px();
	lept4P_i[2] = lept4P.Py();
	lept4P_i[3] = lept4P.Pz();
	
	/// Load internal containers
	four_particle_Ps_i[0] = lept1P_i;
	four_particle_Ps_i[1] = lept2P_i;
	four_particle_Ps_i[2] = lept3P_i;
	four_particle_Ps_i[3] = lept4P_i;
	
	four_particle_IDs_i[0] = lept1Id;
	four_particle_IDs_i[1] = lept2Id;
	four_particle_IDs_i[2] = lept3Id;
	four_particle_IDs_i[3] = lept4Id;
	
	return computeMEs( four_particle_Ps_i, four_particle_IDs_i );
}

void rochcor2012::momcor_data( TLorentzVector& mu, float charge, int runopt, float& qter){
  
  float ptmu = mu.Pt();

  float muphi = mu.Phi();
  float mueta = mu.Eta(); // same with mu.Eta() in Root

  float px = mu.Px();
  float py = mu.Py();
  float pz = mu.Pz();
  float e = mu.E();
  
  int mu_phibin = phibin(muphi);
  int mu_etabin = etabin(mueta);

  float Mf = 0.0;
  float Af = 0.0;

  if(runopt==0){
    Mf = (dcor_bf[mu_phibin][mu_etabin]+mptsys_da_dm[mu_phibin][mu_etabin]*dcor_bfer[mu_phibin][mu_etabin])/(dpavg[mu_phibin][mu_etabin]+dmavg[mu_phibin][mu_etabin]);
    Af = ((dcor_ma[mu_phibin][mu_etabin]+mptsys_da_da[mu_phibin][mu_etabin]*dcor_maer[mu_phibin][mu_etabin]) - Mf*(dpavg[mu_phibin][mu_etabin]-dmavg[mu_phibin][mu_etabin]));     
  }else if(runopt==1){
    Mf = (dcor_bfD[mu_phibin][mu_etabin]+mptsys_da_dm[mu_phibin][mu_etabin]*dcor_bfDer[mu_phibin][mu_etabin])/(dpavgD[mu_phibin][mu_etabin]+dmavgD[mu_phibin][mu_etabin]);
    Af = ((dcor_maD[mu_phibin][mu_etabin]+mptsys_da_da[mu_phibin][mu_etabin]*dcor_maDer[mu_phibin][mu_etabin]) - Mf*(dpavgD[mu_phibin][mu_etabin]-dmavgD[mu_phibin][mu_etabin]));     
  }
  
  float cor = 1.0/(1.0 + 2.0*Mf + charge*Af*ptmu);
  
  px *= cor;
  py *= cor;
  pz *= cor;
  e  *= cor;
  
  //after Z pt correction
  float gscler = dgscl_stat;
  float gscl = (genm_smr/drecm);
  
  px *= (gscl + gscler_da_dev*gscler);
  py *= (gscl + gscler_da_dev*gscler);
  pz *= (gscl + gscler_da_dev*gscler);
  e  *= (gscl + gscler_da_dev*gscler);
  
  float momscl = sqrt(px*px + py*py)/ptmu;
  qter *= momscl;
  
  mu.SetPxPyPzE(px,py,pz,e);
  
}

TLorentzVector doCalEnergy(double BeamEnergy,
                           TLorentzVector Particle1,
                           TLorentzVector Particle2,
                           double nucleusMass,
                           double Particle2Mass,
                           double Particle3Mass)
{
    double E_Particle1   = Particle1.E();
    double p_Particle1_x = Particle1.Px();
    double p_Particle1_y = Particle1.Py();
    double p_Particle1_z = Particle1.Pz();
    double p_Particle1   = sqrt(TMath::Power(p_Particle1_x,2.0) +
                                    TMath::Power(p_Particle1_y,2.0) +
                                    TMath::Power(p_Particle1_z,2.0));
    double phi   = Particle2.Phi();
    double theta = Particle2.Theta();
    double b     = 2.0 * ( p_Particle1_x * cos(phi) * sin(theta) +
                           p_Particle1_y * sin(phi) * sin(theta) +
                           p_Particle1_z * cos(theta) -
                           BeamEnergy * cos(theta)
                         );
    double c     = p_Particle1 * p_Particle1 + BeamEnergy * BeamEnergy - 2.0 * BeamEnergy * p_Particle1_z;
    double d     = BeamEnergy + nucleusMass - E_Particle1;
    double e     = TMath::Power(Particle3Mass,2.0) - TMath::Power(Particle2Mass,2.0) - d * d + c;
    double Delta = 16.0 * TMath::Power(d,2.0) * (TMath::Power(e,2.0) +
                                                 TMath::Power(b * Particle2Mass,2.0) -
                                                 TMath::Power(d * Particle2Mass * 2.0,2.0));
    
    TLorentzVector NewParticle(0.0,0.0,0.0,0.0);
    if(Delta>0.)
    {
       double sol2     = (2.0 * e * b + sqrt(Delta)) / (2.0 * (4.0 * TMath::Power(d,2.0) - TMath::Power(b,2.0)));
       double newpxcal = sol2 * cos(phi) * sin(theta);
       double newpycal = sol2 * sin(phi) * sin(theta);
       double newpzcal = sol2 * cos(theta);
       double energy   = sqrt(TMath::Power(sol2,2.0) + TMath::Power(Particle2Mass,2.0));
      
       TLorentzVector NewParticle2(newpxcal,newpycal,newpzcal,energy);
       NewParticle = NewParticle2;
    }
  
    return NewParticle;
}

void rochcor::momcor_data( TLorentzVector& mu, float charge, float sysdev, int runopt){
  
  float ptmu = mu.Pt();

  float muphi = mu.Phi();
  float mueta = mu.Eta(); // same with mu.Eta() in Root

  float px = mu.Px();
  float py = mu.Py();
  float pz = mu.Pz();
  float e = mu.E();
  
  int mu_phibin = phibin(muphi);
  int mu_etabin = etabin(mueta);
  
  //float mptsys1 = sran.Gaus(0.0,sysdev);
  
  float dm = (dcor_bf[mu_phibin][mu_etabin] + mptsys_da_dm[mu_phibin][mu_etabin]*dcor_bfer[mu_phibin][mu_etabin])/dmavg[mu_phibin][mu_etabin];
  float da = dcor_ma[mu_phibin][mu_etabin] + mptsys_da_da[mu_phibin][mu_etabin]*dcor_maer[mu_phibin][mu_etabin];
  
  float cor = 1.0/(1.0 + dm + charge*da*ptmu);
  
  px *= cor;
  py *= cor;
  pz *= cor;
  e  *= cor;
  
  //after Z pt correction
  float gscler = 0.0;
  
  gscler = TMath::Sqrt( TMath::Power(dgscl_stat,2) + TMath::Power(dgscl_syst,2) );
  
  float gscl = (genm_smr/drecm);
  
  px *= (gscl + gscler_da_dev*gscler);
  py *= (gscl + gscler_da_dev*gscler);
  pz *= (gscl + gscler_da_dev*gscler);
  e  *= (gscl + gscler_da_dev*gscler);
  
  mu.SetPxPyPzE(px,py,pz,e);
  
}

//This is the function to determine the interaction point in case of a cylinder (case1), for BDXmini.
//The cylinder is with the axis along y(vertical), center at x=0,y=0,z=ldet, radius is R,height is h
double KinUtils::findInteractionPointCylinder1(const TLorentzVector &chi,double ldet,double h,double R,TVector3 &vin,TVector3 &vout,TVector3 &vhit){

	double tIN,tOUT,t2,t3,t,L;

	double px=chi.Px();
	double py=chi.Py();
	double pz=chi.Pz();

	double delta=pz*pz*ldet*ldet-(pz*pz+px*px)*(ldet*ldet-R*R);
	if (delta<0){
		cout<<"KinUtils::findInteractionPointCylinder1 error, the delta is: "<<delta<<endl;
		return 0;
	}

	//entry point
	tIN = (pz*ldet - sqrt(delta))/(px*px+pz*pz);
	t2 = (pz*ldet + sqrt(delta))/(px*px+pz*pz);

	t3 = (h/2)/py;


	if ((t3>0)&&(t3<t2)&&(t3>tIN)){
		tOUT=t3;
	}else{
		tOUT=t2;
	}

	t=Rand.Uniform(tIN,tOUT);

	vin.SetXYZ(tIN*px,tIN*py,tIN*pz);
	vout.SetXYZ(tOUT*px,tOUT*py,tOUT*pz);

	vhit.SetXYZ(t*px,t*py,t*pz);


	L = (vout - vin).Mag();



	return L;
}

void DalitzChiSq2::setap3(TLorentzVector vec){
	ap3.v.SetXYZM(vec.Px(),vec.Py(),vec.Pz(),vec.M());
	//em & ep
    	    double cos12 = (ap1.v.Px()*ap2.v.Px() + ap1.v.Py()*ap2.v.Py() + ap1.v.Pz()*ap2.v.Pz())/(ap1.v.P()*ap2.v.P());
    	   //em & gm
   	    double cos23 = (ap2.v.Px()*ap3.v.Px() + ap2.v.Py()*ap3.v.Py() + ap2.v.Pz()*ap3.v.Pz())/(ap3.v.P()*ap2.v.P());
            //ep & gm
            double cos13 = (ap1.v.Px()*ap3.v.Px() + ap1.v.Py()*ap3.v.Py() + ap1.v.Pz()*ap3.v.Pz())/(ap1.v.P()*ap3.v.P());
		
		double beta1 = ap1.v.P()/ap1.v.E();
		double beta2 = ap2.v.P()/ap2.v.E();
	    

	    double z12 = 2*((1/(beta1*beta2)) -cos12);
	    double z23 = 2*((1/beta2)-cos23);
	    double z13 = 2*((1/beta1)-cos13);
	double topterm =   M*M - 2*m_e*m_e - ( z12/(TMath::Sin(ap1.theta) * TMath::Sin(ap2.theta) *ap1.x_m*ap2.x_m));
	double bottomterm = ( (z13/(TMath::Sin(ap1.theta) * ap1.x_m)) + (z23/(TMath::Sin(ap2.theta) * ap2.x_m)) );

	ap3.x_m=topterm/bottomterm;
	
}

void print_se(TLorentzVector e_vec, TLorentzVector p_vec){
    double_t se = 2.0 * e_vec.E() * (abs(p_vec.Px()) + p_vec.E()) + deut_mass*deut_mass;
    cout << "se: " << se <<endl;
}

void print_vec(TLorentzVector v){
    cout <<"("<< v.Px() << " , " << v.Py() << " , " << v.Pz() << " , " << v.E() << ")" <<endl;
}

//.........这里部分代码省略.........
  genevent->set_event_number(_eventcount);

  // Set the PDF information
  HepMC::PdfInfo pdfinfo;
  pdfinfo.set_scalePDF(event->Q2);
  genevent->set_pdf_info(pdfinfo); 

  // We would also like to save:
  //
  // event->t;
  // event->x;
  // event->y;
  // event->s;
  // event->W;
  // event->xpom;
  // (event->polarization == transverse ? 0 : 1);
  // (event->diffractiveMode == coherent ? 0 : 1);
  // 
  // but there doesn't seem to be a good place to do so 
  // within the HepMC event information?
  //
  // t, W and Q^2 form a minial set of good variables for diffractive events
  // Maybe what I do is record the input particles to the event at the HepMC
  // vertices and reconstruct the kinematics from there? 
  
  // Create HepMC vertices and add final state particles to them

  // First, the emitter(electron)-virtual photon vertex:

  HepMC::GenVertex* egammavtx = new HepMC::GenVertex(CLHEP::HepLorentzVector(0.0,0.0,0.0,0.0));
  genevent->add_vertex(egammavtx); 

  egammavtx->add_particle_in( 
			 new HepMC::GenParticle( CLHEP::HepLorentzVector(eIn->Px(),
									 eIn->Py(),
									 eIn->Pz(),
									 eIn->E()), 
						 event->particles[0].pdgId, 
						 3 ) 
			  );

  HepMC::GenParticle *hgamma =  new HepMC::GenParticle( CLHEP::HepLorentzVector(gamma->Px(),
									 gamma->Py(),
									 gamma->Pz(),
									 gamma->E()), 
						event->particles[3].pdgId, 
						3 ); 

  egammavtx->add_particle_out(hgamma);

  egammavtx->add_particle_out( 
			 new HepMC::GenParticle( CLHEP::HepLorentzVector(eOut->Px(),
									 eOut->Py(),
									 eOut->Pz(),
									 eOut->E()), 
						 event->particles[2].pdgId, 
						 1 ) 
			  );

  // Next, the hadron-pomeron vertex:

  HepMC::GenVertex* ppomvtx = new HepMC::GenVertex(CLHEP::HepLorentzVector(0.0,0.0,0.0,0.0));
  genevent->add_vertex(ppomvtx); 


  ppomvtx->add_particle_in( 

void PlotTheta( TString inputfilename, TString outputfilename = "output.root"){
//    infile= new TFile("../PATGrid.SM.10k.root","READ");
    infile = new TFile(inputfilename, "READ");
    tree = (TTree*)infile->Get("Event");
    outputFile = new TFile(outputfilename, "RECREATE");
    outTree = new TTree("MyTree","Untersuchung der RekoObjekte");

    //TH2::SetDefaultSumw2();

    histogram__CosThetaDiff = new TH1D("histogram__CosThetaDiff", "Differenz CosTheta gen-reko", 400, -2, 2);
    histogram__CosTheta_GenReko = new TH2D("histogram__CosTheta_GenReko", "Reko-cos(theta) gegen Gen-cos(theta)", 50, -1, 1, 50, -1, 1);

    histogram__gen_A = new TH2D("histogram__gen_A", "histogram__gen_A", 5, -1, 1, 5, -1, 1);
    histogram__gen_N = new TH2D("histogram__gen_N", "histogram__gen_N", 5, -1, 1, 5, -1, 1);

    histogram__gen_LL = new TH2D("histogram__gen_LL", "histogram__gen_LL", 5, -1, 1, 5, -1, 1);
    histogram__gen_LR = new TH2D("histogram__gen_LR", "histogram__gen_LR", 5, -1, 1, 5, -1, 1);
    histogram__gen_RR = new TH2D("histogram__gen_RR", "histogram__gen_RR", 5, -1, 1, 5, -1, 1);
    histogram__gen_RL = new TH2D("histogram__gen_RL", "histogram__gen_RL", 5, -1, 1, 5, -1, 1);

    histogram__gen_Correlation = new TH2D("histogram__gen_Correlation", "histogram__gen_Correlation", 5, -1, 1, 5, -1, 1);

    histogram__A = new TH2D("histogram__A", "histogram__A", 5, -1, 1, 5, -1, 1);
    histogram__N = new TH2D("histogram__N", "histogram__N", 5, -1, 1, 5, -1, 1);

    histogram__Correlation = new TH2D("histogram__Correlation", "histogram__Correlation", 5, -1, 1, 5, -1, 1);
    histogram__Correlation_L15_B50_T1 = new TH2D("histogram__Correlation_L15_B50_T1", "histogram__Correlation_L15_B50_T1", 5, -1, 1, 5, -1, 1);
    histogram__A_L15_B50_T1 = new TH2D("histogram__A_L15_B50_T1", "histogram__A_L15_B50_T1", 5, -1, 1, 5, -1, 1);
    histogram__N_L15_B50_T1 = new TH2D("histogram__N_L15_B50_T1", "histogram__N_L15_B50_T1", 5, -1, 1, 5, -1, 1);
    histogram__Correlation_L20 = new TH2D("histogram__Correlation_L20", "histogram__Correlation_L20", 5, -1, 1, 5, -1, 1);
    histogram__A_L20 = new TH2D("histogram__A_L20", "histogram__A_L20", 5, -1, 1, 5, -1, 1);
    histogram__N_L20 = new TH2D("histogram__N_L20", "histogram__N_L20", 5, -1, 1, 5, -1, 1);
    histogram__Correlation_L20_B40 = new TH2D("histogram__Correlation_L20_B40", "histogram__Correlation_L20_B40", 5, -1, 1, 5, -1, 1);
    histogram__A_L20_B40 = new TH2D("histogram__A_L20_B40", "histogram__A_L20_B40", 5, -1, 1, 5, -1, 1);
    histogram__N_L20_B40 = new TH2D("histogram__N_L20_B40", "histogram__N_L20_B40", 5, -1, 1, 5, -1, 1);
    
    histogram__Correlation_L20_B30_T1 = new TH2D("histogram__Correlation_L20_B30_T1", "histogram__Correlation_L20_B30_T1", 5, -1, 1, 5, -1, 1);
    histogram__A_L20_B30_T1 = new TH2D("histogram__A_L20_B30_T1", "histogram__A_L20_B30_T1", 5, -1, 1, 5, -1, 1);
    histogram__N_L20_B30_T1 = new TH2D("histogram__N_L20_B30_T1", "histogram__N_L20_B30_T1", 5, -1, 1, 5, -1, 1);

    histogram__Correlation_L20_B40_T1 = new TH2D("histogram__Correlation_L20_B40_T1", "histogram__Correlation_L20_B40_T1", 5, -1, 1, 5, -1, 1);
    histogram__A_L20_B40_T1 = new TH2D("histogram__A_L20_B40_T1", "histogram__A_L20_B40_T1", 5, -1, 1, 5, -1, 1);
    histogram__N_L20_B40_T1 = new TH2D("histogram__N_L20_B40_T1", "histogram__N_L20_B40_T1", 5, -1, 1, 5, -1, 1);
    histogram__Correlation_T1 = new TH2D("histogram__Correlation_T1", "histogram__Correlation_T1", 5, -1, 1, 5, -1, 1);
    histogram__A_T1 = new TH2D("histogram__A_T1", "histogram__A_T1", 5, -1, 1, 5, -1, 1);
    histogram__N_T1 = new TH2D("histogram__N_T1", "histogram__N_T1", 5, -1, 1, 5, -1, 1);



    histogram__CosThetaDiff_TTbarPt = new TH2D("histogram__CosThetaDiff_TTbarPt", "histogram__CosThetaDiff_TTbarPt", 100, 0, 1000, 400, -2, 2);


    histogram__LeptonRelIso = new TH1D("histogram__LeptonRelIso", "histogram__LeptonRelIso", 101, 0, 1.01);


    histogram__semilepton_BLeptonMinus = new TH1D("histogram__semilepton_BLeptonMinus","histogram__semilepton_BLeptonMinus", 200, -1, 1);
    histogram__semilepton_BLeptonPlus = new TH1D("histogram__semilepton_BLeptonPlus","histogram__semilepton_BLeptonPlus", 200, -1, 1);



    histogram_nupx_gen_reco = new TH2D(" histogram_nupx_gen_reco", " histogram_nupx_gen_reco", 600, -300, 300, 600, -300, 300);
    histogram_nupy_gen_reco = new TH2D(" histogram_nupy_gen_reco", " histogram_nupy_gen_reco", 600, -300, 300, 600, -300, 300);
    histogram_nupz_gen_reco = new TH2D(" histogram_nupz_gen_reco", " histogram_nupz_gen_reco", 600, -300, 300, 600, -300, 300);

    histogram_nubpx_gen_reco = new TH2D(" histogram_nubpx_gen_reco", " histogram_nubpx_gen_reco", 600, -300, 300, 600, -300, 300);
    histogram_nubpy_gen_reco = new TH2D(" histogram_nubpy_gen_reco", " histogram_nubpy_gen_reco", 600, -300, 300, 600, -300, 300);
    histogram_nubpz_gen_reco = new TH2D(" histogram_nubpz_gen_reco", " histogram_nubpz_gen_reco", 600, -300, 300, 600, -300, 300);

    outTree->Branch("EventIsGood", &EventIsGood, "Event ist rekonstruiert/I");
    outTree->Branch("numberOfJets", &numberOfJets, "Anzahl der Jets/I");
    outTree->Branch("numberOfGoodJets", &numberOfGoodJets, "Anzahl der guten Jets/I");

    outTree->Branch("CosThetaDiff" ,&CosThetaDiff ,"Differenz im cosTheta Reko zu Gen/D");
    outTree->Branch("CosThetaPlus" ,&CosThetaPlus ,"cosTheta LeptonPlus/D");
    outTree->Branch("CosThetaMinus" ,&CosThetaMinus ,"cosTheta LeptonMinus/D");
    outTree->Branch("RekoCosThetaPlus" ,&RekoCosThetaPlus ,"cosTheta RekoLeptonPlus/D");
    outTree->Branch("RekoCosThetaMinus" ,&RekoCosThetaMinus ,"cosTheta RekoLeptonMinus/D");

    outTree->Branch("CosLeptonAngleD", &CosLeptonAngleD, "CosinusLeptonWinkel D/D");
    outTree->Branch("CosRekoLeptonAngleD", &CosRekoLeptonAngleD, "CosinusRekoLeptonWinkel D/D");

    outTree->Branch("TTbar_Pt", &TTbar_Pt, "Pt des TTbarsystems Generator/D");
    outTree->Branch("RekoTTbar_Pt", &RekoTTbar_Pt, "Pt des TTbarsystems Reko/D");
    outTree->Branch("TTbar_M", &TTbar_M, "Masse des TTbarsystems Generator/D");
    outTree->Branch("RekoTTbar_M", &RekoTTbar_M, "Masse des TTbarsystems Reko/D");
    outTree->Branch("Top_Pt", &Top_Pt, "Pt des Tops Generator/D");
    outTree->Branch("Top_M", &Top_M, "M des Tops Generator/D");
    outTree->Branch("AntiTop_Pt", &AntiTop_Pt, "Pt des AntiTops Generator/D");
    outTree->Branch("AntiTop_M", &AntiTop_M, "M des AntiTops Generator/D");
    outTree->Branch("RekoTop_Pt", &RekoTop_Pt, "Pt des Tops Reko/D");
    outTree->Branch("RekoAntiTop_Pt", &RekoAntiTop_Pt, "Pt des AntiTops Reko/D");
    outTree->Branch("RekoTop_M", &RekoTop_M, "M des Tops Reko/D");
    outTree->Branch("RekoAntiTop_M", &RekoAntiTop_M, "M des AntiTops Reko/D");
    outTree->Branch("Nu_Px", &Nu_Px, "Px des Neutrinos Generator/D");
    outTree->Branch("Nu_Py", &Nu_Py, "Py des Neutrinos Generator/D");
    outTree->Branch("Nu_Pz", &Nu_Pz, "Pz des Neutrinos Generator/D");
    outTree->Branch("AntiNu_Px", &AntiNu_Px, "Px des AntiNeutrinos Generator/D");
    outTree->Branch("AntiNu_Py", &AntiNu_Py, "Py des AntiNeutrinos Generator/D");
    outTree->Branch("AntiNu_Pz", &AntiNu_Pz, "Pz des AntiNeutrinos Generator/D");
    outTree->Branch("RekoNu_Px", &RekoNu_Px, "Px des Neutrinos Reko/D");
//.........这里部分代码省略.........