}
for (Int_t iSideDigit=0; iSideDigit<fSideDigits->GetEntries(); iSideDigit++) {
- AliMFTDigit *newDigit = (AliMFTDigit*) fSideDigits->At(iSideDigit);
- new ((*pSDigList[sDigit.GetPlane()])[pSDigList[sDigit.GetPlane()]->GetEntries()]) AliMFTDigit(*newDigit);
+ AliMFTDigit *newSDigit = (AliMFTDigit*) fSideDigits->At(iSideDigit);
+ new ((*pSDigList[sDigit.GetPlane()])[pSDigList[sDigit.GetPlane()]->GetEntries()]) AliMFTDigit(*newSDigit);
}
- fSideDigits->Clear();
+ fSideDigits->Delete();
}
+
+ // ------------ In case we should simulate a rectangular pattern of pixel...
+
+ for (Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
+ if (fSegmentation->GetPlane(iPlane)->HasPixelRectangularPatternAlongY()) {
+ Int_t nSDigits = pSDigList[iPlane]->GetEntries();
+ for (Int_t iSDigit=0; iSDigit<nSDigits; iSDigit++) {
+ AliMFTDigit *mySDig = (AliMFTDigit*) (pSDigList[iPlane]->At(iSDigit));
+ if (mySDig->GetPixelX()%2 == mySDig->GetPixelY()%2) { // both pair or both odd
+ Int_t xPixelNew = mySDig->GetPixelX();
+ Int_t yPixelNew = mySDig->GetPixelY()+1;
+ if (fSegmentation->DoesPixelExist(mySDig->GetDetElemID(), xPixelNew, yPixelNew)) {
+ AliMFTDigit newSDigit;
+ newSDigit.SetEloss(0.);
+ newSDigit.SetDetElemID(mySDig->GetDetElemID());
+ newSDigit.SetPlane(iPlane);
+ newSDigit.SetPixID(xPixelNew, yPixelNew, 0);
+ newSDigit.SetPixWidth(fSegmentation->GetPixelSizeX(newSDigit.GetDetElemID()),
+ fSegmentation->GetPixelSizeY(newSDigit.GetDetElemID()),
+ fSegmentation->GetPixelSizeZ(newSDigit.GetDetElemID()));
+ newSDigit.SetPixCenter(fSegmentation->GetPixelCenterX(newSDigit.GetDetElemID(), xPixelNew),
+ fSegmentation->GetPixelCenterY(newSDigit.GetDetElemID(), yPixelNew),
+ fSegmentation->GetPixelCenterZ(newSDigit.GetDetElemID(), 0));
+ new ((*pSDigList[iPlane])[pSDigList[iPlane]->GetEntries()]) AliMFTDigit(newSDigit);
+ }
+ }
+ else { // pair-odd
+ Int_t xPixelNew = mySDig->GetPixelX();
+ Int_t yPixelNew = mySDig->GetPixelY()-1;
+ if (fSegmentation->DoesPixelExist(mySDig->GetDetElemID(), xPixelNew, yPixelNew)) {
+ AliMFTDigit newSDigit;
+ newSDigit.SetEloss(0.);
+ newSDigit.SetDetElemID(mySDig->GetDetElemID());
+ newSDigit.SetPlane(iPlane);
+ newSDigit.SetPixID(xPixelNew, yPixelNew, 0);
+ newSDigit.SetPixWidth(fSegmentation->GetPixelSizeX(newSDigit.GetDetElemID()),
+ fSegmentation->GetPixelSizeY(newSDigit.GetDetElemID()),
+ fSegmentation->GetPixelSizeZ(newSDigit.GetDetElemID()));
+ newSDigit.SetPixCenter(fSegmentation->GetPixelCenterX(newSDigit.GetDetElemID(), xPixelNew),
+ fSegmentation->GetPixelCenterY(newSDigit.GetDetElemID(), yPixelNew),
+ fSegmentation->GetPixelCenterZ(newSDigit.GetDetElemID(), 0));
+ new ((*pSDigList[iPlane])[pSDigList[iPlane]->GetEntries()]) AliMFTDigit(newSDigit);
+ }
+ }
+ }
+ }
+ }
+
+ //------------------------------------------------------------------------
AliDebug(1,"Stop Hits2SDigitsLocal");
// default constructor
for (Int_t iPlane=0; iPlane<fNMaxPlanes; iPlane++) {
- fHistNClustersPerEvent[iPlane] = 0;
- fHistNPixelsPerCluster[iPlane] = 0;
- fHistClusterSizeX[iPlane] = 0;
- fHistClusterSizeY[iPlane] = 0;
- fClusterScatterPlotXY[iPlane] = 0;
+ fHistNClustersPerEvent[iPlane] = 0;
+ fHistNPixelsPerCluster[iPlane] = 0;
+ fHistClusterSizeX[iPlane] = 0;
+ fHistClusterSizeY[iPlane] = 0;
+ fHistClusterRadialPosition[iPlane] = 0;
+ fClusterScatterPlotXY[iPlane] = 0;
}
}
fHistNPixelsPerCluster[iPlane] -> Fill(cluster->GetSize());
fHistClusterSizeX[iPlane] -> Fill(cluster->GetErrX()*1.e4); // converted in microns
fHistClusterSizeY[iPlane] -> Fill(cluster->GetErrY()*1.e4); // converted in microns
- fClusterScatterPlotXY[iPlane] -> Fill(cluster->GetX(), cluster->GetY());
+ fHistClusterRadialPosition[iPlane] -> Fill(TMath::Sqrt(cluster->GetX()*cluster->GetX()+cluster->GetY()*cluster->GetY()));
+ fClusterScatterPlotXY[iPlane] -> Fill(cluster->GetX(), cluster->GetY());
}
}
//------------------------------------------------------------
- Int_t rMax = Int_t(10.*(fMFT->GetSegmentation()->GetPlane(iPlane)->GetRMaxSupport()));
+ Int_t rMax = Int_t(10.*(fMFT->GetSegmentation()->GetPlane(iPlane)->GetRMaxSupport())); // expressed in mm
+
+ fHistClusterRadialPosition[iPlane] = new TH1D(Form("fHistClusterRadialPosition_Pl%02d",iPlane),
+ Form("Cluster radial position (Plane%02d)",iPlane),
+ rMax, 0, Double_t(rMax)/10.);
+
fClusterScatterPlotXY[iPlane] = new TH2D(Form("fClusterScatterPlotXY_Pl%02d",iPlane),
Form("Cluster scatter plot (Plane%02d)",iPlane),
2*rMax+1, (-rMax-0.5)/10., (rMax+0.5)/10., 2*rMax+1, (-rMax-0.5)/10., (rMax+0.5)/10.);
- fClusterScatterPlotXY[iPlane] -> Sumw2();
+ fHistClusterRadialPosition[iPlane] -> SetXTitle("R [cm]");
+ fClusterScatterPlotXY[iPlane] -> SetXTitle("X [cm]");
+ fClusterScatterPlotXY[iPlane] -> SetYTitle("Y [cm]");
+
+ fHistClusterRadialPosition[iPlane] -> Sumw2();
+ fClusterScatterPlotXY[iPlane] -> Sumw2();
}
AliInfo("Writing QA histos...");
+ // ---- equalize radial clusters distribution ----------------------
+
+ for (Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
+ for (Int_t iBin=1; iBin<=fHistClusterRadialPosition[iPlane]->GetNbinsX(); iBin++) {
+ Double_t rMin = fHistClusterRadialPosition[iPlane]->GetBinLowEdge(iBin); // in cm
+ Double_t rMax = fHistClusterRadialPosition[iPlane]->GetBinWidth(iBin) + rMin; // in cm
+ Double_t area = 100.*TMath::Pi()*(rMax*rMax - rMin*rMin); // in mm^2
+ Double_t density = fHistClusterRadialPosition[iPlane]->GetBinContent(iBin)/area;
+ fHistClusterRadialPosition[iPlane]->SetBinContent(iBin, density);
+ fHistClusterRadialPosition[iPlane]->SetBinError(iBin, fHistClusterRadialPosition[iPlane]->GetBinError(iBin)/area);
+ }
+ fHistClusterRadialPosition[iPlane] -> SetBinContent(1, fEv); // "scaler" bin
+ }
+
+ // -----------------------------------------------------------------
+
fFileOut->cd();
for (Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
- fHistNClustersPerEvent[iPlane] -> Write();
- fHistNPixelsPerCluster[iPlane] -> Write();
- fHistClusterSizeX[iPlane] -> Write();
- fHistClusterSizeY[iPlane] -> Write();
- fClusterScatterPlotXY[iPlane] -> Write();
+ fHistNClustersPerEvent[iPlane] -> Write();
+ fHistNPixelsPerCluster[iPlane] -> Write();
+ fHistClusterSizeX[iPlane] -> Write();
+ fHistClusterSizeY[iPlane] -> Write();
+ fHistClusterRadialPosition[iPlane] -> Write();
+ fClusterScatterPlotXY[iPlane] -> Write();
}
fFileOut -> Close();
static const Int_t fNMaxPlanes = AliMFTConstants::fNMaxPlanes;
TH1D *fHistNClustersPerEvent[fNMaxPlanes], *fHistNPixelsPerCluster[fNMaxPlanes];
- TH1D *fHistClusterSizeX[fNMaxPlanes], *fHistClusterSizeY[fNMaxPlanes];
+ TH1D *fHistClusterSizeX[fNMaxPlanes], *fHistClusterSizeY[fNMaxPlanes], *fHistClusterRadialPosition[fNMaxPlanes];
TH2D *fClusterScatterPlotXY[fNMaxPlanes];
AliLoader *fMFTLoader;
};
-//======================================================================================================
+//====================================================================================================================================================
#endif
static const Int_t fNMaxPlanes = 20;
- static const Int_t fNMaxDigitsPerCluster = 12; ///< max number of digits per cluster
+ static const Int_t fNMaxDigitsPerCluster = 50; ///< max number of digits per cluster
static const Double_t fCutForAvailableDigits; ///<
static const Double_t fCutForAttachingDigits; ///<
Double_t GetPixelWidthX() const { return fPixelWidthX; }
Double_t GetPixelWidthY() const { return fPixelWidthY; }
Double_t GetPixelWidthZ() const { return fPixelWidthZ; }
-
+
protected:
static const Double_t fElossPerElectron;
fEquivalentSiliconBeforeBack(0),
fActiveElements(0),
fReadoutElements(0),
- fSupportElements(0)
+ fSupportElements(0),
+ fHasPixelRectangularPatternAlongY(kFALSE)
{
// default constructor
fEquivalentSiliconBeforeBack(0),
fActiveElements(new TClonesArray("THnSparseC")),
fReadoutElements(new TClonesArray("THnSparseC")),
- fSupportElements(new TClonesArray("THnSparseC"))
+ fSupportElements(new TClonesArray("THnSparseC")),
+ fHasPixelRectangularPatternAlongY(kFALSE)
{
// constructor
fEquivalentSiliconBeforeBack(plane.fEquivalentSiliconBeforeBack),
fActiveElements(new TClonesArray("THnSparseC")),
fReadoutElements(new TClonesArray("THnSparseC")),
- fSupportElements(new TClonesArray("THnSparseC"))
+ fSupportElements(new TClonesArray("THnSparseC")),
+ fHasPixelRectangularPatternAlongY(plane.fHasPixelRectangularPatternAlongY)
{
// copy constructor
// base class assignement
TNamed::operator=(plane);
- fPlaneNumber = plane.fPlaneNumber;
- fZCenter = plane.fZCenter;
- fRMinSupport = plane.fRMinSupport;
- fRMax = plane.fRMax;
- fRMaxSupport = plane.fRMaxSupport;
- fPixelSizeX = plane.fPixelSizeX;
- fPixelSizeY = plane.fPixelSizeY;
- fThicknessActive = plane.fThicknessActive;
- fThicknessSupport = plane.fThicknessSupport;
- fThicknessReadout = plane.fThicknessReadout;
- fZCenterActiveFront = plane.fZCenterActiveFront;
- fZCenterActiveBack = plane.fZCenterActiveBack;
- fEquivalentSilicon = plane.fEquivalentSilicon;
- fEquivalentSiliconBeforeFront = plane.fEquivalentSiliconBeforeFront;
- fEquivalentSiliconBeforeBack = plane.fEquivalentSiliconBeforeBack;
- *fActiveElements = *plane.fActiveElements;
- *fReadoutElements = *plane.fReadoutElements;
- *fSupportElements = *plane.fSupportElements;
+ fPlaneNumber = plane.fPlaneNumber;
+ fZCenter = plane.fZCenter;
+ fRMinSupport = plane.fRMinSupport;
+ fRMax = plane.fRMax;
+ fRMaxSupport = plane.fRMaxSupport;
+ fPixelSizeX = plane.fPixelSizeX;
+ fPixelSizeY = plane.fPixelSizeY;
+ fThicknessActive = plane.fThicknessActive;
+ fThicknessSupport = plane.fThicknessSupport;
+ fThicknessReadout = plane.fThicknessReadout;
+ fZCenterActiveFront = plane.fZCenterActiveFront;
+ fZCenterActiveBack = plane.fZCenterActiveBack;
+ fEquivalentSilicon = plane.fEquivalentSilicon;
+ fEquivalentSiliconBeforeFront = plane.fEquivalentSiliconBeforeFront;
+ fEquivalentSiliconBeforeBack = plane.fEquivalentSiliconBeforeBack;
+ *fActiveElements = *plane.fActiveElements;
+ *fReadoutElements = *plane.fReadoutElements;
+ *fSupportElements = *plane.fSupportElements;
+ fHasPixelRectangularPatternAlongY = plane.fHasPixelRectangularPatternAlongY;
}
return *this;
Double_t pixelSizeY,
Double_t thicknessActive,
Double_t thicknessSupport,
- Double_t thicknessReadout) {
+ Double_t thicknessReadout,
+ Bool_t hasPixelRectangularPatternAlongY) {
AliDebug(1, Form("Initializing Plane Structure for Plane %s", GetName()));
fThicknessSupport = thicknessSupport;
fThicknessReadout = thicknessReadout;
+ fHasPixelRectangularPatternAlongY = hasPixelRectangularPatternAlongY;
+
fZCenterActiveFront = fZCenter - 0.5*fThicknessSupport - 0.5*fThicknessActive;
fZCenterActiveBack = fZCenter + 0.5*fThicknessSupport + 0.5*fThicknessActive;
Double_t pixelSizeY,
Double_t thicknessActive,
Double_t thicknessSupport,
- Double_t thicknessReadout);
+ Double_t thicknessReadout,
+ Bool_t hasPixelRectangularPatternAlongY);
Bool_t CreateStructure();
Double_t GetEquivalentSiliconBeforeBack() const { return fEquivalentSiliconBeforeBack; }
Int_t GetNumberOfChips(Option_t *opt);
+ Bool_t HasPixelRectangularPatternAlongY() { return fHasPixelRectangularPatternAlongY; }
private:
TClonesArray *fActiveElements, *fReadoutElements, *fSupportElements;
+ Bool_t fHasPixelRectangularPatternAlongY;
+
ClassDef(AliMFTPlane, 1)
};
// constructor
Float_t zCenter, rMin, rMax, pixelSizeX, pixelSizeY, thicknessActive, thicknessSupport, thicknessReadout;
- Float_t equivalentSilicon, equivalentSiliconBeforeFront, equivalentSiliconBeforeBack;
+ Float_t equivalentSilicon, equivalentSiliconBeforeFront, equivalentSiliconBeforeBack, hasPixelRectangularPatternAlongY;
TFile *geomFile = new TFile(nameGeomFile);
TNtuple *geomNtuple = (TNtuple*) geomFile->Get("AliMFTGeometry");
geomNtuple -> SetBranchAddress("equivalentSilicon", &equivalentSilicon);
geomNtuple -> SetBranchAddress("equivalentSiliconBeforeFront", &equivalentSiliconBeforeFront);
geomNtuple -> SetBranchAddress("equivalentSiliconBeforeBack", &equivalentSiliconBeforeBack);
+ if (geomNtuple -> GetBranch("hasPixelRectangularPatternAlongY")) {
+ geomNtuple -> SetBranchAddress("hasPixelRectangularPatternAlongY", &hasPixelRectangularPatternAlongY);
+ }
+ else hasPixelRectangularPatternAlongY = 0.;
Int_t nPlanes = geomNtuple->GetEntries();
zCenter = TMath::Abs(zCenter);
AliMFTPlane *plane = new AliMFTPlane(Form("MFTPlane_%02d", iPlane), Form("MFTPlane_%02d", iPlane));
- plane -> Init(iPlane, zCenter, rMin, rMax, pixelSizeX, pixelSizeY, thicknessActive, thicknessSupport, thicknessReadout);
+
+ plane -> Init(iPlane,
+ zCenter,
+ rMin,
+ rMax,
+ pixelSizeX,
+ pixelSizeY,
+ thicknessActive,
+ thicknessSupport,
+ thicknessReadout,
+ (hasPixelRectangularPatternAlongY>0.5));
+
plane -> SetEquivalentSilicon(equivalentSilicon);
plane -> SetEquivalentSiliconBeforeFront(equivalentSiliconBeforeFront);
plane -> SetEquivalentSiliconBeforeBack(equivalentSiliconBeforeBack);
Bool_t AliMFTSegmentation::Hit2PixelID(Double_t xHit, Double_t yHit, Int_t detElemID, Int_t &xPixel, Int_t &yPixel) {
+ // xPixel and yPixel start from 0
+
THnSparseC *detElem = GetDetElem(detElemID);
if ( xHit<detElem->GetAxis(0)->GetXmin() ||
//====================================================================================================================================================
+Bool_t AliMFTSegmentation::DoesPixelExist(Int_t detElemID, Int_t xPixel, Int_t yPixel) {
+
+ THnSparseC *detElem = GetDetElem(detElemID);
+
+ if (xPixel>=0 && xPixel<detElem->GetAxis(0)->GetNbins() && yPixel>=0 && yPixel<detElem->GetAxis(1)->GetNbins()) return kTRUE;
+ else return kFALSE;
+
+}
+
+//====================================================================================================================================================
+
Int_t GetNPlanes() const { return fMFTPlanes->GetEntries(); }
AliMFTPlane* GetPlane(Int_t iPlane) const { if (iPlane>=0 && iPlane<fMFTPlanes->GetEntries()) return (AliMFTPlane*) fMFTPlanes->At(iPlane); else return NULL; }
+
+ Bool_t DoesPixelExist(Int_t detElemID, Int_t xPixel, Int_t yPixel);
protected:
Int_t nMassBin = 100,
Double_t massMin = 0.,
Double_t massMax = 10.,
+ Bool_t useCutOnOffsetChi2 = kFALSE,
+ Int_t maxNWrongClusters = 999,
const Char_t *outDir = ".",
Bool_t singleMuonAnalysis = kTRUE,
Bool_t muonPairAnalysis = kTRUE,
Int_t firstEvent = -1,
Int_t lastEvent = -1,
Int_t myRandom = 0,
- Int_t maxNWrongClusters = 999,
Double_t ptMinSingleMuons = 0.0) {
gROOT -> LoadMacro("./AliMuonForwardTrackAnalysis.cxx+");
myAnalysis->SetInputDir(readDir);
myAnalysis->SetOutputDir(outDir);
myAnalysis->SetMassRange(nMassBin, massMin, massMax);
- myAnalysis->SetPtDimuRange(10, 0., 5.);
myAnalysis->SetSingleMuonAnalysis(singleMuonAnalysis);
myAnalysis->SetMuonPairAnalysis(muonPairAnalysis);
myAnalysis->SetOption(option);
- myAnalysis->SetMatchTrigger(kTRUE);
myAnalysis->SetMaxNWrongClustersMC(maxNWrongClusters);
myAnalysis->SetPtMinSingleMuons(ptMinSingleMuons);
+ myAnalysis->UseCutOnOffsetChi2(useCutOnOffsetChi2);
- myAnalysis->UseCutOnOffsetChi2(kFALSE); // cut on the single muons
-
+ myAnalysis->SetPtDimuRange(10, 0., 5.);
+ myAnalysis->SetMatchTrigger(kTRUE);
myAnalysis->UseBransonForCut(kFALSE);
myAnalysis->UseBransonForKinematics(kFALSE);
AliMuonForwardTrackPair::AliMuonForwardTrackPair():
TObject(),
fMuonForwardTracks(0),
- fKinemMC(0,0,0,0)
+ fKinemMC(0,0,0,0),
+ fKinem(0,0,0,0),
+ fIsKinemSet(kFALSE)
{
// default constructor
AliMuonForwardTrackPair::AliMuonForwardTrackPair(AliMuonForwardTrack *track0, AliMuonForwardTrack *track1):
TObject(),
fMuonForwardTracks(0),
- fKinemMC(0,0,0,0)
+ fKinemMC(0,0,0,0),
+ fKinem(0,0,0,0),
+ fIsKinemSet(kFALSE)
{
fMuonForwardTracks = new TClonesArray("AliMuonForwardTrack", 2);
AliMuonForwardTrackPair::AliMuonForwardTrackPair(const AliMuonForwardTrackPair& trackPair):
TObject(trackPair),
fMuonForwardTracks(trackPair.fMuonForwardTracks),
- fKinemMC(trackPair.fKinemMC)
+ fKinemMC(trackPair.fKinemMC),
+ fKinem(trackPair.fKinem),
+ fIsKinemSet(trackPair.fIsKinemSet)
{
// copy constructor
Clear();
fMuonForwardTracks = trackPair.fMuonForwardTracks;
+ fKinemMC = trackPair.fKinemMC;
+ fKinem = trackPair.fKinem;
+ fIsKinemSet = trackPair.fIsKinemSet;
return *this;
//====================================================================================================================================================
-Double_t AliMuonForwardTrackPair::GetMass(Double_t z, Int_t nClusters) {
-
- Int_t idCluster[2] = {0};
- if (nClusters>0) {
- idCluster[0] = ((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetNMFTClusters() - nClusters;
- idCluster[1] = ((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetNMFTClusters() - nClusters;
- }
- if (idCluster[0]<0) idCluster[0] = 0;
- if (idCluster[1]<0) idCluster[1] = 0;
-
- Double_t momentum[2] = {0};
-
- AliMUONTrackParam *param0 = ((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetTrackParamAtMFTCluster(idCluster[0]);
- AliMUONTrackParam *param1 = ((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetTrackParamAtMFTCluster(idCluster[1]);
-
- AliDebug(2, Form("MFT before extrap: 1st muon = (%f, %f, %f) 2nd muon = (%f, %f, %f)",
- param0->Px(), param0->Py(), param0->Pz(),
- param1->Px(), param1->Py(), param1->Pz()));
-
- if (TMath::Abs(z)<1e6) {
- AliDebug(2, Form("Extrapolating 1st muon from z = %f to z = %f", param0->GetZ(), z));
- AliMUONTrackExtrap::ExtrapToZCov(param0, z);
- AliDebug(2, Form("Extrapolating 2nd muon from z = %f to z = %f", param1->GetZ(), z));
- AliMUONTrackExtrap::ExtrapToZCov(param1, z);
- }
-
- AliDebug(2, Form("MFT after extrap: 1st muon = (%f, %f, %f) 2nd muon = (%f, %f, %f)",
- param0->Px(), param0->Py(), param0->Pz(),
- param1->Px(), param1->Py(), param1->Pz()));
-
- momentum[0] = (param0->P());
- momentum[1] = (param1->P());
-
- Double_t mMu = TDatabasePDG::Instance()->GetParticle("mu-")->Mass();
-
- TLorentzVector dimu;
-
- dimu.SetE(TMath::Sqrt(mMu*mMu + momentum[0]*momentum[0]) + TMath::Sqrt(mMu*mMu + momentum[1]*momentum[1]));
-
- dimu.SetPx(param0->Px() + param1->Px());
- dimu.SetPy(param0->Py() + param1->Py());
- dimu.SetPz(param0->Pz() + param1->Pz());
-
- return dimu.M();
-
-}
-
-//====================================================================================================================================================
-
Double_t AliMuonForwardTrackPair::GetMassWithoutMFT(Double_t x, Double_t y, Double_t z, Int_t nClusters) {
Int_t idCluster[2] = {0};
void AliMuonForwardTrackPair::SetKinemMC() {
+ if ( !(((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetMCTrackRef()) ||
+ !(((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetMCTrackRef()) ) return;
+
AliDebug(2, Form("MC: 1st muon = (%f, %f, %f) 2nd muon = (%f, %f, %f)",
((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetMCTrackRef()->Px(),
((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetMCTrackRef()->Py(),
//====================================================================================================================================================
+void AliMuonForwardTrackPair::SetKinem(Double_t z, Int_t nClusters) {
+
+// if (!fMuonForwardTracks) return kFALSE;
+// if (!fMuonForwardTracks->At(0) || !fMuonForwardTracks->At(1)) return kFALSE;
+
+ Int_t idCluster[2] = {0};
+ if (nClusters>0) {
+ idCluster[0] = ((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetNMFTClusters() - nClusters;
+ idCluster[1] = ((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetNMFTClusters() - nClusters;
+ }
+ if (idCluster[0]<0) idCluster[0] = 0;
+ if (idCluster[1]<0) idCluster[1] = 0;
+
+ Double_t momentum[2] = {0};
+
+ AliMUONTrackParam *param0 = ((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetTrackParamAtMFTCluster(idCluster[0]);
+ AliMUONTrackParam *param1 = ((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetTrackParamAtMFTCluster(idCluster[1]);
+
+ AliDebug(2, Form("MFT before extrap: 1st muon = (%f, %f, %f) 2nd muon = (%f, %f, %f)",
+ param0->Px(), param0->Py(), param0->Pz(),
+ param1->Px(), param1->Py(), param1->Pz()));
+
+ if (TMath::Abs(z)<1e6) {
+ AliDebug(2, Form("Extrapolating 1st muon from z = %f to z = %f", param0->GetZ(), z));
+ AliMUONTrackExtrap::ExtrapToZCov(param0, z);
+ AliDebug(2, Form("Extrapolating 2nd muon from z = %f to z = %f", param1->GetZ(), z));
+ AliMUONTrackExtrap::ExtrapToZCov(param1, z);
+ }
+
+ AliDebug(2, Form("MFT after extrap: 1st muon = (%f, %f, %f) 2nd muon = (%f, %f, %f)",
+ param0->Px(), param0->Py(), param0->Pz(),
+ param1->Px(), param1->Py(), param1->Pz()));
+
+ momentum[0] = (param0->P());
+ momentum[1] = (param1->P());
+
+ Double_t mMu = TDatabasePDG::Instance()->GetParticle("mu-")->Mass();
+
+ fKinem.SetE(TMath::Sqrt(mMu*mMu + momentum[0]*momentum[0]) + TMath::Sqrt(mMu*mMu + momentum[1]*momentum[1]));
+ fKinem.SetPx(param0->Px() + param1->Px());
+ fKinem.SetPy(param0->Py() + param1->Py());
+ fKinem.SetPz(param0->Pz() + param1->Pz());
+
+ fIsKinemSet = kTRUE;
+
+ // return fKinem.M();
+
+}
+
+//====================================================================================================================================================
+
Bool_t AliMuonForwardTrackPair::IsResonance() {
Bool_t result = kFALSE;
}
void SetKinemMC();
+ void SetKinem(Double_t z, Int_t nClusters=-1);
+ Bool_t IsKinemSet() { return fIsKinemSet; }
Double_t GetWeightedOffset(Double_t x, Double_t y, Double_t z);
- Double_t GetMass(Double_t z, Int_t nClusters=-1);
Double_t GetMassWithoutMFT(Double_t x, Double_t y, Double_t z, Int_t nClusters=-1);
Double_t GetMassMC() { return fKinemMC.M(); }
Double_t GetRapidityMC() { return fKinemMC.Rapidity(); }
Double_t GetPtMC() { return fKinemMC.Pt(); }
+ Double_t GetMass() { return fKinem.M(); }
+ Double_t GetRapidity() { return fKinem.Rapidity(); }
+ Double_t GetPt() { return fKinem.Pt(); }
Bool_t IsResonance();
protected:
TClonesArray *fMuonForwardTracks;
- TLorentzVector fKinemMC;
+ TLorentzVector fKinemMC, fKinem;
+ Bool_t fIsKinemSet;
ClassDef(AliMuonForwardTrackPair,1)
const Float_t equivalentSiliconBeforeFront[nPlanes] = { 0.e-4, 0.e-4, 0.e-4, 0.e-4, 0.e-4}; // expressed in cm
const Float_t equivalentSiliconBeforeBack[nPlanes] = { 250.e-4, 250.e-4, 250.e-4, 250.e-4, 250.e-4}; // expressed in cm
- TNtuple *geomMFT = new TNtuple("AliMFTGeometry", "ALICE MFT Geometry", "zCenter:rMin:rMax:pixelSizeX:pixelSizeY:thicknessActive:thicknessSupport:thicknessReadout:equivalentSilicon:equivalentSiliconBeforeFront:equivalentSiliconBeforeBack");
+ const Float_t hasPixelRectangularPatternAlongY[nPlanes] = {0., 0., 0., 0., 0.};
+
+ TNtuple *geomMFT = new TNtuple("AliMFTGeometry", "ALICE MFT Geometry", "zCenter:rMin:rMax:pixelSizeX:pixelSizeY:thicknessActive:thicknessSupport:thicknessReadout:equivalentSilicon:equivalentSiliconBeforeFront:equivalentSiliconBeforeBack:hasPixelRectangularPatternAlongY");
for (Int_t iPlane=0; iPlane<nPlanes; iPlane++) geomMFT -> Fill(zCenter[iPlane],
rMin[iPlane],
thicknessReadout[iPlane],
equivalentSilicon[iPlane],
equivalentSiliconBeforeFront[iPlane],
- equivalentSiliconBeforeBack[iPlane]);
+ equivalentSiliconBeforeBack[iPlane],
+ hasPixelRectangularPatternAlongY[iPlane]);
TFile *fileGeomMFT = new TFile("AliMFTGeometry.root", "recreate");
geomMFT -> Write();
// GRP from local OCDB\r
reco->SetSpecificStorage("GRP/GRP/Data",Form("local://%s",gSystem->pwd()));\r
\r
- // MUON Tracker\r
- // reco->SetSpecificStorage("MUON/Align/Data","local:///$OCDB/simulation/2008/v4-15-Release/Residual");\r
+ // MUON Tracker -> local:///$OCDB should reflect the content of alien://folder=/alice\r
+ reco->SetSpecificStorage("MUON/Align/Data", "local:///$OCDB/simulation/2008/v4-15-Release/Residual");\r
+ reco->SetSpecificStorage("MUON/Calib/RecoParam", "local:///$OCDB/simulation/2008/v4-15-Release/Full");\r
\r
reco->SetRunReconstruction("MUON MFT");\r
reco->SetRunLocalReconstruction("MUON MFT");\r
simulator->SetRunQA("ALL");\r
simulator->SetRunHLT("");\r
\r
- // MUON Tracker\r
- // simulator->SetSpecificStorage("MUON/Align/Data", "local:///$OCDB/simulation/2008/v4-15-Release/Ideal");\r
+ // MUON Tracker -> local:///$OCDB should reflect the content of alien://folder=/alice\r
+ simulator->SetSpecificStorage("MUON/Align/Data", "local:///$OCDB/simulation/2008/v4-15-Release/Ideal");\r
+ simulator->SetSpecificStorage("MUON/Calib/Gains","local:///$OCDB/simulation/2008/v4-15-Release/Ideal");\r
\r
// The rest\r
TStopwatch timer;\r