AlgShowerSSList Class Reference

#include <AlgShowerSSList.h>

Inheritance diagram for AlgShowerSSList:

List of all members.

Public Member Functions
	AlgShowerSSList ()
virtual	~AlgShowerSSList ()
virtual void	RunAlg (AlgConfig &ac, CandHandle &ch, CandContext &cx)
virtual void	Trace (const char *c) const

---

Constructor & Destructor Documentation

AlgShowerSSList::AlgShowerSSList (   )

Definition at line 49 of file AlgShowerSSList.cxx. { }

AlgShowerSSList::~AlgShowerSSList (   )  [virtual]

Definition at line 54 of file AlgShowerSSList.cxx. { }

---

Member Function Documentation

void AlgShowerSSList::RunAlg (  AlgConfig &  ac, CandHandle &  ch, CandContext &  cx )  [virtual]

We begin with a nested loop over all CandSlices, and over all CandClusters within each CandSlice, generating list of all clusters in the U and V planes. Next, we execute a nested loop over U-view and V-view CandClusters, checking for matching 2D clusters. The following criteria are used: Beginning planes for the two 2D 'long' CandClusters (defined by PlaneScale) must be within DiffViewPlaneMatch., or DiffViewPlaneMatchShort if one or more clusters is short. Beginning times for the two 2D CandClusters must be within DiffViewTimeMatch. The fractional overlap of U and V clusters must exceed DiffViewPlaneCoverage, if both clusters are long. The total pulse heights of the U and V clusters must agree to within DiffViewPulseHeightCut if at least one of the clusters is long. (Note: This should probably be changed to a charge ratio) If a U/V cluster set is found which satisfies these requirements, a second loop over U clusters is performed, searching for alternative U/V pairs involving the original V cluster satisfying the requirements above, but higher in total energy. If no better alternative is found a CandShower is constructed. In this way, a given CandCluster will only be used once, in the shower with the highest total energy. Reimplemented from AlgShowerSRList. Definition at line 60 of file AlgShowerSSList.cxx. References CandHandle::AddDaughterLink(), Registry::Get(), AlgFactory::GetAlgHandle(), CandSubShowerSRHandle::GetAvgDev(), CandRecoHandle::GetBegPlane(), CandContext::GetCandRecord(), CandRecoHandle::GetCandSlice(), CandStripHandle::GetCharge(), CandSubShowerSRHandle::GetClusterID(), CandContext::GetDataIn(), CandHandle::GetDaughterIterator(), CandRecoHandle::GetEndPlane(), CandShowerHandle::GetEnergy(), CandSubShowerSRHandle::GetEnergy(), AlgFactory::GetInstance(), CandSubShowerSRHandle::GetMaxU(), CandSubShowerSRHandle::GetMaxV(), CandSubShowerSRHandle::GetMinU(), CandSubShowerSRHandle::GetMinV(), CandContext::GetMom(), CandSubShowerSRHandle::GetPlaneView(), CandSubShowerSRHandle::GetSlope(), CandStripHandle::GetTPos(), RecMinos::GetVldContext(), CandRecoHandle::GetVtxU(), CandRecoHandle::GetVtxV(), CandShowerSR::MakeCandidate(), MSG, PITCHINMETRES, CandHandle::Print(), and CandRecoHandle::SetCandSlice(). { assert(cx.GetDataIn()); if (!(cx.GetDataIn()->InheritsFrom("TObjArray"))) { return; } const CandSliceListHandle *slicelist = 0; const CandClusterListHandle *clusterlist = 0; const CandSubShowerSRListHandle *subshowerlist = 0; const TObjArray *cxin = dynamic_cast<const TObjArray *>(cx.GetDataIn()); for (Int_t i=0; i<=cxin->GetLast(); i++) { TObject *tobj = cxin->At(i); if (tobj->InheritsFrom("CandSliceListHandle")) { slicelist = dynamic_cast<CandSliceListHandle*>(tobj); } if (tobj->InheritsFrom("CandClusterListHandle")) { clusterlist = dynamic_cast<CandClusterListHandle*>(tobj); } if (tobj->InheritsFrom("CandSubShowerSRListHandle")) { subshowerlist = dynamic_cast<CandSubShowerSRListHandle*>(tobj); } } if (!slicelist || !clusterlist || !subshowerlist) { MSG("AlgShowerSS",Msg::kError) << "CandSliceListHandle, CandClusterListHandle or CandSubShowerSRListHandle missing\n"; return; } CandContext cxx(this,cx.GetMom()); Double_t MinNonPhysicsEnergy = 0.2; //in GeV Int_t MaxPlaneGap = 2; Int_t NPeakFindingBins = 32; //get config for AlgSubShowerSS const char *charShowerSSAlgConfig = 0; ac.Get("ShowerSSAlgConfig",charShowerSSAlgConfig); ac.Get("MaxPlaneGap",MaxPlaneGap); ac.Get("NPeakFindingBins",NPeakFindingBins); ac.Get("MinNonPhysicsEnergy",MinNonPhysicsEnergy); //Get singleton instance of AlgFactory AlgFactory &af = AlgFactory::GetInstance(); AlgHandle ah = af.GetAlgHandle("AlgShowerSS",charShowerSSAlgConfig); const CandRecord *candrec = cx.GetCandRecord(); assert(candrec); const VldContext *vldcptr = candrec->GetVldContext(); assert(vldcptr); VldContext vldc = *vldcptr; UgliGeomHandle ugh(vldc); CandSliceHandleItr sliceItr(slicelist->GetDaughterIterator()); Int_t nslice = 0; Int_t nsubshower = 0; while (CandSliceHandle *slice = sliceItr()) { //slice loop MSG("AlgShowerSS",Msg::kDebug) << "Slice: " << nslice << endl; // Select SubShowers within this slice CandSubShowerSRHandleItr subshowerItr(subshowerlist->GetDaughterIterator()); // iterate over subshowers, and make a shower out of // all subshowers in each slice, in contiguous planes int plane[500] = {0}; int firstPlane = 500; int lastPlane = 0; nsubshower = 0; while (CandSubShowerSRHandle *subshower = subshowerItr()) { if (*subshower->GetCandSlice()==*slice) { MSG("AlgShowerSS",Msg::kDebug) << "SubShower: " << nsubshower << " in slice " << nslice << endl; for(int i=subshower->GetBegPlane();i<=subshower->GetEndPlane();i++){ MSG("AlgShowerSS",Msg::kDebug) << "Plane: " << i << endl; plane[i] = 1; if(i<firstPlane) firstPlane = i; if(i>lastPlane) lastPlane = i; } } nsubshower++; } MSG("AlgShowerSS",Msg::kDebug) << "First Plane = " << firstPlane << endl; MSG("AlgShowerSS",Msg::kDebug) << "Last Plane = " << lastPlane << endl; //separate up into plane chunks according to MaxPlaneGap Int_t nshower = 0; while(firstPlane<=lastPlane){ //while loop over planes Int_t begPlane = 500; Int_t endPlane = 0; Int_t nGaps = 0; for(int i=firstPlane;i<=lastPlane;i++){ if(plane[i]==0) { if(begPlane!=500) nGaps+=1; } else { if(begPlane==500) begPlane = i; endPlane = i; } plane[i] = 0; if(nGaps>MaxPlaneGap) { MSG("AlgShowerSS",Msg::kDebug) << "begPlane = " << begPlane << endl; MSG("AlgShowerSS",Msg::kDebug) << "endPlane = " << endPlane << endl; break; } } firstPlane = endPlane+1; //could have multiple showers per plane chunk //pull out all subshowers in plane chunk and look for TPos peaks: //loop through subshower list again and histogram tposvtx //in order to look for peaks // (add an extra bin at each end outside of expected range, // in order to be able to detect peaks at (far) detector edges) TH1F *vtxUHistAll = new TH1F("vtxUHistAll","vtxUHistAll",NPeakFindingBins+2, -4. - 8./Float_t(NPeakFindingBins), +4. + 8./Float_t(NPeakFindingBins)); TH1F *vtxVHistAll = new TH1F("vtxVHistAll","vtxVHistAll",NPeakFindingBins+2, -4. - 8./Float_t(NPeakFindingBins), +4. + 8./Float_t(NPeakFindingBins)); TH1F *vtxUHist = new TH1F("vtxUHist","vtxUHist",NPeakFindingBins+2, -4. - 8./Float_t(NPeakFindingBins), +4. + 8./Float_t(NPeakFindingBins)); TH1F *vtxVHist = new TH1F("vtxVHist","vtxVHist",NPeakFindingBins+2, -4. - 8./Float_t(NPeakFindingBins), +4. + 8./Float_t(NPeakFindingBins)); subshowerItr.ResetFirst(); nsubshower = 0; while (CandSubShowerSRHandle *subshower = subshowerItr()) { if (*subshower->GetCandSlice()==*slice) { if(subshower->GetBegPlane()>=begPlane && subshower->GetBegPlane()<=endPlane){ if(subshower->GetPlaneView()==PlaneView::kU || subshower->GetPlaneView()==PlaneView::kX) { Float_t weightedCentre = subshower->GetVtxU() + subshower->GetSlope() * PITCHINMETRES * Float_t(subshower->GetEndPlane() - subshower->GetBegPlane() + 1)/2.; if(weightedCentre<-4. || weightedCentre>4.) weightedCentre = subshower->GetVtxU(); MSG("AlgShowerSS",Msg::kDebug) << "Filling vtxUHist with Subshower " << nsubshower << " with vtvU = " << subshower->GetVtxU() << " and energy = " << subshower->GetEnergy() << " and slope = " << subshower->GetSlope() << " and AvgDev = " << subshower->GetAvgDev() << " and weightedCentre = " << weightedCentre << endl; CandStripHandleItr stripssItr(subshower->GetDaughterIterator()); while (CandStripHandle *stp = stripssItr()) { Float_t val = stp->GetCharge(CalDigitType::kSigCorr) - vtxUHistAll->GetBinContent(vtxUHistAll->FindBin(stp->GetTPos())); if(val>0) vtxUHistAll->Fill(stp->GetTPos(),val); //vtxUHistAll->Fill(weightedCentre); if( subshower->GetClusterID()==ClusterType::kEMLike || subshower->GetClusterID()==ClusterType::kHadLike || subshower->GetClusterID()==ClusterType::kTrkLike){ val = stp->GetCharge(CalDigitType::kSigCorr) - vtxUHist->GetBinContent(vtxUHist->FindBin(stp->GetTPos())); if(val>0) vtxUHist->Fill(stp->GetTPos(),val); //vtxUHist->Fill(weightedCentre); } } } else if(subshower->GetPlaneView()==PlaneView::kV || subshower->GetPlaneView()==PlaneView::kY){ Float_t weightedCentre = subshower->GetVtxV() + subshower->GetSlope() * PITCHINMETRES * Float_t(subshower->GetEndPlane() - subshower->GetBegPlane() + 1)/2; if(weightedCentre<-4. || weightedCentre>4.) weightedCentre = subshower->GetVtxV(); MSG("AlgShowerSS",Msg::kDebug) << "Filling vtxVHist with Subshower " << nsubshower << " with vtvV = " << subshower->GetVtxV() << " and energy = " << subshower->GetEnergy() << " and slope = " << subshower->GetSlope() << " and AvgDev = " << subshower->GetAvgDev() << " and weightedCentre = " << weightedCentre << endl; CandStripHandleItr stripssItr(subshower->GetDaughterIterator()); while (CandStripHandle *stp = stripssItr()) { //fill with max charge for this tpos Float_t val = stp->GetCharge(CalDigitType::kSigCorr) - vtxVHistAll->GetBinContent(vtxVHistAll->FindBin(stp->GetTPos())); if(val>0) vtxVHistAll->Fill(stp->GetTPos(),val); //vtxVHistAll->Fill(weightedCentre); if( subshower->GetClusterID()==ClusterType::kEMLike || subshower->GetClusterID()==ClusterType::kHadLike || subshower->GetClusterID()==ClusterType::kTrkLike){ val = stp->GetCharge(CalDigitType::kSigCorr) - vtxVHist->GetBinContent(vtxVHist->FindBin(stp->GetTPos())); if(val>0) vtxVHist->Fill(stp->GetTPos(),val); //vtxVHist->Fill(weightedCentre); } } } } } nsubshower++; } if(false){ TCanvas *can = new TCanvas(); can->Divide(2,2); can->cd(1); vtxUHistAll->Draw(); can->cd(2); vtxVHistAll->Draw(); can->cd(3); vtxUHist->Draw(); can->cd(4); vtxVHist->Draw(); char nomus[256]; sprintf(nomus,"peakPlot_slice%i.eps",nslice); can->Print(nomus); delete can; } TSpectrum *specUAll = new TSpectrum(10); specUAll->Search(vtxUHistAll,1); TSpectrum *specU = new TSpectrum(10); specU->Search(vtxUHist,1); TSpectrum *specVAll = new TSpectrum(10); specVAll->Search(vtxVHistAll,1); TSpectrum *specV = new TSpectrum(10); specV->Search(vtxVHist,1); Int_t nUPeaksAll = specUAll->GetNPeaks(); Int_t nUPeaks = specU->GetNPeaks(); Float_t *peakUPos = specU->GetPositionX(); Float_t *peakUVal = specU->GetPositionY(); Int_t nVPeaksAll = specVAll->GetNPeaks(); Int_t nVPeaks = specV->GetNPeaks(); Float_t *peakVPos = specV->GetPositionX(); Float_t *peakVVal = specV->GetPositionY(); //if no peaks in either view if(nUPeaksAll==0 || nVPeaksAll==0) { //can't form a 3D shower in these conditions! MSG("AlgShowerSS",Msg::kDebug) << "nUPeaks = " << nUPeaks << " nVPeaks = " << nVPeaks << " Can't form a 3D shower" << endl; delete vtxUHist; delete vtxVHist; delete specU; delete specV; delete vtxUHistAll; delete vtxVHistAll; delete specUAll; delete specVAll; continue; } //one or both views have no good "physics" subshower: if(nUPeaks==0 || nVPeaks==0){ MSG("AlgShowerSS",Msg::kDebug) << "Have something in both views but " << "nUPeaks = " << nUPeaks << "and nVPeaks = " << nVPeaks << endl; //try to form anyway with whatever's there if(nUPeaks==0) nUPeaks = 1; if(nVPeaks==0) nVPeaks = 1; } //Try to form some 3D showers, yeah! //have one or more peak => maybe one or //more shower in this plane range //get boundaries that define each shower in each view: Float_t *peakUPosLow = new Float_t[nUPeaks]; Float_t *peakUPosUpp = new Float_t[nUPeaks]; MSG("AlgShowerSS",Msg::kDebug) << "Number of U Peaks found = " << nUPeaks << " with positions, heights " << "and bounds: " << endl; for (int i=0;i<nUPeaks;i++){ if(i==0) peakUPosLow[i] = -4; else peakUPosLow[i] = peakUPosUpp[i-1]; if(i==nUPeaks-1) peakUPosUpp[i] = 4; else { Int_t bin = vtxUHist->GetXaxis()->FindBin(peakUPos[i]); Int_t bestBin = vtxUHist->GetNbinsX(); Float_t threeBinAve = 1000000; //sigcor! while(bin<vtxUHist->GetNbinsX() && vtxUHist->GetBinCenter(bin)<peakUPos[i+1]) { Float_t ave = (vtxUHist->GetBinContent(bin-1) + vtxUHist->GetBinContent(bin) + vtxUHist->GetBinContent(bin+1)); if(ave<threeBinAve) { threeBinAve = ave; bestBin = bin; } bin++; } peakUPosUpp[i] = vtxUHist->GetBinCenter(bestBin); } MSG("AlgShowerSS",Msg::kDebug) << peakUPos[i] << " " << peakUVal[i] << " [" << peakUPosLow[i] << "," << peakUPosUpp[i] << "]" << endl; } Float_t *peakVPosLow = new Float_t[nVPeaks]; Float_t *peakVPosUpp = new Float_t[nVPeaks]; MSG("AlgShowerSS",Msg::kDebug) << "Number of V Peaks found = " << nVPeaks << " with positions, heights " << "and bounds: " << endl; for (int i=0;i<nVPeaks;i++){ if(i==0) peakVPosLow[i] = -4; else peakVPosLow[i] = peakVPosUpp[i-1]; if(i==nVPeaks-1) peakVPosUpp[i] = 4; else { Int_t bin = vtxVHist->GetXaxis()->FindBin(peakVPos[i]); Int_t bestBin = vtxVHist->GetNbinsX(); Float_t threeBinAve = 1000000; //gev! while(bin<vtxVHist->GetNbinsX() && vtxVHist->GetBinCenter(bin)<peakVPos[i+1]) { Float_t ave = (vtxVHist->GetBinContent(bin-1) + vtxVHist->GetBinContent(bin) + vtxVHist->GetBinContent(bin+1)); if(ave<threeBinAve) { threeBinAve = ave; bestBin = bin; } bin++; } peakVPosUpp[i] = vtxVHist->GetBinCenter(bestBin); } MSG("AlgShowerSS",Msg::kDebug) << peakVPos[i] << " " << peakVVal[i] << " [" << peakVPosLow[i] << "," << peakVPosUpp[i] << "]" << endl; } //for matching up U/V views of showers: Float_t *totUEnergy = new Float_t[nUPeaks]; Float_t *approxVfromTime = new Float_t[nUPeaks]; Float_t *maxPhysU = new Float_t[nUPeaks]; Float_t *minPhysU = new Float_t[nUPeaks]; Int_t *maxPlaneU = new Int_t[nUPeaks]; Int_t *minPlaneU = new Int_t[nUPeaks]; Int_t *nSubShowersU = new Int_t[nUPeaks]; for(int i=0;i<nUPeaks;i++) { totUEnergy[i] = 0; approxVfromTime[i] = 0; maxPhysU[i] = -1; maxPlaneU[i] = -1; minPhysU[i] = 999; minPlaneU[i] = -1; nSubShowersU[i] = 0; } Float_t *totVEnergy = new Float_t[nVPeaks]; Float_t *approxUfromTime = new Float_t[nVPeaks]; Float_t *maxPhysV = new Float_t[nVPeaks]; Float_t *minPhysV = new Float_t[nVPeaks]; Int_t *maxPlaneV = new Int_t[nVPeaks]; Int_t *minPlaneV = new Int_t[nVPeaks]; Int_t *nSubShowersV = new Int_t[nVPeaks]; for(int i=0;i<nVPeaks;i++) { totVEnergy[i] = 0; approxUfromTime[i] = 0; maxPhysV[i] = -1; maxPlaneV[i] = -1; minPhysV[i] = 999; minPlaneV[i] = -1; nSubShowersV[i] = 0; } //loop through subshowerlist and add up energy of showers subshowerItr.ResetFirst(); while (CandSubShowerSRHandle *subshower = subshowerItr()) { if (*subshower->GetCandSlice()==*slice) { if(subshower->GetBegPlane()>=begPlane && subshower->GetEndPlane()<=endPlane) { if(subshower->GetPlaneView()==PlaneView::kU){ for(int i=0;i<nUPeaks;i++){ Float_t weightedCentre = subshower->GetVtxU() + subshower->GetSlope() * PITCHINMETRES * Float_t(subshower->GetEndPlane() - subshower->GetBegPlane() + 1)/2.; if(weightedCentre>peakUPosLow[i] && weightedCentre<=peakUPosUpp[i]) { totUEnergy[i]+=subshower->GetEnergy(); if(nUPeaks>1 && subshower->GetClusterID()==ClusterType::kEMLike || subshower->GetClusterID()==ClusterType::kHadLike || subshower->GetClusterID()==ClusterType::kTrkLike){ approxVfromTime[i]+=subshower->GetVtxV(); nSubShowersU[i] += 1; for(int j=subshower->GetBegPlane();j<=subshower->GetEndPlane();j++){ Float_t testMinU = subshower->GetMinU(j); Float_t testMaxU = subshower->GetMaxU(j); if(testMaxU>maxPhysU[i]) { maxPhysU[i] = testMaxU; maxPlaneU[i] = j; } if(testMinU<minPhysU[i]) { minPhysU[i] = testMinU; minPlaneU[i] = j; } } } } } } else if(subshower->GetPlaneView()==PlaneView::kV){ for(int i=0;i<nVPeaks;i++){ Float_t weightedCentre = subshower->GetVtxV() + subshower->GetSlope() * PITCHINMETRES * Float_t(subshower->GetEndPlane() - subshower->GetBegPlane() + 1)/2.; if(weightedCentre>peakVPosLow[i] && weightedCentre<=peakVPosUpp[i]) { totVEnergy[i]+=subshower->GetEnergy(); if(nVPeaks>1 && subshower->GetClusterID()==ClusterType::kEMLike || subshower->GetClusterID()==ClusterType::kHadLike || subshower->GetClusterID()==ClusterType::kTrkLike){ approxUfromTime[i]+=subshower->GetVtxU(); nSubShowersV[i] += 1; for(int j=subshower->GetBegPlane();j<=subshower->GetEndPlane();j++){ Float_t testMinV = subshower->GetMinV(j); Float_t testMaxV = subshower->GetMaxV(j); if(testMaxV>maxPhysV[i]) { maxPhysV[i] = testMaxV; maxPlaneV[i] = j; } if(testMinV<minPhysV[i]) { minPhysV[i] = testMinV; minPlaneV[i] = j; } } } } } } } } } //Test to see if any peaks should be combined: //Look for adjacent peaks for which min - max "physics" strip difference is ~1 for(int i=0;i<nUPeaks-1;i++){ if(minPlaneU[i+1]==-1 || maxPlaneU[i]==-1 || TMath::Abs(minPlaneU[i+1] - maxPlaneU[i]) > 3 ) //within +/- 1 plane in view continue; if(minPhysU[i+1] - maxPhysU[i] < 1.5*STRIPWIDTHINMETRES ) {//combine peaks //add energy to second peak and //remove energy from first peak totUEnergy[i+1] += totUEnergy[i]; totUEnergy[i] = 0; //add other-view vertex estimate sum approxVfromTime[i+1] += approxVfromTime[i]; approxVfromTime[i] = 0; nSubShowersU[i+1] += nSubShowersU[i]; nSubShowersU[i] = 0; //shift min U of second peak down and //set min U plane of second peak minPhysU[i+1] = minPhysU[i]; minPlaneU[i+1] = minPlaneU[i]; //set min/max limits of first peak to be the same maxPhysU[i] = minPhysU[i]; maxPlaneU[i] = minPlaneU[i]; //set peak lower limit of second peak to that of first peakUPosLow[i+1] = peakUPosLow[i]; //set peak lower and upper limit of first peak to be the same peakUPosUpp[i] = peakUPosLow[i]; } } for(int i=0;i<nVPeaks-1;i++){ if(minPlaneV[i+1]==-1 || maxPlaneV[i]==-1 || TMath::Abs(minPlaneV[i+1] - maxPlaneV[i]) > 3 ) //within +/- 1 plane in view continue; if(minPhysV[i+1] - maxPhysV[i] < 1.5*STRIPWIDTHINMETRES ) {//combine peaks //add energy to second peak and //remove energy from first peak totVEnergy[i+1] += totVEnergy[i]; totVEnergy[i] = 0; //add other-view vertex estimate sum approxUfromTime[i+1] += approxUfromTime[i]; approxUfromTime[i] = 0; nSubShowersV[i+1] += nSubShowersV[i]; nSubShowersV[i] = 0; //shift min V of second peak down and //set min V plane of second peak minPhysV[i+1] = minPhysV[i]; minPlaneV[i+1] = minPlaneV[i]; //set min/max limits of first peak to be the same maxPhysV[i] = minPhysV[i]; maxPlaneV[i] = minPlaneV[i]; //set peak lower limit of second peak to that of first peakVPosLow[i+1] = peakVPosLow[i]; //set peak lower and upper limit of first peak to be the same peakVPosUpp[i] = peakVPosLow[i]; } } //get best order for matching subshower groups in U/V Int_t *bestOrderU = new Int_t[nUPeaks]; Int_t *bestOrderV = new Int_t[nVPeaks]; for(int i=0;i<nUPeaks;i++){ if(totUEnergy[i]>0){ Int_t nbigger = 0; for(int j=0;j<nUPeaks;j++){ if(totUEnergy[i]<totUEnergy[j]) nbigger+=1; } bestOrderU[nbigger] = i; } else { // if there is no energy in this "peak" Int_t nsmaller = 0; // put zero peaks in reverse order for(int j=0;j<i;j++){ // purely to ensure that all indices if(totUEnergy[j]==0) nsmaller+=1; // of bestOrderX are filled } // with legitimate values (prevents bestOrderU[nUPeaks-1-nsmaller] = i; // against a rare case seg fault) } } for(int i=0;i<nVPeaks;i++){ if(totVEnergy[i]>0){ Int_t nbigger = 0; for(int j=0;j<nVPeaks;j++){ if(totVEnergy[i]<totVEnergy[j]) nbigger+=1; } bestOrderV[nbigger] = i; } else { Int_t nsmaller = 0; for(int j=0;j<i;j++){ if(totVEnergy[j]==0) nsmaller+=1; } bestOrderV[nVPeaks-1-nsmaller] = i; } } Int_t nShowers = nUPeaks; if(nShowers>nVPeaks) nShowers = nVPeaks; if(false){ //check energy ordering is ok: for(int i=0;i<nShowers-1;i++){ //loop over showers Float_t EUi = totUEnergy[bestOrderU[i]]; Float_t EVi = totVEnergy[bestOrderV[i]]; Float_t EUii = totUEnergy[bestOrderU[i+1]]; Float_t EVii = totVEnergy[bestOrderV[i+1]]; if(EUi==0 || EUii==0 || EVi==0 || EVii==0) continue; //can we swap the order and still get two good energy matches? Float_t asym = TMath::Abs(EUi - EVi)/(EUi + EVi) + TMath::Abs(EUii - EVii)/(EUii + EVii); Float_t asymSwap = TMath::Abs(EUi - EVii)/(EUi + EVii) + TMath::Abs(EUii - EVi)/(EUii + EVi); //if difference between two asymmetry estimates is <10% //see if vertex estimate from timing helps if( TMath::Abs(asym - asymSwap)/asym < 0.1) { if(nSubShowersU[bestOrderU[i]]==0 || nSubShowersV[bestOrderV[i]]==0 || nSubShowersU[bestOrderU[i+1]]==0 || nSubShowersV[bestOrderV[i+1]]==0) continue; approxVfromTime[bestOrderU[i]] /= nSubShowersU[bestOrderU[i]]; approxUfromTime[bestOrderV[i]] /= nSubShowersV[bestOrderV[i]]; approxVfromTime[bestOrderU[i+1]] /= nSubShowersU[bestOrderU[i+1]]; approxUfromTime[bestOrderV[i+1]] /= nSubShowersV[bestOrderV[i+1]]; //sum distance from both showers Float_t vertexDistance = TMath::Sqrt(TMath::Power(approxUfromTime[bestOrderV[i]] - peakUPos[bestOrderU[i]],2.) + TMath::Power(approxVfromTime[bestOrderU[i]] - peakVPos[bestOrderV[i]],2)) + TMath::Sqrt(TMath::Power(approxUfromTime[bestOrderV[i+1]] - peakUPos[bestOrderU[i+1]],2.) + TMath::Power(approxVfromTime[bestOrderU[i+1]] - peakVPos[bestOrderV[i+1]],2)); Float_t vertexDistanceSwap = TMath::Sqrt(TMath::Power(approxUfromTime[bestOrderV[i+1]] - peakUPos[bestOrderU[i]],2.) + TMath::Power(approxVfromTime[bestOrderU[i]] - peakVPos[bestOrderV[i+1]],2)) + TMath::Sqrt(TMath::Power(approxUfromTime[bestOrderV[i]] - peakUPos[bestOrderU[i+1]],2.) + TMath::Power(approxVfromTime[bestOrderU[i+1]] - peakVPos[bestOrderV[i]],2)); if(vertexDistanceSwap<vertexDistance && vertexDistanceSwap<0.5) { //check that best vertex is reasonable //swap! //arbitrary choice which view to switch: Int_t tempIndex = bestOrderU[i]; bestOrderU[i] = bestOrderU[i+1]; bestOrderU[i+1] = tempIndex; //start checks again so that clusters //can propagate if better matches are found i=0; } } } } //get true number of physics peaks again //for final test nUPeaks = specU->GetNPeaks(); nVPeaks = specV->GetNPeaks(); //loop over number of expected showers: for(int i=0;i<nShowers;i++){ TObjArray newshower; Int_t nsubshower = 0; subshowerItr.ResetFirst(); while (CandSubShowerSRHandle *subshower = subshowerItr()) { if (*subshower->GetCandSlice()==*slice) { if(subshower->GetBegPlane()>=begPlane && subshower->GetEndPlane()<=endPlane) { if(subshower->GetPlaneView()==PlaneView::kU){ Float_t weightedCentre = subshower->GetVtxU() + subshower->GetSlope() * PITCHINMETRES * Float_t(subshower->GetEndPlane() - subshower->GetBegPlane() + 1)/2.; if(weightedCentre>peakUPosLow[bestOrderU[i]] && weightedCentre<=peakUPosUpp[bestOrderU[i]]) { newshower.Add(subshower); MSG("AlgShowerSS",Msg::kDebug) << "Adding SubShower: " << nsubshower << " in slice " << nslice << " to newshower " << i << endl; } } if(subshower->GetPlaneView()==PlaneView::kV){ Float_t weightedCentre = subshower->GetVtxV() + subshower->GetSlope() * PITCHINMETRES * Float_t(subshower->GetEndPlane() - subshower->GetBegPlane() + 1)/2.; if(weightedCentre>peakVPosLow[bestOrderV[i]] && weightedCentre<=peakVPosUpp[bestOrderV[i]]) { newshower.Add(subshower); MSG("AlgShowerSS",Msg::kDebug) << "Adding SubShower: " << nsubshower << " in slice " << nslice << " to newshower " << i << endl; } } } } nsubshower+=1; } cxx.SetDataIn(&newshower); CandShowerSRHandle showerhandle = CandShowerSR::MakeCandidate(ah,cxx); if((showerhandle.GetEnergy()>0 && nUPeaks>0 && nVPeaks>0) || showerhandle.GetEnergy()>MinNonPhysicsEnergy) { showerhandle.SetCandSlice(slice); ch.AddDaughterLink(showerhandle); } else { MSG("AlgShowerSS",Msg::kDebug) << "Shower " << nshower << " has zero energy." << " - Not adding shower to list." << endl; } } delete [] peakUPosUpp; delete [] peakUPosLow; delete [] peakVPosUpp; delete [] peakVPosLow; delete [] totUEnergy; delete [] totVEnergy; delete [] approxVfromTime; delete [] approxUfromTime; delete [] nSubShowersU; delete [] nSubShowersV; delete [] maxPhysU; delete [] maxPhysV; delete [] minPhysU; delete [] minPhysV; delete [] maxPlaneU; delete [] maxPlaneV; delete [] minPlaneU; delete [] minPlaneV; delete [] bestOrderU; delete [] bestOrderV; delete vtxUHist; delete vtxVHist; delete vtxUHistAll; delete vtxVHistAll; delete specU; delete specV; delete specUAll; delete specVAll; } nslice++; } }

void AlgShowerSSList::Trace (  const char *  c  )  const [virtual]

Reimplemented from AlgShowerSRList. Definition at line 755 of file AlgShowerSSList.cxx. { }

---

The documentation for this class was generated from the following files:

• AlgShowerSSList.h
• AlgShowerSSList.cxx

---