AlgShowerAtNu Class Reference

#include <AlgShowerAtNu.h>

Inheritance diagram for AlgShowerAtNu:

List of all members.

Public Member Functions
	AlgShowerAtNu ()
	~AlgShowerAtNu ()
void	RunAlg (AlgConfig &ac, CandHandle &ch, CandContext &cx)
void	Trace (const char *c) const
Private Attributes
TObjArray	fHitArr [500]

---

Constructor & Destructor Documentation

AlgShowerAtNu::AlgShowerAtNu (   )

Definition at line 34 of file AlgShowerAtNu.cxx. { }

AlgShowerAtNu::~AlgShowerAtNu (   )

Definition at line 39 of file AlgShowerAtNu.cxx. { }

---

Member Function Documentation

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

Implements AlgBase. Definition at line 44 of file AlgShowerAtNu.cxx. References CandHandle::AddDaughterLink(), CandRecoHandle::CalibrateMIP(), CandRecoHandle::CalibrateSigMapped(), UgliStripHandle::ClearFiber(), fHitArr, FloatErr, ShowerAtNu::GetAssocTrackAt(), Calibrator::GetAttenCorrectedTpos(), ObjTrackAtNu::GetBegPlane(), ObjShowerAtNu::GetBegPlane(), ShowerAtNu::GetBegVtxFlag(), CandContext::GetCandRecord(), ShowerAtNu::GetCandSliceHandle(), HitAtNu::GetCandStripHandle(), CandStripHandle::GetCharge(), HitAtNu::GetCharge(), CandContext::GetDataIn(), VldContext::GetDetector(), CandShowerAtNuHandle::GetDirCosErrU(), CandShowerAtNuHandle::GetDirCosErrV(), CandRecoHandle::GetDirCosU(), CandRecoHandle::GetDirCosV(), CandRecoHandle::GetDirCosZ(), CandShowerAtNuHandle::GetDirTimeScore(), CandShowerAtNuHandle::GetDirVtxShwScore(), Registry::GetDouble(), PlexStripEndId::GetEncoded(), ObjTrackAtNu::GetEndPlane(), ObjShowerAtNu::GetEndPlane(), ShowerAtNu::GetEndVtxFlag(), CandShowerHandle::GetEnergy(), UgliStripHandle::GetHalfLength(), ObjAtNu::GetHitAt(), ObjAtNu::GetHitLast(), Registry::GetInt(), CandShowerAtNuHandle::GetMaxPlane(), CandShowerAtNuHandle::GetMinPlane(), Calibrator::GetMIP(), CandStripHandle::GetNDigit(), CandShowerAtNuHandle::GetNPlanes(), CandShowerAtNuHandle::GetNStrips(), HitAtNu::GetPlane(), CandStripHandle::GetPlaneView(), HitAtNu::GetPlaneView(), CandShowerAtNuHandle::GetReseedFlag(), ObjShowerAtNu::GetReseedFlag(), HitAtNu::GetShwFlag(), UgliGeomHandle::GetSteelPlnHandle(), CandStripHandle::GetStripEndId(), UgliGeomHandle::GetStripHandle(), CandStripHandle::GetTime(), CandRecoHandle::GetTimeOffset(), CandRecoHandle::GetTimeSlope(), HitAtNu::GetTPos(), RecMinos::GetVldContext(), CandRecoHandle::GetVtxPlane(), CandShowerAtNuHandle::GetVtxR(), CandShowerAtNuHandle::GetVtxShw(), CandRecoHandle::GetVtxT(), CandRecoHandle::GetVtxU(), CandRecoHandle::GetVtxV(), CandRecoHandle::GetVtxZ(), UgliSteelPlnHandle::GetZ0(), HitAtNu::GetZPos(), Calibrator::Instance(), UgliSteelPlnHandle::IsValid(), MSG, Calibrator::ReInitialise(), CandRecoHandle::SetCandSlice(), CandShowerAtNuHandle::SetDirCosErrU(), CandShowerAtNuHandle::SetDirCosErrV(), CandRecoHandle::SetDirCosU(), CandRecoHandle::SetDirCosV(), CandRecoHandle::SetDirCosZ(), CandShowerAtNuHandle::SetDirTimeScore(), CandShowerAtNuHandle::SetDirVtxShwScore(), CandShowerHandle::SetEnergy(), ValueErr< T >::SetError(), CandShowerAtNuHandle::SetMaxPlane(), CandShowerAtNuHandle::SetMinPlane(), CandShowerAtNuHandle::SetNPlanes(), CandShowerAtNuHandle::SetNStrips(), CandShowerAtNuHandle::SetReseedFlag(), CandRecoHandle::SetTimeOffset(), CandRecoHandle::SetTimeSlope(), CandShowerHandle::SetU(), CandShowerHandle::SetV(), CandRecoHandle::SetVtxPlane(), CandShowerAtNuHandle::SetVtxR(), CandShowerAtNuHandle::SetVtxShw(), CandRecoHandle::SetVtxT(), CandRecoHandle::SetVtxU(), CandRecoHandle::SetVtxV(), CandRecoHandle::SetVtxZ(), and UgliStripHandle::WlsPigtail(). { MSG("AlgShowerAtNu", Msg::kDebug) << "AlgShowerAtNu::RunAlg(...)" << endl; // Get CandShowerAtNuHandle CandShowerAtNuHandle& shower = dynamic_cast<CandShowerAtNuHandle&>(ch); // Unpack AlgConfig Double_t fFibreIndex; Int_t fAtNuAnaOnOff; fFibreIndex = ac.GetDouble("FibreIndex"); fAtNuAnaOnOff = ac.GetInt("AtNuAnaOnOff"); MSG("AlgShowerAtNu", Msg::kDebug) << " FibreIndex=" << fFibreIndex << endl; MSG("AlgShowerAtNu", Msg::kDebug) << " AtNuAnaOnOff=" << fAtNuAnaOnOff << endl; // Unpack CandContext const TObjArray *arr = dynamic_cast<const TObjArray*>(cx.GetDataIn()); ShowerAtNu* shwu = (ShowerAtNu*)(arr->At(0)); ShowerAtNu* shwv = (ShowerAtNu*)(arr->At(1)); shower.SetCandSlice(shwu->GetCandSliceHandle()); CandRecord* candrec = (CandRecord*)(cx.GetCandRecord()); VldContext *vldc = (VldContext*)(candrec->GetVldContext()); UgliGeomHandle ugh(*vldc); // Reset Calibrator Calibrator& cal = Calibrator::Instance(); cal.ReInitialise(*vldc); // Get Geometry Int_t atmosflag,modtype,modnum,modpln,modstr; switch(candrec->GetVldContext()->GetDetector()){ case Detector::kFar: atmosflag=1; modtype=1; modnum=2; modpln=248; modstr=192; break; case Detector::kNear: atmosflag=0; modtype=2; modnum=1; modpln=282; modstr=96; break; case Detector::kCalib: atmosflag=0; modtype=3; modnum=1; modpln=60; modstr=24; break; default: atmosflag=0; modtype=-1; modnum=0; modpln=0; modstr=0; break; } Int_t i,j,k; Int_t vuw; Int_t pln,bpln,epln,npln; Double_t begZsm1=-9999.9,endZsm1=-9999.9; Double_t begZsm2=-9999.9,endZsm2=-9999.9; Int_t begsm1=-1,endsm1=-1,begsm2=-1,endsm2=-1; for(j=0;j<486;j++){ PlexPlaneId myplaneid(vldc->GetDetector(),j,1); UgliSteelPlnHandle myplanehandle = ugh.GetSteelPlnHandle(myplaneid); if(myplanehandle.IsValid()){ if(j<=modpln && begZsm1<0.0){ begsm1=j; begZsm1=myplanehandle.GetZ0(); } if(j<=modpln && myplanehandle.GetZ0()>endZsm1){ endsm1=j; endZsm1=myplanehandle.GetZ0(); } if(j>modpln && begZsm2<0.0){ begsm2=j; begZsm2=myplanehandle.GetZ0(); } if(j>modpln && myplanehandle.GetZ0()>endZsm2){ endsm2=j; endZsm2=myplanehandle.GetZ0(); } } } MSG("AlgShowerAtNu", Msg::kDebug) << " *** detector geometry *** " << endl << " SM1 : pln " << begsm1 << " -> " << endsm1 << " Z " << begZsm1 << " -> " << endZsm1 << endl << " SM2 : pln " << begsm2 << " -> " << endsm2 << " Z " << begZsm2 << " -> " << endZsm2 << endl; if(shwu->GetBegPlane()<shwv->GetBegPlane()) bpln=shwu->GetBegPlane(); else bpln=shwv->GetBegPlane(); if(shwu->GetEndPlane()>shwv->GetEndPlane()) epln=shwu->GetEndPlane(); else epln=shwv->GetEndPlane(); // *********************************************** // * C A L C U L A T E C O - O R D I N A T E S * // *********************************************** MSG("AlgShowerAtNu", Msg::kDebug) << " *** calculate co-ordinates *** " << endl; npln = 1+epln-bpln; Double_t* Z = new Double_t[npln]; Double_t* U = new Double_t[npln]; Double_t* V = new Double_t[npln]; Double_t* CT = new Double_t[npln]; Double_t* Q = new Double_t[npln]; Double_t* Qm = new Double_t[npln]; Double_t* Qp = new Double_t[npln]; Double_t* CTm = new Double_t[npln]; Double_t* CTp = new Double_t[npln]; Int_t* plnvuw = new Int_t[npln]; Int_t* plnnum = new Int_t[npln]; for(i=0;i<npln;i++){ plnvuw[i]=-1; plnnum[i]=0; U[i]=0.0; V[i]=0.0; Z[i]=0.0; Q[i]=0.0; CT[i]=0.0; Qm[i]=0.0; Qp[i]=0.0; CTm[i]=0.0; CTp[i]=0.0; } // loop over hits and fill arrays for(i=0;i<1+shwu->GetHitLast();i++){ HitAtNu* hit = (HitAtNu*)(shwu->GetHitAt(i)); CandStripHandle* strip = (CandStripHandle*)(hit->GetCandStripHandle()); pln = hit->GetPlane()-bpln; fHitArr[pln].Add(hit); if(plnvuw[pln]<0){ plnvuw[pln]=hit->GetPlaneView(); Z[pln]=hit->GetZPos(); } U[pln]+=hit->GetTPos()*hit->GetCharge(); Q[pln]+=hit->GetCharge(); shower.AddDaughterLink(*strip); plnnum[pln]++; } for(i=0;i<1+shwv->GetHitLast();i++){ HitAtNu* hit = (HitAtNu*)(shwv->GetHitAt(i)); CandStripHandle* strip = (CandStripHandle*)(hit->GetCandStripHandle()); pln = hit->GetPlane()-bpln; fHitArr[pln].Add(hit); if(plnvuw[pln]<0){ plnvuw[pln]=hit->GetPlaneView(); Z[pln]=hit->GetZPos(); } V[pln]+=hit->GetTPos()*hit->GetCharge(); Q[pln]+=hit->GetCharge(); shower.AddDaughterLink(*strip); plnnum[pln]++; } Double_t totph=0.; Int_t nplanes=0,nstrips=0; for(i=0;i<npln;i++){ if(plnvuw[i]>-1){ if(Q[i]>0.0){ U[i]=U[i]/Q[i]; V[i]=V[i]/Q[i]; totph+=Q[i]; nstrips+=plnnum[i]; nplanes++; } } } // estimate transverse co-ordinates Int_t km,kp; Double_t q,sq,sqt; for(pln=0;pln<npln;pln++){ if(plnvuw[pln]>-1){ vuw=plnvuw[pln]; sq=0.0; sqt=0.0; k=0; km=-1; while(pln-k>0 && km<0){ k++; if(plnvuw[pln-k]>-1 && plnvuw[pln-k]!=vuw) km=k; } k=0; kp=-1; while(pln+k<npln-1 && kp<0){ k++; if(plnvuw[pln+k]>-1 && plnvuw[pln+k]!=vuw) kp=k; } if(km>0){ if(vuw==0){ q=1.0/km; sq+=q*Q[pln-km]; sqt+=q*Q[pln-km]*V[pln-km]; } if(vuw==1){ q=1.0/km; sq+=q*Q[pln-km]; sqt+=q*Q[pln-km]*U[pln-km]; } } if(kp>0){ if(vuw==0){ q=1.0/kp; sq+=q*Q[pln+kp]; sqt+=q*Q[pln+kp]*V[pln+kp]; } if(vuw==1){ q=1.0/kp; sq+=q*Q[pln+kp]; sqt+=q*Q[pln+kp]*U[pln+kp]; } } if(sq>0.0){ if(vuw==0){ V[pln]=sqt/sq; } if(vuw==1){ U[pln]=sqt/sq; } } else{ if(vuw==0){ V[pln]=0.0; } if(vuw==1){ U[pln]=0.0; } } } } // get timing information Double_t tpos = 0., fibre = 0.; Double_t ctime = 0., digchg = 0.; Double_t indx=fFibreIndex; for(pln=0;pln<npln;pln++){ for(i=0;i<1+fHitArr[pln].GetLast();i++){ HitAtNu* hit = (HitAtNu*)(fHitArr[pln].At(i)); CandStripHandle* strip = (CandStripHandle*)(hit->GetCandStripHandle()); if( strip->GetPlaneView()==PlaneView::kU || strip->GetPlaneView()==PlaneView::kX ) tpos = -V[pln]; if( strip->GetPlaneView()==PlaneView::kV || strip->GetPlaneView()==PlaneView::kY ) tpos = U[pln]; PlexStripEndId stripid = strip->GetStripEndId(); UgliStripHandle striphandle = ugh.GetStripHandle(stripid); if(strip->GetNDigit(StripEnd::kPositive)>0){ fibre = striphandle.ClearFiber(StripEnd::kPositive)+striphandle.WlsPigtail(StripEnd::kPositive)+striphandle.GetHalfLength()-tpos; ctime = 3.0e8*strip->GetTime(StripEnd::kPositive)-indx*fibre; digchg=strip->GetCharge(StripEnd::kPositive); Qp[pln]+=digchg; CTp[pln]+=digchg*ctime; } if(strip->GetNDigit(StripEnd::kNegative)>0){ fibre = striphandle.ClearFiber(StripEnd::kNegative)+striphandle.WlsPigtail(StripEnd::kNegative)+striphandle.GetHalfLength()+tpos; ctime = 3.0e8*strip->GetTime(StripEnd::kNegative)-indx*fibre; digchg=strip->GetCharge(StripEnd::kNegative); Qm[pln]+=digchg; CTm[pln]+=digchg*ctime; } } } // set map values for(i=0;i<npln;i++){ pln=bpln+i; if(plnvuw[i]>-1){ shower.SetU(pln,U[i]); shower.SetV(pln,V[i]); if(Qm[i]>0.0){ CTm[i]=CTm[i]/Qm[i]; } if(Qp[i]>0.0){ CTp[i]=CTp[i]/Qp[i]; } if(Qm[i]+Qp[i]>0){ CT[i]=(Qm[i]*CTm[i]+Qp[i]*CTp[i])/(Qm[i]+Qp[i]); } } } MSG("AlgShowerAtNu", Msg::kDebug) << " *** shower hits *** " << endl; for(pln=0;pln<npln;pln++){ if( plnvuw[pln]>-1 ){ MSG("AlgShowerAtNu",Msg::kDebug) << " pln=" << pln << " vuw=" << plnvuw[pln] << " z=" << Z[pln] << " u=" << U[pln] << " v=" << V[pln] << endl; } } // *********************************** // * D E T E R M I N E V E R T E X * // *********************************** MSG("AlgShowerAtNu", Msg::kDebug) << " *** determine vertex *** " << endl; Bool_t vtxshw=0; Double_t trkscr=0.0,shwscr=0.0,dirscr=0.0; Double_t vtxu=0.0,vtxv=0.0,vtxx=0.0,vtxy=0.0,vtxz=0.0; Int_t vtxpln=-1,trkbpln=-1,trkepln=-1,shwbpln=-1,shwepln=-1; // determine vertex plane Double_t totchgpln=0.0,totchg=0.0; for(i=0;i<npln;i++){ for(j=0;j<1+fHitArr[i].GetLast();j++){ HitAtNu* hit = (HitAtNu*)(fHitArr[i].At(j)); if(hit->GetShwFlag()>1){ totchg+=hit->GetCharge(); totchgpln+=hit->GetCharge()*hit->GetPlane(); if(shwbpln<0||hit->GetPlane()<shwbpln) shwbpln=hit->GetPlane(); if(shwepln<0||hit->GetPlane()>shwepln) shwepln=hit->GetPlane(); } } } if(totchg>0.0){ vtxpln=(Int_t)(totchgpln/totchg); } else{ vtxpln=(Int_t)(0.5*(bpln+epln)); shwbpln=bpln; shwepln=epln; } TrackAtNu* vtxtrkU = (TrackAtNu*)(shwu->GetAssocTrackAt(0)); TrackAtNu* vtxtrkV = (TrackAtNu*)(shwv->GetAssocTrackAt(0)); if(vtxtrkU && vtxtrkV){ vtxshw=1; if(vtxtrkU->GetBegPlane()<vtxtrkV->GetBegPlane()) trkbpln=vtxtrkU->GetBegPlane(); else trkbpln=vtxtrkV->GetBegPlane(); if(vtxtrkU->GetEndPlane()>vtxtrkV->GetEndPlane()) trkepln=vtxtrkU->GetEndPlane(); else trkepln=vtxtrkV->GetEndPlane(); trkscr = (trkbpln+trkepln-2.0*vtxpln)/(1.0+trkepln-trkbpln); if(trkscr>0.0){ trkscr = (shwbpln+shwepln-2.0*trkbpln)/(1.0+shwepln-shwbpln); } if(trkscr<0.0){ trkscr = (shwbpln+shwepln-2.0*trkepln)/(1.0+shwepln-shwbpln); } if( ( shwu->GetBegVtxFlag() && shwv->GetBegVtxFlag()) && !( shwu->GetEndVtxFlag() && shwv->GetEndVtxFlag()) ) vtxpln=trkbpln; if( ( shwu->GetEndVtxFlag() && shwv->GetEndVtxFlag()) && !( shwu->GetBegVtxFlag() && shwv->GetBegVtxFlag()) ) vtxpln=trkepln; } // determine vertex position pln=vtxpln-bpln; if( pln>-1 && pln<npln && plnvuw[pln]>-1 ){ vtxu=U[pln]; vtxv=V[pln]; vtxz=Z[pln]; } else{ k=0; km=-1; while(pln-k>0 && km<0){ k++; if(pln-k<npln && plnvuw[pln-k]>-1) km=k; } k=0; kp=-1; while(pln+k<npln-1 && kp<0){ k++; if(pln+k>-1 && plnvuw[pln+k]>-1) kp=k; } if(km>0 && kp>0){ vtxu=(km*U[pln+kp]+kp*U[pln-km])/(km+kp); vtxv=(km*V[pln+kp]+kp*V[pln-km])/(km+kp); vtxz=(km*Z[pln+kp]+kp*Z[pln-km])/(km+kp); } else{ if(km>0){ vtxu=U[pln-km]; vtxv=V[pln-km]; vtxz=Z[pln-km]; } if(kp>0){ vtxu=U[pln+kp]; vtxv=V[pln+kp]; vtxz=Z[pln+kp]; } } } vtxx=0.7071*(vtxu-vtxv); vtxy=0.7071*(vtxu+vtxv); MSG("AlgShowerAtNu", Msg::kVerbose) << " ... vertex: " << " plane = " << vtxpln << " position = (" << vtxu << "," << vtxv << "," << vtxz << ")" << endl; // ************************************* // * C O N T A I N M E N T S T U F F * // ************************************* MSG("AlgShowerAtNu", Msg::kDebug) << " *** containment stuff *** " << endl; Double_t vtxr=0.0; if(atmosflag){ // closest distance to edge Double_t xm,xp,ym,yp,um,up,vm,vp; Double_t rmin; rmin=4.0; up=4.0-vtxu; if(up<rmin) rmin=up; um=4.0+vtxu; if(um<rmin) rmin=um; vp=4.0-vtxv; if(vp<rmin) rmin=vp; vm=4.0+vtxv; if(vm<rmin) rmin=vm; xp=4.0-vtxx; if(xp<rmin) rmin=xp; xm=4.0+vtxx; if(xm<rmin) rmin=xm; yp=4.0-vtxy; if(yp<rmin) rmin=yp; ym=4.0+vtxy; if(ym<rmin) rmin=ym; vtxr=rmin; } //***************************************** //* D E T E R M I N E D I R E C T I O N * //***************************************** MSG("AlgShowerAtNu", Msg::kDebug) << " *** determine direction *** " << endl; // measure overall timeslope Double_t ctimeslope=0.0,ctimeoffset=0.0,ctimeaverage=0.0,ctimescatter=0.0; Double_t sw,swz,swt,swzz,swzt,swtt; Double_t du,dv,dz; Double_t z,t,w; Double_t diru,dirv,dirz,dirs; Double_t direrru,direrrv; Double_t offset; Int_t flag,npts; npts=0; sw=0.0; swz=0.0; swzz=0.0; swt=0.0; swtt=0.0; swzt=0.0; for(i=0;i<npln;i++){ if( plnvuw[i]>-1 ){ z=Z[i]; t=CT[i]; w=Q[i]; swz+=w*z; swzz+=w*z*z; swt+=w*t; swtt+=w*t*t; swzt+=w*z*t; sw+=w; npts++; } } if(npts>1){ ctimeslope=(sw*swzt-swz*swt)/(sw*swzz-swz*swz); ctimeoffset=(swt*swzz-swz*swzt)/(sw*swzz-swz*swz); ctimeaverage=(swt/sw); if( (swtt/sw)-(swt/sw)*(swt/sw)>0.0 && (swzz/sw)-(swz/sw)*(swz/sw)>0.0 ){ ctimescatter=((swzt/sw)-(swz/sw)*(swt/sw))/((sqrt((swzz/sw)-(swz/sw)*(swz/sw)))*(sqrt((swtt/sw)-(swt/sw)*(swt/sw)))); } } offset=ctimeaverage; shwscr=ctimescatter; if(vtxshw) dirscr=trkscr; else dirscr=shwscr; npts=0; diru=0.0; direrru=0.0; if(1+shwu->GetEndPlane()-shwu->GetBegPlane()>1){ sw=0.0; swz=0.0; swt=0.0; swzz=0.0; swzt=0.0; swtt=0.0; for(i=0;i<npln;i++){ flag=0; for(j=0;j<1+fHitArr[i].GetLast();j++){ HitAtNu* hit = (HitAtNu*)(fHitArr[i].At(j)); if(hit->GetShwFlag()>1 && hit->GetPlaneView()==0){ du=hit->GetTPos()-vtxu; dz=hit->GetZPos()-vtxz; w=hit->GetCharge(); swz+=w*dz; swt+=w*du; swzz+=w*dz*dz; swzt+=w*dz*du; swtt+=w*du*du; sw+=w; flag=1; } } if(flag) npts++; } if(npts>1){ diru=(sw*swzt-swz*swt)/(sw*swzz-swz*swz); if( (swtt/sw)-(swt/sw)*(swt/sw)>0.0 && (swzz/sw)-(swz/sw)*(swz/sw)>0.0 ){ direrru=((swzt/sw)-(swz/sw)*(swt/sw))/((sqrt((swzz/sw)-(swz/sw)*(swz/sw)))*(sqrt((swtt/sw)-(swt/sw)*(swt/sw)))); } } } npts=0; dirv=0.0; direrrv=0.0; if(1+shwv->GetEndPlane()-shwv->GetBegPlane()>1){ sw=0.0; swz=0.0; swt=0.0; swzz=0.0; swzt=0.0; swtt=0.0; for(i=0;i<npln;i++){ flag=0; for(j=0;j<1+fHitArr[i].GetLast();j++){ HitAtNu* hit = (HitAtNu*)(fHitArr[i].At(j)); if(hit->GetShwFlag()>1 && hit->GetPlaneView()==1){ dv=hit->GetTPos()-vtxv; dz=hit->GetZPos()-vtxz; w=hit->GetCharge(); swz+=w*dz; swt+=w*dv; swzz+=w*dz*dz; swzt+=w*dz*dv; swtt+=w*dv*dv; sw+=w; flag=1; } } if(flag) npts++; } if(npts>1){ dirv=(sw*swzt-swz*swt)/(sw*swzz-swz*swz); if( (swtt/sw)-(swt/sw)*(swt/sw)>0.0 && (swzz/sw)-(swz/sw)*(swz/sw)>0.0 ){ direrrv=((swzt/sw)-(swz/sw)*(swt/sw))/((sqrt((swzz/sw)-(swz/sw)*(swz/sw)))*(sqrt((swtt/sw)-(swt/sw)*(swt/sw)))); } } } dirs=sqrt(diru*diru+dirv*dirv+1.0); diru=diru/dirs; dirv=dirv/dirs; dirz=1.0/dirs; MSG("AlgShowerAtNu", Msg::kVerbose) << " ... direction: " << " (" << diru << "," << dirv << "," << dirz << ")" << endl; // *********************************** // * C A L C U L A T E E N E R G Y * // *********************************** MSG("AlgShowerAtNu", Msg::kDebug) << " *** calculate energy *** " << endl; // do the SIGMAP/SIGMIP calibration PlexStripEndId plexseid; FloatErr sigcorr,sigmip,sigmap; FloatErr lpos; for(pln=0;pln<npln;pln++){ for(i=0;i<1+fHitArr[pln].GetLast();i++){ HitAtNu* hit = (HitAtNu*)(fHitArr[pln].At(i)); CandStripHandle* strip = (CandStripHandle*)(hit->GetCandStripHandle()); lpos = 0.0; if( strip->GetPlaneView()==PlaneView::kU || strip->GetPlaneView()==PlaneView::kX ) lpos = V[pln]; if( strip->GetPlaneView()==PlaneView::kV || strip->GetPlaneView()==PlaneView::kY ) lpos = U[pln]; lpos.SetError(1.0); if(strip->GetNDigit(StripEnd::kPositive)>0){ plexseid=strip->GetStripEndId(StripEnd::kPositive); sigcorr=strip->GetCharge(StripEnd::kPositive,CalDigitType::kSigCorr); sigmap=cal.GetAttenCorrectedTpos(sigcorr,lpos,plexseid); sigmip=cal.GetMIP(sigmap,plexseid); shower.CalibrateSigMapped(plexseid.GetEncoded(),sigmap); shower.CalibrateMIP(plexseid.GetEncoded(),sigmip); } if(strip->GetNDigit(StripEnd::kNegative)>0){ plexseid=strip->GetStripEndId(StripEnd::kNegative); sigcorr=strip->GetCharge(StripEnd::kNegative,CalDigitType::kSigCorr); sigmap=cal.GetAttenCorrectedTpos(sigcorr,lpos,plexseid); sigmip=cal.GetMIP(sigmap,plexseid); shower.CalibrateSigMapped(plexseid.GetEncoded(),sigmap); shower.CalibrateMIP(plexseid.GetEncoded(),sigmip); } } } // determine overall direction Double_t dir=1.0; if(dirscr>0.0) dir=1.0; if(dirscr<0.0) dir=-1.0; // calculate energy from pulse height Double_t energy = 0.; if(totph>0.0){ energy = totph; // NEEDS TUNING energy = 0.3 + totph/150.0; } // ******************************************* // * S E T S H O W E R V A R I A B L E S * // ******************************************* // set reseed flag Bool_t reseedflag=0; if(shwu->GetReseedFlag()||shwv->GetReseedFlag()) reseedflag=1; // set shower variables shower.SetMinPlane(bpln); shower.SetMaxPlane(epln); shower.SetNPlanes(nplanes); shower.SetNStrips(nstrips); shower.SetTimeSlope((dir)/(3.0e8)); shower.SetTimeOffset((offset)/(3.0e8)); shower.SetDirTimeScore(shwscr); shower.SetDirVtxShwScore(dirscr); shower.SetVtxShw(vtxshw); shower.SetVtxU(vtxu); shower.SetVtxV(vtxv); shower.SetVtxZ(vtxz); shower.SetVtxR(vtxr); shower.SetVtxT((offset)/(3.0e8)); shower.SetVtxPlane(vtxpln); shower.SetDirCosU(dir*diru); shower.SetDirCosV(dir*dirv); shower.SetDirCosZ(dir*dirz); shower.SetDirCosErrU(dir*direrru); shower.SetDirCosErrV(dir*direrrv); shower.SetEnergy(energy); shower.SetReseedFlag(reseedflag); MSG("AlgShowerAtNu", Msg::kDebug) << " *** CREATING NEW SHOWER *** " << endl << " MinPlane = " << shower.GetMinPlane() << endl << " MaxPlane = " << shower.GetMaxPlane() << endl << " NPlanes = " << shower.GetNPlanes() << endl << " NStrips = " << shower.GetNStrips() << endl << " TimeSlope = " << shower.GetTimeSlope() << endl << " TimeOffset = " << shower.GetTimeOffset() << endl << " DirTimeScore = " << shower.GetDirTimeScore() << endl << " DirVtxShwScore= " << shower.GetDirVtxShwScore() << endl << " VtxShw = " << shower.GetVtxShw() << endl << " VtxU = " << shower.GetVtxU() << endl << " VtxV = " << shower.GetVtxV() << endl << " VtxZ = " << shower.GetVtxZ() << endl << " VtxR = " << shower.GetVtxR() << endl << " VtxT = " << shower.GetVtxT() << endl << " VtxPln = " << shower.GetVtxPlane() << endl << " DirCosU = " << shower.GetDirCosU() << endl << " DirCosV = " << shower.GetDirCosV() << endl << " DirCosZ = " << shower.GetDirCosZ() << endl << " DirCosErrU = " << shower.GetDirCosErrU() << endl << " DirCosErrV = " << shower.GetDirCosErrV() << endl << " Energy = " << shower.GetEnergy() << endl << " ReseedFlag = " << shower.GetReseedFlag() << endl; for(i=0;i<npln;i++){ fHitArr[i].Clear(); } delete [] Z; delete [] U; delete [] V; delete [] CT; delete [] Q; delete [] Qm; delete [] Qp; delete [] CTm; delete [] CTp; delete [] plnvuw; delete [] plnnum; return; }

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

Reimplemented from AlgBase. Definition at line 617 of file AlgShowerAtNu.cxx. { }

---

Member Data Documentation

TObjArray AlgShowerAtNu::fHitArr[500] [private]

Definition at line 18 of file AlgShowerAtNu.h. Referenced by RunAlg().

---

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

• AlgShowerAtNu.h
• AlgShowerAtNu.cxx

---