#include "AliMUONPairLight.h"
#include "TString.h"
ClassImp(AliMUONPairLight)
AliMUONPairLight::AliMUONPairLight() :
TObject(),
fMu0(),
fMu1(),
fCreationProcess(-999),
fIsCorrelated(kFALSE),
fCauseOfCorrelation (-1),
fIsFeedDown(kFALSE)
{
;
}
AliMUONPairLight::AliMUONPairLight(AliMUONPairLight &dimuCopy)
: TObject(dimuCopy),
fMu0(dimuCopy.fMu0),
fMu1(dimuCopy.fMu1),
fCreationProcess(dimuCopy.fCreationProcess),
fIsCorrelated(dimuCopy.fIsCorrelated),
fCauseOfCorrelation (dimuCopy.fCauseOfCorrelation),
fIsFeedDown(dimuCopy.fIsFeedDown)
{
;
}
AliMUONPairLight::~AliMUONPairLight(){
}
AliMUONPairLight& AliMUONPairLight::operator=(const AliMUONPairLight& dimuCopy)
{
if (this == &dimuCopy) return *this;
TObject::operator=(dimuCopy);
fMu0 = dimuCopy.fMu0;
fMu1 = dimuCopy.fMu1;
fCreationProcess = dimuCopy.fCreationProcess;
fIsCorrelated = dimuCopy.fIsCorrelated;
fCauseOfCorrelation = dimuCopy.fCauseOfCorrelation;
fIsFeedDown = dimuCopy.fIsFeedDown;
return *this;
}
Bool_t AliMUONPairLight::IsAResonance(){
if (!fIsCorrelated) return kFALSE;
Int_t nparents0 = fMu0.GetNParents();
Int_t nparents1 = fMu1.GetNParents();
Int_t minP = TMath::Min(nparents0, nparents1);
for (Int_t i = 0 ; i < minP; i++) {
if (fMu0.IsMotherAResonance(nparents0-1-i) && fMu1.IsMotherAResonance(nparents1-1-i) &&
fMu0.GetParentPythiaLine(nparents0-1-i)==fMu1.GetParentPythiaLine(nparents1-1-i)) {
if (nparents0-1-i) SetFeedDown(nparents0-1-i);
return kTRUE;
}
}
return kFALSE;
}
AliMUONTrackLight* AliMUONPairLight::GetMuon(Int_t index) {
if (index==0) return &fMu0;
else if (index==1) return &fMu1;
else{ printf ("Index can be either 0 or 1\n"); return 0;}
}
Int_t AliMUONPairLight::GetMuonMotherPDG(Int_t imuon, Int_t mother) {
if (imuon==0) return fMu0.GetParentPDGCode(mother);
else if (imuon==1) return fMu1.GetParentPDGCode(mother);
else { printf ("Index must be only 0 or 1\n"); return -999; }
}
void AliMUONPairLight::SetProcess(){
AliMUONTrackLight *mu1 = &fMu0;
AliMUONTrackLight *mu2 = &fMu1;
Int_t npar1 = mu1->GetNParents();
Int_t npar2 = mu2->GetNParents();
for (Int_t imoth1 = npar1-1; imoth1>=0; imoth1--) {
Int_t lineMo1 = mu1->GetParentPythiaLine(imoth1);
for (Int_t imoth2 = npar2-1; imoth2>=0; imoth2--) {
Int_t lineMo2 = mu2->GetParentPythiaLine(imoth2);
if(lineMo1 == lineMo2) {
if(mu1->IsDiquark(mu1->GetParentPDGCode(imoth1)))return;
this->SetCorrelated(kTRUE);
this->SetCauseOfCorrelation(mu1->GetParentPDGCode(imoth1));
if(!IsAResonance()) fCreationProcess = 3;
else fCreationProcess = -1;
return;
}
}
}
if(this->IsDimuonFromCorrPiK()){
this->SetCorrelated(kTRUE);
this->SetCauseOfCorrelation(mu1->GetParentPDGCode(0));
fCreationProcess = -1;
}
Int_t flavPar1 = mu1->GetParentFlavour(0);
Int_t flavPar2 = mu2->GetParentFlavour(0);
for (Int_t imoth1 = 0; imoth1 < 4; imoth1++) {
Int_t lineMo1 = mu1->GetQuarkPythiaLine(imoth1);
for (Int_t imoth2 = 0; imoth2 < 4; imoth2++) {
Int_t lineMo2 = mu2->GetQuarkPythiaLine(imoth2);
if(lineMo1 == lineMo2 && mu1->GetQuarkPDGCode(imoth1) == 21) {
if(flavPar1 == flavPar2){
this->SetCorrelated(kTRUE);
if(GetCauseOfCorrelation() == -1)
this->SetCauseOfCorrelation(mu1->GetQuarkPDGCode(imoth1));
fCreationProcess = 1;
return;
}
}
}
}
Int_t line1 = mu1->GetQuarkPythiaLine(2);
Int_t line2 = mu2->GetQuarkPythiaLine(2);
Int_t line6or7[2] = {-1, -1};
Int_t flavourLine6or7[2] = {-1, -1};
for (Int_t imoth1 = 3; imoth1>=0; imoth1--) {
Int_t lineMo1 = mu1->GetQuarkPythiaLine(imoth1);
Int_t flavour1 = TMath::Abs(mu1->GetQuarkPDGCode(imoth1));
if(lineMo1 == 6 || lineMo1 == 7){
line6or7[0] = imoth1;
flavourLine6or7[0] = flavour1;
}
for (Int_t imoth2 = 3; imoth2>=0; imoth2--) {
Int_t lineMo2 = mu2->GetQuarkPythiaLine(imoth2);
Int_t flavour2 = TMath::Abs(mu2->GetQuarkPDGCode(imoth2));
if(lineMo2 == 6 || lineMo2 == 7){
line6or7[1] = imoth2;
flavourLine6or7[1] = flavour2;
}
if((line6or7[0] > 0 && line6or7[1] > 0) &&
(flavourLine6or7[0] == 4 || flavourLine6or7[0] == 5) &&
(flavourLine6or7[1] == 4 || flavourLine6or7[1] == 5) &&
(flavPar1 == flavPar2)){
this->SetCorrelated(kTRUE);
fCreationProcess = 0;
return;
}
}
}
Int_t line2or3[2] = {-1, -1};
Int_t flavourLine2or3[2] = {-1, -1};
for (Int_t imoth1 = 3; imoth1>=0; imoth1--) {
Int_t lineMo1 = mu1->GetQuarkPythiaLine(imoth1);
Int_t flavour1 = TMath::Abs(mu1->GetQuarkPDGCode(imoth1));
if(lineMo1 == 2 || lineMo1 == 3){
line2or3[0] = imoth1;
flavourLine2or3[0] = flavour1;
}
for (Int_t imoth2 = 3; imoth2>=0; imoth2--) {
Int_t lineMo2 = mu2->GetQuarkPythiaLine(imoth2);
Int_t flavour2 = TMath::Abs(mu2->GetQuarkPDGCode(imoth2));
if(lineMo2 == 2 || lineMo2 == 3){
line2or3[1] = imoth2;
flavourLine2or3[1] = flavour2;
}
if(((line6or7[0] > 0 && (flavourLine6or7[0] == 4 || flavourLine6or7[0] == 5)) &&
(line2or3[1] > 0 && (flavourLine2or3[1] == 21 || flavourLine2or3[1] < 10))) ||
((line6or7[1] > 0 && (flavourLine6or7[1] == 4 || flavourLine6or7[1] == 5)) &&
(line2or3[0] > 0 && (flavourLine2or3[0] == 21 || flavourLine2or3[0] < 10)))){
if(flavPar1 == flavPar2){
this->SetCorrelated(kTRUE);
fCreationProcess = 2;
return;
}
}
}
}
if(line1 == line2 && (line1 == 2 || line1 == 3)){
if((TMath::Abs(mu1->GetQuarkPDGCode(1)) == 4 && TMath::Abs(mu2->GetQuarkPDGCode(1)) == 4) ||
(TMath::Abs(mu1->GetQuarkPDGCode(1)) == 5 && TMath::Abs(mu2->GetQuarkPDGCode(1)) == 5)){
if(flavPar1 == flavPar2){
this->SetCorrelated(kTRUE);
fCreationProcess = 1;
if(GetCauseOfCorrelation() == -1){
this->SetCauseOfCorrelation(mu1->GetQuarkPDGCode(1));
}
return;
}
}
}
line1 = mu1->GetQuarkPythiaLine(1);
line2 = mu2->GetQuarkPythiaLine(1);
if(line1 == line2 && (line1 == 2 || line1 == 3)){
if((TMath::Abs(mu1->GetQuarkPDGCode(0)) == 4 && TMath::Abs(mu2->GetQuarkPDGCode(0)) == 4) ||
(TMath::Abs(mu1->GetQuarkPDGCode(0)) == 5 && TMath::Abs(mu2->GetQuarkPDGCode(0)) == 5)){
if(flavPar1 == flavPar2){
this->SetCorrelated(kTRUE);
fCreationProcess = 1;
if(GetCauseOfCorrelation() == -1){
this->SetCauseOfCorrelation(mu1->GetQuarkPDGCode(1));
}
return;
}
}
}
line1 = mu1->GetQuarkPythiaLine(1);
line2 = mu2->GetQuarkPythiaLine(1);
if(line1 == line2 && (line1 == 6 || line1 == 7)){
if((TMath::Abs(mu1->GetQuarkPDGCode(0)) == 4 && TMath::Abs(mu2->GetQuarkPDGCode(0)) == 4) ||
(TMath::Abs(mu1->GetQuarkPDGCode(0)) == 5 && TMath::Abs(mu2->GetQuarkPDGCode(0)) == 5)){
if(flavPar1 == flavPar2){
this->SetCorrelated(kTRUE);
fCreationProcess = 1;
if(GetCauseOfCorrelation() == -1){
this->SetCauseOfCorrelation(mu1->GetQuarkPDGCode(1));
}
return;
}
}
}
}
void AliMUONPairLight::SetMuons(const AliMUONTrackLight& mu0, const AliMUONTrackLight& mu1){
fMu0 = mu0;
fMu1 = mu1;
this->SetProcess();
}
void AliMUONPairLight::PrintInfo(const Option_t* opt){
TString options(opt);
options.ToUpper();
if(options.Contains("H") || options.Contains("A")){
AliMUONTrackLight *mu1 = &fMu0;
AliMUONTrackLight *mu2 = &fMu1;
printf("========= History =======================\n");
printf("first muon");
mu1->PrintInfo("H");
printf("second muon");
mu2->PrintInfo("H");
printf("=========================================\n");
}
if(options.Contains("F") || options.Contains("A")){
printf("the flags set for this muon pair are:\n");
printf("=====================================\n");
if(this->IsOneTrackNotAMuon()) printf("(*) one rec. track is not a muon\n");
fIsCorrelated ? printf("(*) it is a correlated pair\n") : printf("(*) it is not a correlated pair\n");
if(IsOpenCharm()) printf("(*) correlated open charm: ");
if(IsOpenBeauty()) printf("(*) correlated open beauty: ");
if(IsOpenCharm() || IsOpenBeauty()){
switch(fCreationProcess){
case 0:
printf("pair creation");
break;
case 1:
printf("gluon splitting");
break;
case 2:
printf("flavour excitation");
break;
case 3:
printf("both muons come from same fragmented mother");
break;
}
if(this->GetMuon(0)->GetOscillation() || this->GetMuon(1)->GetOscillation())
printf("... where oscillation occured\n");
else{
if(IsOpenBeauty())
printf(" (no oscillation)\n");
else
printf("\n");
}
}
IsAResonance() ? printf("(*) it is a resonance: %d\n", this->GetMuonMotherPDG(0, fIsFeedDown)) : printf("(*) it is not a resonance\n");
fIsFeedDown ? printf("(*) mother has feed-down: %d --> %d\n", this->GetMuonMotherPDG(0,fMu0.GetNParents()-2), this->GetMuonMotherPDG(0,fMu0.GetNParents()-1)) : printf("(*) no feed-down\n");
printf("=====================================\n");
}
if(options.Contains("K") || options.Contains("A")){
Double_t *vtx = this->GetMuon(0)->GetVertex();
TLorentzVector momRec = this->GetPRec();
TLorentzVector momGen = this->GetPGen();
printf("the dimuon charge is %d\n", this->GetCharge());
printf("primary Vertex: Vx = %1.3f, Vy = %1.3f, Vz = %1.3f\n", vtx[0], vtx[1], vtx[2]);
printf("Generated: Px = %1.3f, Py = %1.3f, Pz = %1.3f\n", momGen.Px(), momGen.Py(), momGen.Pz());
printf("Reconstructed: Px = %1.3f, Py = %1.3f, Pz = %1.3f\n", momRec.Px(), momRec.Py(), momRec.Pz());
printf("Rec. variables: mass %1.3f, pT %1.3f, pseudo-rapidity %1.3f, openingAngle %1.3f (%1.3f degree), theta %1.3f (%1.3f degree), phi %1.3f (%1.3f degree)\n",
momRec.M(), momRec.Pt(), momRec.Eta(),
TMath::Pi()/180.*this->GetOpeningAngle(), this->GetOpeningAngle(),
momRec.Theta(), 180./TMath::Pi() * momRec.Theta(),
momRec.Phi(), 180./TMath::Pi() * momRec.Phi());
}
}
Double_t AliMUONPairLight::GetOpeningAngle() {
TLorentzVector pRecMu0 = fMu0.GetPRec();
TLorentzVector pRecMu1 = fMu1.GetPRec();
TVector3 pRecMu03 = pRecMu0.Vect();
TVector3 pRecMu13 = pRecMu1.Vect();
Double_t scalar = pRecMu03.Dot(pRecMu13);
Double_t modMu0 = pRecMu03.Mag();
Double_t modMu1 = pRecMu13.Mag();
Double_t theta = (TMath::ACos(scalar/(modMu0*modMu1)))*(180./TMath::Pi());
return theta;
}
Bool_t AliMUONPairLight::IsDimuonFromCorrPiK(){
AliMUONTrackLight *mu0 = this->GetMuon(0), *mu1 = this->GetMuon(1);
Bool_t fromSameLine = kFALSE;
if (mu0->IsParentPionOrKaon() &&
mu1->IsParentPionOrKaon() &&
mu1->GetQuarkPythiaLine() == mu0->GetQuarkPythiaLine()
) fromSameLine = kTRUE;
return fromSameLine;
}