diff --git a/reco/detectors/trd/qa/CbmTrdHitProducerQa.cxx b/reco/detectors/trd/qa/CbmTrdHitProducerQa.cxx index a6f0cd1dc99dd9e553050e3ea20b41c033493792..3fd8bbf871284921afe3fa196109a2230f4bc076 100644 --- a/reco/detectors/trd/qa/CbmTrdHitProducerQa.cxx +++ b/reco/detectors/trd/qa/CbmTrdHitProducerQa.cxx @@ -58,11 +58,25 @@ InitStatus CbmTrdHitProducerQa::Init() const int nStations = fNoTrdStations * fNoTrdPerStation; for (int i = 0; i < nStations; i++) { - fvhHitPullsX.push_back(new TH1F(Form("L%iHitPullsX", i), "", 500, -50, 50)); - fvhHitPullsY.push_back(new TH1F(Form("L%iHitPullsY", i), "", 500, -50, 50)); - - histFolder->Add(fvhHitPullsX[i]); - histFolder->Add(fvhHitPullsY[i]); + fvhHitPullX.push_back(new TH1F(Form("L%iHitPullX", i), "", 500, -50, 50)); + fvhHitPullY.push_back(new TH1F(Form("L%iHitPullY", i), "", 500, -50, 50)); + fvhHitPullT.push_back(new TH1F(Form("L%iHitPullT", i), "", 500, -50, 50)); + fvhHitPullX[i]->SetCanExtend(TH1::kAllAxes); + fvhHitPullY[i]->SetCanExtend(TH1::kAllAxes); + fvhHitPullT[i]->SetCanExtend(TH1::kAllAxes); + histFolder->Add(fvhHitPullX[i]); + histFolder->Add(fvhHitPullY[i]); + histFolder->Add(fvhHitPullT[i]); + + fvhHitResX.push_back(new TH1F(Form("L%iHitResX", i), "", 500, -50, 50)); + fvhHitResY.push_back(new TH1F(Form("L%iHitResY", i), "", 500, -50, 50)); + fvhHitResT.push_back(new TH1F(Form("L%iHitResT", i), "", 500, -50, 50)); + fvhHitResX[i]->SetCanExtend(TH1::kAllAxes); + fvhHitResY[i]->SetCanExtend(TH1::kAllAxes); + fvhHitResT[i]->SetCanExtend(TH1::kAllAxes); + histFolder->Add(fvhHitResX[i]); + histFolder->Add(fvhHitResY[i]); + histFolder->Add(fvhHitResT[i]); fvhedEcut.push_back(new TH1F(Form("L%iedEcut", i), Form("dEdx of e- for layer %i, mom. cut", i), 600, 0., 60.)); fvhedEall.push_back(new TH1F(Form("L%iedEall", i), Form("dEdx of e- for layer %i", i), 600, 0., 60.)); @@ -79,8 +93,9 @@ InitStatus CbmTrdHitProducerQa::Init() fvdECanvas[i]->Divide2D(4); fOutFolder.Add(fvdECanvas[i]); - fvPullCanvas.push_back(new CbmQaCanvas(Form("cL%iPull", i), Form("Pull Distribution Layer %i", i), 2 * 400, 400)); - fvPullCanvas[i]->Divide2D(2); + fvPullCanvas.push_back( + new CbmQaCanvas(Form("cL%iPull", i), Form("Pull Distribution Layer %i", i), 3 * 400, 2 * 400)); + fvPullCanvas[i]->Divide2D(6); fOutFolder.Add(fvPullCanvas[i]); } @@ -140,7 +155,7 @@ InitStatus CbmTrdHitProducerQa::Init() void CbmTrdHitProducerQa::Exec(Option_t*) { // Loop over TRD hits - for (Int_t iHit = 0; iHit < fTrdHitCollection->GetEntriesFast(); iHit++) { + for (int iHit = 0; iHit < fTrdHitCollection->GetEntriesFast(); iHit++) { const CbmTrdHit* trdHit = (CbmTrdHit*) fTrdHitCollection->At(iHit); if (NULL == trdHit) continue; @@ -153,16 +168,15 @@ void CbmTrdHitProducerQa::Exec(Option_t*) const CbmTrdPoint* trdPoint = (CbmTrdPoint*) fTrdPointCollection->At(trdDigiMatch->GetMatchedLink().GetIndex()); if (NULL == trdPoint) continue; - const Int_t planeId = trdHit->GetPlaneId(); + const int planeId = trdHit->GetPlaneId(); if (planeId >= fNoTrdStations * fNoTrdPerStation) { cout << GetName() << ": Warning, TRD plane out of bounds, skipping hit." << endl; continue; } - const Int_t partID = (((CbmMCTrack*) fMCTrackArray->At(trdPoint->GetTrackID()))->GetPdgCode()); - const Float_t momentum = - TMath::Sqrt((trdPoint->GetPx() * trdPoint->GetPx()) + (trdPoint->GetPy() * trdPoint->GetPy()) - + (trdPoint->GetPz() * trdPoint->GetPz())); + const int partID = (((CbmMCTrack*) fMCTrackArray->At(trdPoint->GetTrackID()))->GetPdgCode()); + const float momentum = TMath::Sqrt((trdPoint->GetPx() * trdPoint->GetPx()) + (trdPoint->GetPy() * trdPoint->GetPy()) + + (trdPoint->GetPz() * trdPoint->GetPz())); //electrons if ((TMath::Abs(partID) == 11) && (momentum > fMomCutLower) && (momentum < fMomCutUpper)) { @@ -177,24 +191,29 @@ void CbmTrdHitProducerQa::Exec(Option_t*) if (TMath::Abs(partID) == 211) { fvhpidEall[planeId]->Fill((trdHit->GetELoss()) * 1000000); } // compute the Hit pulls for X and Y coordinate - const Float_t hitPosX = trdHit->GetX(); - const Float_t pointPosX = trdPoint->GetX(); - const Float_t hitPullX = (hitPosX - pointPosX); + const float hitPosX = trdHit->GetX(); + const float pointPosX = (trdPoint->GetX() + trdPoint->GetXOut()) / 2.; + const float hitResX = hitPosX - pointPosX; + const float hitPullX = hitResX / trdHit->GetDx(); + + const float hitPosY = trdHit->GetY(); + const float pointPosY = (trdPoint->GetY() + trdPoint->GetYOut()) / 2.; + const float hitResY = hitPosY - pointPosY; + const float hitPullY = hitResY / trdHit->GetDy(); - const Float_t hitPosY = trdHit->GetY(); - const Float_t pointPosY = trdPoint->GetY(); - const Float_t hitPullY = (hitPosY - pointPosY); + const double hitPosT = trdHit->GetTime(); + const double pointPosT = trdPoint->GetTime(); + const double hitResT = hitPosT - pointPosT; + const double hitPullT = hitResT / trdHit->GetTimeError(); // fill histograms - fvhHitPullsX[planeId]->Fill(hitPullX); - fvhHitPullsY[planeId]->Fill(hitPullY); - - //Float_t hitErrX = 0; - //Float_t hitErrY = 0; - // hitErrX = trdHit->GetDx(); - // hitErrX /= 10000; // micrometers to centimeters - // hitErrY = trdHit->GetDy(); - // hitErrY /= 10000; // micrometers to centimeters + fvhHitPullX[planeId]->Fill(hitPullX); + fvhHitPullY[planeId]->Fill(hitPullY); + fvhHitPullT[planeId]->Fill(hitPullT); + + fvhHitResX[planeId]->Fill(hitResX); + fvhHitResY[planeId]->Fill(hitResY); + fvhHitResT[planeId]->Fill(hitResT); } } // -------------------------------------------------------------------------- @@ -206,7 +225,7 @@ void CbmTrdHitProducerQa::Finish() { WriteHistograms(); } // ---- Write test histograms ------------------------------------------------ void CbmTrdHitProducerQa::WriteHistograms() { - for (UInt_t i = 0; i < fvdECanvas.size(); i++) { + for (int i = 0; i < fvdECanvas.size(); i++) { fvdECanvas[i]->cd(1); fvhedEcut[i]->DrawCopy("", ""); @@ -220,12 +239,24 @@ void CbmTrdHitProducerQa::WriteHistograms() fvhpidEall[i]->DrawCopy("", ""); } - for (UInt_t i = 0; i < fvPullCanvas.size(); i++) { + for (int i = 0; i < fvPullCanvas.size(); i++) { fvPullCanvas[i]->cd(1); - fvhHitPullsX[i]->DrawCopy("", ""); + fvhHitPullX[i]->DrawCopy("", ""); fvPullCanvas[i]->cd(2); - fvhHitPullsY[i]->DrawCopy("", ""); + fvhHitPullY[i]->DrawCopy("", ""); + + fvPullCanvas[i]->cd(3); + fvhHitPullT[i]->DrawCopy("", ""); + + fvPullCanvas[i]->cd(4); + fvhHitResX[i]->DrawCopy("", ""); + + fvPullCanvas[i]->cd(5); + fvhHitResY[i]->DrawCopy("", ""); + + fvPullCanvas[i]->cd(6); + fvhHitResT[i]->DrawCopy("", ""); } FairSink* sink = FairRootManager::Instance()->GetSink(); diff --git a/reco/detectors/trd/qa/CbmTrdHitProducerQa.h b/reco/detectors/trd/qa/CbmTrdHitProducerQa.h index 4709063aa9d77f053ff79dbac89da544017d0f9c..53750b3d3724eba5907fba93dc76f2b045c7253a 100644 --- a/reco/detectors/trd/qa/CbmTrdHitProducerQa.h +++ b/reco/detectors/trd/qa/CbmTrdHitProducerQa.h @@ -51,17 +51,17 @@ public: virtual void Finish(); /* Set the momentum cuts */ - void SetMomentumCuts(Float_t CutLower, Float_t CutHigher) + void SetMomentumCuts(float CutLower, float CutHigher) { fMomCutLower = CutLower; fMomCutUpper = CutHigher; } /* Set number of TRD stations */ - void SetNumberStations(Int_t nStations) { fNoTrdStations = nStations; } + void SetNumberStations(int nStations) { fNoTrdStations = nStations; } /** Set number of layers per station **/ - void SetLayersPerStations(Int_t nLayers) { fNoTrdPerStation = nLayers; } + void SetLayersPerStations(int nLayers) { fNoTrdPerStation = nLayers; } private: TFolder fOutFolder; /// output folder with histos and canvases @@ -75,17 +75,21 @@ private: TClonesArray* fMCTrackArray = nullptr; /** Number of TRD stations **/ - Int_t fNoTrdStations = 4; + int fNoTrdStations = 4; /** Number of layers per station **/ - Int_t fNoTrdPerStation = 1; + int fNoTrdPerStation = 1; /* Write test histograms */ void WriteHistograms(); /* Test histograms*/ - std::vector<TH1F*> fvhHitPullsX; - std::vector<TH1F*> fvhHitPullsY; + std::vector<TH1F*> fvhHitPullX; + std::vector<TH1F*> fvhHitPullY; + std::vector<TH1F*> fvhHitPullT; + std::vector<TH1F*> fvhHitResX; + std::vector<TH1F*> fvhHitResY; + std::vector<TH1F*> fvhHitResT; std::vector<TH1F*> fvhedEcut; std::vector<TH1F*> fvhedEall; @@ -97,8 +101,8 @@ private: std::vector<CbmQaCanvas*> fvPullCanvas; /* Momentum cuts for energy distributions */ - Float_t fMomCutLower = 1.; - Float_t fMomCutUpper = 7.; + float fMomCutLower = 1.; + float fMomCutUpper = 7.; ClassDef(CbmTrdHitProducerQa, 3) };