#include <CondensedNtpModuleAtm.h>
Inheritance diagram for CondensedNtpModuleAtm:

Public Member Functions | |
| CondensedNtpModuleAtm () | |
| virtual | ~CondensedNtpModuleAtm () |
| void | BeginJob () |
| JobCResult | Ana (const MomNavigator *mom) |
| const Registry & | DefaultConfig () const |
| void | Config (const Registry &r) |
| void | Help () |
| void | EndJob () |
Private Member Functions | |
| void | FillEventInformation (ANtpRecoNtpManipulator *ntpManipulator, NtpSREvent *ntpEvent) |
| void | FillMCInformation (NtpMCTruth *ntpMCTruth, NtpMCStdHep *ntpMCStdHep) |
| void | FillTrackInformation (NtpSRTrack *ntpTrack, TClonesArray *stripArray, double azimuth, double ra, double dec) |
| void | FillTrackTimingVariables (NtpSRTrack *ntpTrack, TClonesArray *stripArray) |
| void | FillTrackMagneticBendingVariables (NtpSRTrack *ntpTrack, TClonesArray *stripArray) |
| void | ResetTreeVariables () |
| Double_t | GetTimeWeight (Float_t ph) |
| Double_t | ArrivalTime_Weight (Double_t npe) const |
| Double_t | ArrivalTime_Uncertainty (Double_t npe) const |
| void | WeightedAverage (const Int_t n, const Double_t *x, const Double_t *w, Double_t *parm, Double_t *eparm) |
| void | LinearFit_Weighted (const Int_t n, const Double_t *x, const Double_t *y, const Double_t *w, Double_t *parm, Double_t *eparm) |
| void | FitMinimizingResidual (Int_t &n, Double_t *x, Double_t *y, Double_t *w, Double_t *param, Double_t *eparam, bool findSlope=true) |
| void | FindLinearFit (Double_t *x, Double_t *y, Double_t *weight, Int_t nPoints, Double_t &a1, Double_t &a2, Double_t &chiSq) |
| void | FindChi2 (Double_t *x, Double_t *y, Double_t *weight, Int_t nPoints, Double_t &a1, Double_t &a2, Double_t &chiSq) |
| void | FindStraightLineDeviation (Double_t *x, Double_t *y, Int_t nPoints, Double_t &a1, Double_t &a2, Double_t &deviation) |
| void | FindMeanAndRMS (double *x, double *weight, int nPoints, double &mean, double &rms) |
| Float_t | CheapBaseline (float cosz) |
Private Attributes | |
| std::string | fFileName |
| std::string | fRockMapFileName |
| std::string | fTreeName |
| std::string | fBeamTreeName |
| std::string | fBeamPath |
| TFile * | fNtpFile |
| TTree * | fNtuple |
| TTree * | fFailTree |
| TTree * | fBeamTree |
| Int_t | fDetector |
| Int_t | fGoodTrack |
| Int_t | fDataType |
| Int_t | fEndPlane |
| Double_t | fHornCurrent |
| Double_t | fCurrentPOT |
| Double_t | fBeamHPos |
| Double_t | fBeamVPos |
| Double_t | fTargetPos |
| Double_t | fBeamTimeSec |
| Int_t | fIsMultiple |
| Int_t | fValidPlaneFail |
| Int_t | fFailDeMux |
| Int_t | fMajorRelease |
| Double_t | fPathCosZenDeg [1500] |
| Double_t | fPathAzimuthDeg [1500] |
| Double_t | fColDenCosZen [1500] |
| Double_t | fDensity [1500] |
| bool | fNewMC |
| ANtpEventInfo * | fEventInfo |
| ANtpHeaderInfo * | fHeaderInfo |
| ANtpBeamInfo * | fBeamInfo |
| ANtpTrackInfoAtm * | fTrackInfo |
| ANtpTruthInfoAtm * | fTruthInfo |
| ANtpInfoObjectFillerBeam * | fInfoFiller |
| Int_t | fCtr |
Definition at line 41 of file CondensedNtpModuleAtm.h.
| CondensedNtpModuleAtm::CondensedNtpModuleAtm | ( | ) |
Definition at line 70 of file CondensedNtpModuleAtm.cxx.
References Msg::kDebug, and MSG.
00070 : 00071 fFileName("analysisNtuple.root"), 00072 fTreeName("analysisNtuple"), 00073 fNtpFile(0), 00074 fNtuple(0), 00075 fBeamTree(0), 00076 fDetector(1), 00077 fEndPlane(485), 00078 fEventInfo(0), 00079 fHeaderInfo(0), 00080 fTrackInfo(0), 00081 fTruthInfo(0), 00082 fCtr(0) 00083 { 00084 MSG("JobC", Msg::kDebug) << "CondensedNtpModuleAtm::Constructor" << endl; 00085 00086 fHeaderInfo = new ANtpHeaderInfo(); 00087 fBeamInfo = new ANtpBeamInfo(); 00088 fEventInfo = new ANtpEventInfo(); 00089 fTrackInfo = new ANtpTrackInfoAtm(); 00090 fTruthInfo = new ANtpTruthInfoAtm(); 00091 00092 fInfoFiller = new ANtpInfoObjectFillerBeam(); 00093 00094 }
| CondensedNtpModuleAtm::~CondensedNtpModuleAtm | ( | ) | [virtual] |
Definition at line 97 of file CondensedNtpModuleAtm.cxx.
References fBeamInfo, fEventInfo, fHeaderInfo, fTrackInfo, fTruthInfo, Msg::kDebug, and MSG.
00098 { 00099 00100 MSG("JobC", Msg::kDebug) << "CondensedNtpModuleAtm::Destructor" << endl; 00101 00102 if(fHeaderInfo) delete fHeaderInfo; 00103 if(fBeamInfo) delete fBeamInfo; 00104 if(fEventInfo) delete fEventInfo; 00105 if(fTrackInfo) delete fTrackInfo; 00106 if(fTruthInfo) delete fTruthInfo; 00107 00108 }
| JobCResult CondensedNtpModuleAtm::Ana | ( | const MomNavigator * | mom | ) | [virtual] |
Implement this for read only access to the MomNavigator
Reimplemented from JobCModule.
Definition at line 190 of file CondensedNtpModuleAtm.cxx.
References det, ANtpHeaderInfo::detector, fBeamInfo, fCtr, fDataType, fDetector, fEndPlane, fEventInfo, fFailDeMux, fFailTree, fGoodTrack, fHeaderInfo, ANtpInfoObjectFillerBeam::FillBeamInformation(), FillEventInformation(), ANtpInfoObjectFiller::FillHeaderInformation(), FillMCInformation(), FillTrackInformation(), FillTrackMagneticBendingVariables(), FillTrackTimingVariables(), fInfoFiller, fIsMultiple, fMajorRelease, fNewMC, fNtuple, fTrackInfo, fValidPlaneFail, ANtpRecoNtpManipulator::GetDmxStatusInfo(), ANtpEventManipulator::GetEvent(), ANtpRecoNtpManipulator::GetEventManipulator(), MomNavigator::GetFragment(), ANtpRecoNtpManipulator::GetSnarlEventSummary(), ANtpRecoNtpManipulator::GetStripArray(), VldContext::GetTimeStamp(), RecRecordImp< T >::GetVldContext(), ANtpEventInfo::index, NtpSRDmxStatus::ismultimuon, SimFlag::kData, Msg::kDebug, Detector::kFar, SimFlag::kMC, Detector::kNear, JobCResult::kPassed, Msg::kWarning, MSG, NtpSREventSummary::nevent, NtpSRDmxStatus::nonphysicalfail, NtpSREventSummary::ntrack, ANtpTrackInfo::planes, ResetTreeVariables(), ANtpHeaderInfo::run, ANtpEventManipulator::SetEventInSnarl(), JobCResult::SetFailed(), ANtpRecoNtpManipulator::SetPrimaryShowerCriteria(), ANtpRecoNtpManipulator::SetPrimaryTrackCriteria(), ANtpInfoObjectFiller::SetStripArray(), JobCResult::SetWarning(), ANtpHeaderInfo::snarl, ANtpTrackInfo::totalStrips, ANtpEventInfo::totalStrips, NtpSRDmxStatus::ustrayplanes, NtpSRDmxStatus::uvalidplanes, NtpSRDmxStatus::validplanesfail, NtpSRDmxStatus::vertexplanefail, NtpSRDmxStatus::vstrayplanes, NtpSRDmxStatus::vvalidplanes, and ANtpHeaderInfo::year.
00191 { 00192 MSG("CondensedNtpModuleAtm", Msg::kDebug) << "start ana" << endl; 00193 00194 JobCResult result(JobCResult::kPassed); 00195 00196 //get the records from MOM 00197 assert(mom); 00198 00199 NtpSRRecord *record = dynamic_cast<NtpSRRecord *>(mom->GetFragment("NtpSRRecord")); 00200 NtpStRecord *stRecord = dynamic_cast<NtpStRecord *>(mom->GetFragment("NtpStRecord")); 00201 00202 if(!record && !stRecord){ 00203 MSG("CondensedNtpModuleAtm", Msg::kWarning) << "Could not get NtpSR or NtpSt Record from MOM" << endl; 00204 result.SetWarning().SetFailed(); 00205 return result; 00206 }//end if no NtpSRRecord 00207 00208 NtpMCRecord *mcRecord = 0; 00209 NtpTHRecord *thRecord = 0; 00210 00211 mcRecord = dynamic_cast<NtpMCRecord *>(mom->GetFragment("NtpMCRecord")); 00212 thRecord = dynamic_cast<NtpTHRecord *>(mom->GetFragment("NtpTHRecord")); 00213 00214 //get the record for the beam tree 00215 NtpBDLiteRecord *bdRecord = dynamic_cast<NtpBDLiteRecord *>(mom->GetFragment("NtpBDLiteRecord")); 00216 if(!bdRecord) MSG("CondensedNtpModuleAtm", Msg::kDebug) << "NO BEAM RECORD" << endl; 00217 00218 //fill the header information 00219 Detector::Detector_t det = Detector::kFar; 00220 if(fDetector<1) det = Detector::kNear; 00221 //get a sim flag 00222 SimFlag::SimFlag_t dataType = SimFlag::kData; 00223 if(fDataType) dataType = SimFlag::kMC; 00224 00225 //instantiate a NtpHelper object to help you get the info you want 00226 ANtpRecoNtpManipulator *ntpManipulator = 0; 00227 if(record) 00228 ntpManipulator = new ANtpRecoNtpManipulator(record, mcRecord, thRecord); 00229 else if(stRecord) 00230 ntpManipulator = new ANtpRecoNtpManipulator(stRecord); 00231 00232 VldTimeStamp timeStamp = stRecord ? stRecord->GetVldContext()->GetTimeStamp() : record->GetVldContext()->GetTimeStamp(); 00233 00234 MSG("CondensedNtpModuleAtm", Msg::kDebug) << timeStamp << endl; 00235 00236 // TStopwatch stopwatch; 00237 // stopwatch.Start(); 00238 00239 fInfoFiller->SetStripArray(ntpManipulator->GetStripArray()); 00240 fInfoFiller->FillHeaderInformation(ntpManipulator, det, dataType, fHeaderInfo); 00241 00242 // stopwatch.Stop(); 00243 // MSG("CondensedNtpModuleAtm", Msg::kDebug) << "took " << stopwatch.RealTime() 00244 // << " to fill header" << endl; 00245 // stopwatch.Reset(); 00246 00247 // stopwatch.Start(); 00248 00249 if(fHeaderInfo->year > 2004){ 00250 if(bdRecord) fInfoFiller->FillBeamInformation(bdRecord, fBeamInfo); 00251 else fInfoFiller->FillBeamInformation(timeStamp, det, dataType, fBeamInfo); 00252 } 00253 00254 // stopwatch.Stop(); 00255 // MSG("CondensedNtpModuleAtm", Msg::kDebug) << "took " << stopwatch.RealTime() 00256 // << " to fill beam" << endl; 00257 // stopwatch.Reset(); 00258 00259 //set up which flags you want to use to determine the primary shower or track 00260 //a value of 0 for a flag means it will not be used 00261 ntpManipulator->SetPrimaryTrackCriteria(0,1,0); // nplanes, length, total pulse height 00262 ntpManipulator->SetPrimaryShowerCriteria(0,1); // nplanes, total pulse height 00263 00264 //check that the current event passed demuxing 00265 bool passedDeMux = true; 00266 fGoodTrack = true; 00267 fIsMultiple = 0; 00268 fValidPlaneFail = 0; 00269 fFailDeMux = 0; 00270 00271 MSG("CondensedNtpModuleAtm", Msg::kDebug) << "run,snarl = " << fHeaderInfo->run 00272 << "," << fHeaderInfo->snarl << " " 00273 << ntpManipulator->GetSnarlEventSummary().ntrack 00274 << " " << passedDeMux << endl; 00275 00276 if(fHeaderInfo->detector == Detector::kFar){ 00277 00278 if(ntpManipulator->GetDmxStatusInfo().ismultimuon){ 00279 passedDeMux = false; 00280 fIsMultiple = 1; 00281 } 00282 if(ntpManipulator->GetDmxStatusInfo().validplanesfail){ 00283 passedDeMux = false; 00284 fValidPlaneFail = 1; 00285 } 00286 if(fMajorRelease < 18){ 00287 00288 if(ntpManipulator->GetDmxStatusInfo().vertexplanefail 00289 || ntpManipulator->GetDmxStatusInfo().nonphysicalfail) passedDeMux = false; 00290 00291 UShort_t uValidPlanes = ntpManipulator->GetDmxStatusInfo().uvalidplanes; 00292 UShort_t uStrayPlanes = ntpManipulator->GetDmxStatusInfo().ustrayplanes; 00293 UShort_t vValidPlanes = ntpManipulator->GetDmxStatusInfo().vvalidplanes; 00294 UShort_t vStrayPlanes = ntpManipulator->GetDmxStatusInfo().vstrayplanes; 00295 Float_t ufrac = 0., vfrac = 0.; 00296 00297 if(uValidPlanes>0) 00298 ufrac = (1.*uStrayPlanes)/(1.*uValidPlanes); 00299 else ufrac = 1.; 00300 if(vValidPlanes>0) 00301 vfrac = (1.*vStrayPlanes)/(1.*vValidPlanes); 00302 else vfrac = 1.; 00303 if(ufrac>=0.45) passedDeMux = false; 00304 else if( vfrac>=0.45 ) passedDeMux = false; 00305 else if(uStrayPlanes>=3&&uValidPlanes<20) passedDeMux = false; 00306 else if(vStrayPlanes>=3&&vValidPlanes<20) passedDeMux = false; 00307 else if(uStrayPlanes>=4&&uValidPlanes>=20 00308 && uValidPlanes<=40) passedDeMux = false; 00309 else if(vStrayPlanes>=4&&vValidPlanes>=20 00310 && vValidPlanes<=40) passedDeMux = false; 00311 else if(ufrac>=0.1 && uValidPlanes>=40) passedDeMux = false; 00312 else if(vfrac>=0.1 && vValidPlanes>=40) passedDeMux = false; 00313 00314 }//end if before switch over to cambridge demuxing 00315 }//end if far detector 00316 00317 //if the event isnt demuxed well, then forget it also forget it 00318 //if there arent any tracks in the event 00319 if(!passedDeMux || ntpManipulator->GetSnarlEventSummary().ntrack==0 || fEndPlane<248){ 00320 fFailTree->Fill(); 00321 return result; 00322 } 00323 00324 //set up some pointers to typical NtpSR objects 00325 NtpSRTrack *ntpTrack = 0; 00326 NtpSREvent *ntpEvent = 0; 00327 NtpMCTruth *ntpMCTruth = 0; 00328 NtpMCStdHep *ntpMCStdHep = 0; 00329 NtpTHEvent *ntpTHEvent = 0; 00330 00331 //loop over the events in the snarl 00332 for(Int_t i = 0; i < ntpManipulator->GetSnarlEventSummary().nevent; ++i){ 00333 00334 MSG("CondensedNtpModuleAtm", Msg::kDebug) << "on event " << i << endl; 00335 00336 ResetTreeVariables(); 00337 00338 fEventInfo->index = i; 00339 00340 // stopwatch.Start(); 00341 //get event i. the call to ANtpRecoNtpManipulator::GetEventInSnarl sets the 00342 //NtpSREvent data member in that object to the current event so that 00343 //you can later call for the primary track, shower, etc. 00344 ntpManipulator->GetEventManipulator()->SetEventInSnarl(i); 00345 ntpEvent = ntpManipulator->GetEventManipulator()->GetEvent(); 00346 if(ntpEvent) FillEventInformation(ntpManipulator,ntpEvent); 00347 // stopwatch.Stop(); 00348 // MSG("CondensedNtpModuleAtm", Msg::kInfo) << "took " << stopwatch.RealTime() 00349 // << " to fill event" << endl; 00350 // stopwatch.Reset(); 00351 00352 // stopwatch.Start(); 00353 //get the primary track for the event - if no track is present it 00354 //returns 0 00355 fGoodTrack = false; 00356 ntpTrack = ntpManipulator->GetEventManipulator()->GetPrimaryTrack(); 00357 if(ntpTrack){ 00358 fGoodTrack = true; 00359 FillTrackInformation(ntpTrack, ntpManipulator->GetStripArray(), 00360 ntpManipulator->GetCosmicRayInfo().azimuth, 00361 ntpManipulator->GetCosmicRayInfo().ra, 00362 ntpManipulator->GetCosmicRayInfo().dec); 00363 FillTrackMagneticBendingVariables(ntpTrack, ntpManipulator->GetStripArray()); 00364 if(fDetector>0&&fGoodTrack) FillTrackTimingVariables(ntpTrack, ntpManipulator->GetStripArray()); 00365 } 00366 // stopwatch.Stop(); 00367 // MSG("CondensedNtpModuleAtm", Msg::kInfo) << "took " << stopwatch.RealTime() 00368 // << " to fill track" << endl; 00369 // stopwatch.Reset(); 00370 00371 // stopwatch.Start(); 00372 // Get best neu match from Truthhelper 00373 ntpTHEvent = ntpManipulator->GetMCManipulator()->GetNtpTHEvent(i); 00374 MSG("CondensedNtpModuleAtm", Msg::kDebug) << "start filling mc info " << ntpTHEvent->neumc << endl; 00375 if(ntpTHEvent) ntpMCTruth = ntpManipulator->GetMCManipulator()->GetNtpMCTruth(ntpTHEvent->neumc); 00376 if(ntpMCTruth){ 00377 //get the NtpMCStdHep object - use the fact that i put the values i want in the parent 00378 //place, which means i want the 0th entry 00379 if(fNewMC) ntpMCStdHep = ntpManipulator->GetMCManipulator()->GetNtpMCStdHep(0); 00380 FillMCInformation(ntpMCTruth, ntpMCStdHep); 00381 } 00382 // stopwatch.Stop(); 00383 // MSG("CondensedNtpModuleAtm", Msg::kInfo) << "took " << stopwatch.RealTime() 00384 // << " to fill truth" << endl; 00385 // stopwatch.Reset(); 00386 00387 if(fEventInfo->totalStrips<fTrackInfo->totalStrips) fGoodTrack = false; 00388 if(fTrackInfo->planes<1) fGoodTrack = false; 00389 if(fTrackInfo->totalStrips<1) fGoodTrack = false; 00390 00391 MSG("CondensedNtpModuleAtm", Msg::kDebug) << "strips = " << fEventInfo->totalStrips 00392 << "/" << fTrackInfo->totalStrips << " planes = " 00393 << fTrackInfo->planes << " good track = " 00394 << (int)fGoodTrack << endl; 00395 00396 //got all the necessary sr info, now fill the ntuple 00397 if(!fGoodTrack) continue; 00398 00399 fNtuple->Fill(); 00400 }//end loop over events for this snarl 00401 00402 ++fCtr; 00403 00404 delete ntpManipulator; 00405 00406 return result; 00407 }
| Double_t CondensedNtpModuleAtm::ArrivalTime_Uncertainty | ( | Double_t | npe | ) | const [private] |
Definition at line 1025 of file CondensedNtpModuleAtm.cxx.
Referenced by ArrivalTime_Weight().
01026 { 01027 // cout << "in ArrivalTime_Uncertainty" << endl; 01028 01029 Double_t sigma = 0.; // sigma is uncertainty in arrival of first p.e. 01030 01031 // a simple monte carlo was done to calculate the uncertainty in the 01032 // arrival time of the first p.e. as a function of # of p.e. 01033 // 01034 // the only assumption made was about the shape of the pulse, a 2 ns risetime 01035 // and an 8 ns decay time of the fluor 01036 // 01037 // rlee feb 2002 01038 // 01039 01040 if(npe<10.) sigma = 0.90+8.0/npe; 01041 else sigma = 5.0*exp(-0.5*log(npe)); 01042 01043 // cout << sigma << endl; 01044 01045 return sigma; 01046 }
| Double_t CondensedNtpModuleAtm::ArrivalTime_Weight | ( | Double_t | npe | ) | const [private] |
Definition at line 1013 of file CondensedNtpModuleAtm.cxx.
References ArrivalTime_Uncertainty().
Referenced by GetTimeWeight().
01014 { 01015 // cout << "in ArrivalTime_Weight" << endl; 01016 01017 Double_t weight = 1./ArrivalTime_Uncertainty(npe); 01018 01019 // //cout << weight << endl; 01020 01021 return weight*weight; 01022 }
| void CondensedNtpModuleAtm::BeginJob | ( | ) | [virtual] |
Implement for notification of begin of job
Reimplemented from JobCModule.
Definition at line 111 of file CondensedNtpModuleAtm.cxx.
References fBeamInfo, fBeamPath, fBeamTreeName, fColDenCosZen, fDataType, fDensity, fEventInfo, fFailDeMux, fFailTree, fFileName, fHeaderInfo, fIsMultiple, fNtpFile, fNtuple, fPathAzimuthDeg, fPathCosZenDeg, fRockMapFileName, fTrackInfo, fTreeName, fTruthInfo, fValidPlaneFail, Msg::kDebug, Msg::kInfo, and MSG.
00112 { 00113 // save the current working directory 00114 TDirectory* savedir = gDirectory; 00115 00116 00117 //create the new TFile for holding the electronics tree 00118 // this changes the value of gDirectory 00119 fNtpFile = new TFile(fFileName.c_str(),"RECREATE"); 00120 00121 MSG("CondensedNtpModuleAtm", Msg::kDebug) << fFileName << " " << fTreeName 00122 << " " << fNtpFile->GetName() << endl; 00123 00124 00125 // create TTree, these will attach themselves to the current 00126 //working directory 00127 fNtuple = new TTree(fTreeName.c_str(), "Analysis Tree"); 00128 00129 fNtuple->Branch("header.", "ANtpHeaderInfo", &fHeaderInfo, 64000, 2); 00130 fNtuple->Branch("beam.", "ANtpBeamInfo", &fBeamInfo, 64000, 2); 00131 fNtuple->Branch("event.", "ANtpEventInfo", &fEventInfo, 64000, 2); 00132 fNtuple->Branch("track.", "ANtpTrackInfoAtm", &fTrackInfo, 64000, 2); 00133 if(fDataType==1) fNtuple->Branch("truth.", "ANtpTruthInfoAtm", &fTruthInfo, 64000, 2); 00134 00135 MSG("CondensedNtpModuleAtm", Msg::kDebug) << "got branches" << endl; 00136 00137 fFailTree = new TTree("failures", "Failed Events"); 00138 fFailTree->Branch("header.", "ANtpHeaderInfo", &fHeaderInfo, 64000, 2); 00139 fFailTree->Branch("multiple", &fIsMultiple, "multiple/I"); 00140 fFailTree->Branch("validplane", &fValidPlaneFail, "validplane/I"); 00141 fFailTree->Branch("failDeMux", &fFailDeMux, "failDeMux/I"); 00142 00143 // change back to current working directory before leaving constructor 00144 savedir->cd(); 00145 00146 //load the vertical intensity information into some arrays 00147 //get the slant depths and the number of events for each slant depth 00148 //read in rock data 00149 // http://www.hep.umn.ed/~schubert/far/s2rock/slantdepth.data 00150 // 1-cos(zen); azi(degrees); slantdepth (gmcm**2) * cos(zen); 00151 // density (gm/cm**2) 00152 // N.B., units for density are incorrectly labelled in the file 00153 ifstream densityFile(fRockMapFileName.c_str()); 00154 if( !densityFile.is_open() ){ 00155 cout << "could not open rock map - exit" << endl; 00156 exit(0); 00157 } 00158 00159 int ir = 0; 00160 while( ir<1404 ){ 00161 densityFile >> fPathCosZenDeg[ir] >> fPathAzimuthDeg[ir] >> fColDenCosZen[ir] >> fDensity[ir]; 00162 ++ir; 00163 } 00164 00165 //get the beam info tree 00166 //check to see if you need to calculate the total Protons on Target 00167 MSG("BeamTestModule", Msg::kInfo) << "beam info path = " << fBeamPath.c_str() 00168 << " tree name = " << fBeamTreeName.c_str() << endl; 00169 00170 // if(fDataType < 1){ 00171 00172 // TFile *beamFile = new TFile(fBeamPath.c_str()); 00173 // fBeamTree = dynamic_cast<TTree *>(beamFile->Get(fBeamTreeName.c_str())); 00174 00175 // //build the index for the tree to select the utc timestamps 00176 // fBeamTree->BuildIndex("utc"); 00177 // fBeamTree->SetBranchAddress("intensity", &fCurrentPOT); 00178 // fBeamTree->SetBranchAddress("hornCurrent", &fHornCurrent); 00179 // fBeamTree->SetBranchAddress("hpos", &fBeamHPos); 00180 // fBeamTree->SetBranchAddress("vpos", &fBeamVPos); 00181 // fBeamTree->SetBranchAddress("tgtpos", &fTargetPos); 00182 // } 00183 00184 return; 00185 }
| Float_t CondensedNtpModuleAtm::CheapBaseline | ( | float | cosz | ) | [private] |
Definition at line 1399 of file CondensedNtpModuleAtm.cxx.
Referenced by FillMCInformation().
01400 { 01401 static float Rearth = 6378140.0; /* Equatorial radius of earth */ 01402 static float nuHeight = 1500.0; /* Average Neutrino production height */ 01403 static float altitudeSK = -210.0; /* Altitude of the detector - MINOS is 689ft below sea level*/ 01404 01405 static float x1s = (Rearth+nuHeight)*(Rearth+nuHeight); 01406 static float x2 = Rearth+altitudeSK; 01407 static float x2s = (Rearth+altitudeSK)*(Rearth+altitudeSK); 01408 01409 return sqrt(x1s+x2s*(cosz*cosz-1.0))-x2*cosz; 01410 }
| void CondensedNtpModuleAtm::Config | ( | const Registry & | r | ) | [virtual] |
Return the actual configuration. If your module directly pulls its configuration from the fConfig Registry, you don't need to override this. Override if you have local config variables.
Reimplemented from JobCModule.
Definition at line 474 of file CondensedNtpModuleAtm.cxx.
References fBeamPath, fBeamTreeName, fDataType, fDetector, fEndPlane, fFileName, fMajorRelease, fNewMC, fRockMapFileName, fTreeName, and Registry::Get().
00475 { 00476 int tmpb; // a temp bool. See comment under DefaultConfig... 00477 int tmpi; // a temp int. 00478 double tmpd; // a temp double. 00479 const char* tmps; // a temp string 00480 00481 if (r.Get("FileName", tmps)) fFileName = tmps; 00482 if (r.Get("RockMapFileName", tmps)) fRockMapFileName = tmps; 00483 if (r.Get("TreeName", tmps)) fTreeName = tmps; 00484 if (r.Get("BeamPath", tmps)) fBeamPath = tmps; 00485 if (r.Get("BeamTreeName", tmps)) fBeamTreeName = tmps; 00486 if (r.Get("EndPlane", tmpi)) fEndPlane = tmpi; 00487 if (r.Get("Detector", tmpi)) fDetector = tmpi; 00488 if (r.Get("DataType", tmpi)) fDataType = tmpi; 00489 if (r.Get("NewMC", tmpb)) fNewMC = tmpb; 00490 if (r.Get("MajorRelease", tmpi)) fMajorRelease = tmpi; 00491 00492 //quiet the compiler warnings 00493 tmpb = 0; 00494 tmpd = 1.; 00495 00496 return; 00497 }
| const Registry & CondensedNtpModuleAtm::DefaultConfig | ( | ) | const [virtual] |
Get the default configuration registry. This should normally be overridden. One useful idiom is to implement it like:
const Registry& MyModule::DefaultConfig() const { static Registry cfg; // never is destroyed if (cfg.Size()) return cfg; // already filled it // set defaults: cfg.Set("TheAnswer",42); cfg.Set("Units","unknown"); return cfg; }
Reimplemented from JobCModule.
Definition at line 446 of file CondensedNtpModuleAtm.cxx.
References Registry::LockValues(), Registry::Set(), and Registry::UnLockValues().
00447 { 00448 00449 int itrue = 1; // work around for lack of bool in registry 00450 int ifalse = 0; // work around for lack of bool in registry 00451 00452 static Registry r; 00453 00454 r.UnLockValues(); 00455 00456 r.Set("FileName", "analysisTree.root"); 00457 r.Set("RockMapFileName", "density.txt"); 00458 r.Set("TreeName", "analysisTree"); 00459 r.Set("BeamPath", "/afs/fnal.gov/files/data/minos/d04/brebel/beamInfo.root"); 00460 r.Set("BeamTreeName", "beam"); 00461 r.Set("EndPlane", 485); 00462 r.Set("Detector", 1); 00463 r.Set("DataType", 0); 00464 r.Set("NewMC", itrue); 00465 r.Set("MajorRelease", 14); 00466 r.LockValues(); 00467 00468 itrue = ifalse; //quiet the compiler warnings 00469 00470 return r; 00471 }
| void CondensedNtpModuleAtm::EndJob | ( | ) | [virtual] |
Implement for notification of end of job
Reimplemented from JobCModule.
Definition at line 420 of file CondensedNtpModuleAtm.cxx.
References fFailTree, fNtpFile, fNtuple, Msg::kDebug, Msg::kInfo, and MSG.
00421 { 00422 00423 MSG("CondensedNtpModuleAtm", Msg::kDebug) << "start end job method" << endl; 00424 00425 //save current working directory 00426 TDirectory *savedir = gDirectory; 00427 00428 //cd to the "directory" of the file 00429 fNtpFile->cd(); 00430 00431 MSG("CondensedNtpModuleAtm", Msg::kInfo) << fNtuple->GetEntries() << endl; 00432 00433 fNtuple->Write(); 00434 fFailTree->Write(); 00435 00436 //go back to original working directory 00437 savedir->cd(); 00438 00439 fNtpFile->Write(); 00440 fNtpFile->Close(); 00441 00442 return; 00443 }
| void CondensedNtpModuleAtm::FillEventInformation | ( | ANtpRecoNtpManipulator * | ntpManipulator, | |
| NtpSREvent * | ntpEvent | |||
| ) | [private] |
Definition at line 935 of file CondensedNtpModuleAtm.cxx.
References fEventInfo, ANtpInfoObjectFiller::FillEventInformation(), fInfoFiller, Msg::kDebug, and MSG.
Referenced by Ana().
00937 { 00938 MSG("CondensedNtpModuleAtm", Msg::kDebug) << "in FillEventInformation" << endl; 00939 00940 fInfoFiller->FillEventInformation(ntpManipulator, 00941 ntpEvent, 00942 fEventInfo); 00943 00944 return; 00945 }
| void CondensedNtpModuleAtm::FillMCInformation | ( | NtpMCTruth * | ntpMCTruth, | |
| NtpMCStdHep * | ntpMCStdHep | |||
| ) | [private] |
Definition at line 949 of file CondensedNtpModuleAtm.cxx.
References ANtpTruthInfoAtm::azimuth, ANtpTruthInfoAtm::baseLine, CheapBaseline(), ANtpInfoObjectFiller::FillMCTruthInformation(), fInfoFiller, fTruthInfo, Msg::kDebug, ANtpTruthInfo::leptonDCosX, ANtpTruthInfo::leptonDCosY, ANtpTruthInfo::leptonDCosZ, MSG, ANtpTruthInfo::nuEnergy, NtpMCTruth::p4neu, NtpMCStdHep::vtx, and ANtpTruthInfoAtm::weight.
Referenced by Ana().
00950 { 00951 MSG("CondensedNtpModuleAtm", Msg::kDebug) << "in FillMCInformation" << endl; 00952 00953 const double rotMatrix[] = { 0.894507011,0,0.447053919, 00954 0, 1, 0, 00955 -0.447053919,0,0.894507011}; 00956 00957 double dcosx_ideal, dcosy_ideal, dcosz_ideal; 00958 00959 fInfoFiller->FillMCTruthInformation(ntpMCTruth, fTruthInfo); 00960 00961 if(fTruthInfo->nuEnergy>0.) fTruthInfo->baseLine = 0.001*CheapBaseline(-ntpMCTruth->p4neu[1] 00962 /TMath::Abs(fTruthInfo->nuEnergy)); 00963 00964 //get the azimuth from AstroUtill 00965 dcosx_ideal = (rotMatrix[0]*fTruthInfo->leptonDCosX 00966 + rotMatrix[1]*fTruthInfo->leptonDCosY 00967 + rotMatrix[2]*fTruthInfo->leptonDCosZ); 00968 dcosy_ideal = (rotMatrix[3]*fTruthInfo->leptonDCosX 00969 + rotMatrix[4]*fTruthInfo->leptonDCosY 00970 + rotMatrix[5]*fTruthInfo->leptonDCosZ); 00971 dcosz_ideal = (rotMatrix[6]*fTruthInfo->leptonDCosX 00972 + rotMatrix[7]*fTruthInfo->leptonDCosY 00973 + rotMatrix[8]*fTruthInfo->leptonDCosZ); 00974 fTruthInfo->azimuth = 180.*TMath::ATan2(-dcosx_ideal,dcosz_ideal)/TMath::Pi(); 00975 if(fTruthInfo->azimuth < 0.) fTruthInfo->azimuth += 360.; 00976 00977 //use the NtpMCStdHep object to figure out the weight for this event based on the 00978 //fact that the new MC picks events at the surface. we plugged the value into 00979 //the vtx Y location in the stdhep common bloc, gminos multiplies it by 10^-3 00980 //because it thought it was converting mm to m 00981 if(ntpMCStdHep){ 00982 fTruthInfo->weight = 1.e3*ntpMCStdHep->vtx[1]; 00983 } 00984 00985 return; 00986 }
| void CondensedNtpModuleAtm::FillTrackInformation | ( | NtpSRTrack * | ntpTrack, | |
| TClonesArray * | stripArray, | |||
| double | azimuth, | |||
| double | ra, | |||
| double | dec | |||
| ) | [private] |
Definition at line 501 of file CondensedNtpModuleAtm.cxx.
References ANtpTrackInfoAtm::alternateChi2, ANtpTrackInfoAtm::azimuth, NtpSRTrackTime::cdtds, NtpSRVertex::dcosx, ANtpTrackInfo::dcosXEnd, ANtpTrackInfo::dcosXVtx, NtpSRVertex::dcosy, ANtpTrackInfo::dcosYEnd, ANtpTrackInfo::dcosYVtx, NtpSRVertex::dcosz, ANtpTrackInfo::dcosZEnd, ANtpTrackInfo::dcosZVtx, ANtpTrackInfoAtm::dec, NtpSRFiducial::dr, NtpSRFiducial::dz, NtpSRTrack::end, ANtpTrackInfo::endX, ANtpTrackInfo::endY, ANtpTrackInfo::endZ, fDensity, fDetector, fEndPlane, fEventInfo, fGoodTrack, NtpSRTrack::fidend, ANtpTrackInfoAtm::fiducialEndDr, ANtpTrackInfoAtm::fiducialEndDz, ANtpTrackInfoAtm::fiducialVtxDr, ANtpTrackInfoAtm::fiducialVtxDz, NtpSRTrack::fidvtx, ANtpInfoObjectFiller::FillTrackInformation(), fInfoFiller, NtpSRTrack::fit, fTrackInfo, ANtpTrackInfoAtm::impactParameter, ANtpTrackInfoAtm::invBeta, Msg::kDebug, PlaneView::kU, PlaneView::kV, ANtpTrackInfo::length, NtpSRTrack::momentum, MSG, NtpSRStrip::ndigit, NtpSRFitTrack::ndof, Munits::ns, NtpSRTrack::nstrip, NtpSRTrackPlane::ntrklike, NtpSRPlane::nu, NtpSRPlane::nv, NtpSRStrip::plane, NtpSRTrack::plane, ANtpEventInfo::planes, ANtpTrackInfoAtm::planesIn10Meter, ANtpTrackInfoAtm::planesIn15Meter, ANtpTrackInfoAtm::planesIn20Meter, ANtpTrackInfoAtm::planesIn25Meter, ANtpTrackInfoAtm::planesIn30Meter, ANtpTrackInfoAtm::planesIn35Meter, ANtpTrackInfoAtm::planesIn40Meter, ANtpTrackInfoAtm::planeUseFraction, NtpSRStrip::planeview, ANtpTrackInfo::pulseHeight, ANtpEventInfo::pulseHeight, NtpSRMomentum::qp, ANtpTrackInfoAtm::ra, ANtpTrackInfoAtm::signalUseFraction, ANtpTrackInfoAtm::slantDepth, NtpSRTrack::stp, NtpSRTrack::stpu, NtpSRTrack::stpv, NtpSRTrack::stpx, NtpSRTrack::stpy, NtpSRStrip::strip, NtpSRTrack::time, ANtpTrackInfo::totalStrips, NtpSRStrip::tpos, ANtpTrackInfoAtm::trackLikePlanes, ANtpTrackInfoAtm::twoEndStripFraction, ANtpTrackInfoAtm::uvAsymmetry, NtpSRTrack::vtx, ANtpTrackInfo::vtxX, ANtpTrackInfo::vtxY, ANtpTrackInfo::vtxZ, NtpSRVertex::x, NtpSRVertex::y, and NtpSRVertex::z.
Referenced by Ana().
00504 { 00505 MSG("CondensedNtpModuleAtm", Msg::kDebug) << "filling track info" << endl; 00506 00507 //fill the basics first 00508 fInfoFiller->FillTrackInformation(ntpTrack, fTrackInfo); 00509 00510 NtpSRStrip *ntpStrip = 0; 00511 00512 Int_t plane = 0, strip = 0, numDoubleEndedStrips = 0, prevPlane = -1; 00513 float xPos = 0., yPos = 0.; 00514 Float_t radius = 0.; 00515 00516 Int_t index = 0; 00517 00518 //stuff to find the slant depth 00519 double vec1[39],vec2[36]; 00520 int index1, index2, index3; 00521 00522 for(int i=0;i<39;++i){ 00523 vec1[i] = i*0.02; 00524 // cout << vec1[i] << endl; 00525 if(i<36){ 00526 vec2[i] = i*10. + 5.; 00527 // cout << vec2[i] << endl; 00528 } 00529 } 00530 00531 // Kashahara picks events according to where they come from 00532 // she uses altitude, not zenith angle 00533 index1 = TMath::BinarySearch(39,vec1,1.-(-fTrackInfo->dcosYVtx))+1; 00534 index2 = TMath::BinarySearch(36,vec2,azimuth)+1; 00535 index3 = (index1)*36 +(index2); 00536 fTrackInfo->slantDepth = fDensity[index3]*0.01; //density is in 10^5g/cm^2, ie km.w.e, 00537 //0.01 converts to m.w.e 10^3m/km/(10^5g/cm^2/km.w.e) 00538 00539 fTrackInfo->azimuth = azimuth; 00540 00541 fTrackInfo->ra = ra; 00542 fTrackInfo->dec = dec; 00543 fTrackInfo->trackLikePlanes = ntpTrack->plane.ntrklike; 00544 Int_t nVPlanes = ntpTrack->plane.nv; 00545 Int_t nUPlanes = ntpTrack->plane.nu; 00546 fTrackInfo->invBeta = ntpTrack->time.cdtds; 00547 fTrackInfo->planesIn10Meter = 0; 00548 fTrackInfo->planesIn15Meter = 0; 00549 fTrackInfo->planesIn20Meter = 0; 00550 fTrackInfo->planesIn25Meter = 0; 00551 fTrackInfo->planesIn30Meter = 0; 00552 fTrackInfo->planesIn35Meter = 0; 00553 fTrackInfo->planesIn40Meter = 0; 00554 00555 //check if the zenith angle indicates the right direction for the near detector. 00556 //if not, swap the ends and adjust the azimuth 00557 if(fTrackInfo->dcosYVtx > 0. && fDetector<1){ 00558 00559 fTrackInfo->vtxX = ntpTrack->end.x; 00560 fTrackInfo->vtxY = ntpTrack->end.y; 00561 fTrackInfo->vtxZ = ntpTrack->end.z; 00562 fTrackInfo->endX = ntpTrack->vtx.x; 00563 fTrackInfo->endY = ntpTrack->vtx.y; 00564 fTrackInfo->endZ = ntpTrack->vtx.z; 00565 fTrackInfo->dcosXVtx = -ntpTrack->end.dcosx; 00566 fTrackInfo->dcosYVtx = -ntpTrack->end.dcosy; 00567 fTrackInfo->dcosZVtx = -ntpTrack->end.dcosz; 00568 fTrackInfo->dcosXEnd = -ntpTrack->vtx.dcosx; 00569 fTrackInfo->dcosYEnd = -ntpTrack->vtx.dcosy; 00570 fTrackInfo->dcosZEnd = -ntpTrack->vtx.dcosz; 00571 00572 const double rotMatrix[] = { 0.894507011,0,0.447053919, 00573 0, 1, 0, 00574 -0.447053919,0,0.894507011}; 00575 00576 double dcosx_ideal, dcosy_ideal, dcosz_ideal; 00577 00578 //get the azimuth from AstroUtill 00579 dcosx_ideal = (rotMatrix[0]*fTrackInfo->dcosXVtx 00580 + rotMatrix[1]*fTrackInfo->dcosYVtx 00581 + rotMatrix[2]*fTrackInfo->dcosZVtx); 00582 dcosy_ideal = (rotMatrix[3]*fTrackInfo->dcosXVtx 00583 + rotMatrix[4]*fTrackInfo->dcosYVtx 00584 + rotMatrix[5]*fTrackInfo->dcosZVtx); 00585 dcosz_ideal = (rotMatrix[6]*fTrackInfo->dcosXVtx 00586 + rotMatrix[7]*fTrackInfo->dcosYVtx 00587 + rotMatrix[8]*fTrackInfo->dcosZVtx); 00588 00589 fTrackInfo->azimuth = 180.*TMath::ATan2(-dcosx_ideal,dcosz_ideal)/TMath::Pi(); 00590 if(fTrackInfo->azimuth < 0.) fTrackInfo->azimuth += 360.; 00591 } 00592 00593 if(ntpTrack->fit.ndof == 0 || ntpTrack->momentum.qp== 0.){ 00594 fGoodTrack = false; 00595 return; 00596 } 00597 00598 //the value of trk.fidvtx.dz is not valid for events with the full detector 00599 //so hack up my own value of fiducialDz 00600 Float_t fullDetEndZ = 0.0594*485. + 1.; 00601 00602 fTrackInfo->fiducialVtxDz = ntpTrack->fidvtx.dz; 00603 if(fEndPlane>248) 00604 fTrackInfo->fiducialVtxDz = TMath::Min(fullDetEndZ - fTrackInfo->vtxZ, fTrackInfo->vtxZ); 00605 fTrackInfo->fiducialVtxDr = ntpTrack->fidvtx.dr; 00606 00607 fTrackInfo->fiducialEndDz = ntpTrack->fidend.dz; 00608 if(fEndPlane>248) 00609 fTrackInfo->fiducialEndDz = TMath::Min(fullDetEndZ - fTrackInfo->endZ, fTrackInfo->endZ); 00610 fTrackInfo->fiducialEndDr = ntpTrack->fidend.dr; 00611 00612 fTrackInfo->uvAsymmetry = 1.; 00613 if(nVPlanes+nUPlanes>0) 00614 fTrackInfo->uvAsymmetry = TMath::Abs((1.*nUPlanes-1.*nVPlanes)/(1.*nVPlanes+1.*nUPlanes)); 00615 00616 fTrackInfo->signalUseFraction = 0.; 00617 if(fEventInfo->pulseHeight>0.) 00618 fTrackInfo->signalUseFraction = fTrackInfo->pulseHeight/(fEventInfo->pulseHeight); 00619 00620 fTrackInfo->planeUseFraction = 0.; 00621 if(fEventInfo->planes>0) 00622 fTrackInfo->planeUseFraction = 1.*fTrackInfo->trackLikePlanes/(1.*fEventInfo->planes); 00623 00624 Int_t numStrips = 0; 00625 Int_t numDigits = 0; 00626 00627 fTrackInfo->alternateChi2 = 0.; 00628 double deltaTPos = 0.; 00629 double stripPosError = 0.041/TMath::Sqrt(12.); 00630 00631 fTrackInfo->impactParameter = 20.0; 00632 fTrackInfo->totalStrips = 0; 00633 //does the pathlength of the track make sense? 00634 if(fTrackInfo->length<50.){ 00635 00636 00637 //loop over all strips in the track 00638 numStrips = (int)ntpTrack->nstrip; 00639 numDoubleEndedStrips = 0; 00640 00641 MSG("CondensedNtpModuleAtm", Msg::kDebug) << "start looping over track strips" << endl; 00642 00643 for(Int_t ns = 0; ns < numStrips; ++ns){ 00644 00645 //get the index for this strip 00646 index = ntpTrack->stp[ns]; 00647 00648 //get the NtpSRStrip object 00649 ntpStrip = dynamic_cast<NtpSRStrip *>(stripArray->At(index)); 00650 plane = (int)ntpStrip->plane; 00651 strip = (int)ntpStrip->strip; 00652 numDigits = (int)ntpStrip->ndigit; 00653 00654 xPos = ntpTrack->stpx[ns]; 00655 yPos = ntpTrack->stpy[ns]; 00656 00657 if((int)ntpStrip->planeview == PlaneView::kU) 00658 deltaTPos = ntpTrack->stpu[ns]-ntpStrip->tpos; 00659 else if((int)ntpStrip->planeview == PlaneView::kV) 00660 deltaTPos = ntpTrack->stpv[ns]-ntpStrip->tpos; 00661 00662 MSG("CondensedNtpModuleAtm", Msg::kDebug) << ntpTrack->stpu[ns] << " " << ntpStrip->tpos 00663 << " " << ntpTrack->stpv[ns] << " " << deltaTPos << endl; 00664 00665 fTrackInfo->alternateChi2 += deltaTPos*deltaTPos/(stripPosError*stripPosError); 00666 00667 if(numDigits >= 2){ 00668 00669 ++numDoubleEndedStrips; 00670 00671 radius = TMath::Sqrt(yPos*yPos + xPos*xPos); 00672 if(radius<fTrackInfo->impactParameter) 00673 fTrackInfo->impactParameter = TMath::Sqrt(yPos*yPos + xPos*xPos); 00674 00675 //get values to determine the number of plances crossed by the track 00676 //in concentric circles in the detector from 0.3 m out to 00677 //4.0 m. the points come in order, so i dont have to test each one 00678 //to make sure it is before or after the previous points 00679 if(radius>0.3 && radius<1.0 && plane != prevPlane) ++fTrackInfo->planesIn10Meter; 00680 if(radius>0.3 && radius<1.5 && plane != prevPlane) ++fTrackInfo->planesIn15Meter; 00681 if(radius>0.3 && radius<2.0 && plane != prevPlane) ++fTrackInfo->planesIn20Meter; 00682 if(radius>0.3 && radius<2.5 && plane != prevPlane) ++fTrackInfo->planesIn25Meter; 00683 if(radius>0.3 && radius<3.0 && plane != prevPlane) ++fTrackInfo->planesIn30Meter; 00684 if(radius>0.3 && radius<3.5 && plane != prevPlane) ++fTrackInfo->planesIn35Meter; 00685 if(radius>0.3 && radius<4.0 && plane != prevPlane) ++fTrackInfo->planesIn40Meter; 00686 00687 ++fTrackInfo->totalStrips; 00688 00689 } // end if the strip had double sided readout 00690 00691 prevPlane = plane; 00692 00693 } // end loop over strips in the track 00694 00695 fTrackInfo->twoEndStripFraction = 0.; 00696 if(numStrips>0){ 00697 fTrackInfo->twoEndStripFraction = (1.*numDoubleEndedStrips)/(1.*numStrips); 00698 fTrackInfo->alternateChi2 /= (1.*numStrips); 00699 } 00700 }//end if the track length is reasonable 00701 else fGoodTrack = false; 00702 00703 MSG("CondensedNtpModuleAtm", Msg::kDebug) << "end fill track info - good track = " 00704 << (int)fGoodTrack << endl; 00705 00706 return; 00707 }
| void CondensedNtpModuleAtm::FillTrackMagneticBendingVariables | ( | NtpSRTrack * | ntpTrack, | |
| TClonesArray * | stripArray | |||
| ) | [private] |
Definition at line 852 of file CondensedNtpModuleAtm.cxx.
References NtpSRTrack::end, fGoodTrack, fTrackInfo, Msg::kDebug, ANtpTrackInfoAtm::meanDistFromLinearFit, MSG, ANtpTrackInfoAtm::netDistFromLinearFit, Munits::ns, NtpSRTrack::nstrip, ANtpTrackInfoAtm::rmsDistFromLinearFit, ANtpTrackInfoAtm::sagitta, NtpSRTrack::stp, NtpSRTrack::stpx, NtpSRTrack::stpy, NtpSRTrack::stpz, NtpSRTrack::vtx, NtpSRVertex::x, NtpSRVertex::y, NtpSRVertex::z, and ANtpTrackInfoAtm::zeroCurvatureChi2.
Referenced by Ana().
00854 { 00855 00856 MSG("CondensedNtpModuleAtm", Msg::kDebug) << "in FillTrackMagneticBendingVariables" 00857 << endl; 00858 00859 NtpSRStrip *ntpStrip = 0; 00860 00861 00862 //find the direction cosines for the straight line between the vertex and end points 00863 //of the track 00864 Float_t deltaX = ntpTrack->end.x - ntpTrack->vtx.x; 00865 Float_t deltaY = ntpTrack->end.y - ntpTrack->vtx.y; 00866 Float_t deltaZ = ntpTrack->end.z - ntpTrack->vtx.z; 00867 Float_t deltaS = TMath::Sqrt(deltaX*deltaX + deltaY*deltaY + deltaZ*deltaZ); 00868 Float_t dcosX = deltaX/deltaS; 00869 Float_t dcosY = deltaY/deltaS; 00870 Float_t dcosZ = deltaZ/deltaS; 00871 Float_t xPos = 0.; 00872 Float_t yPos = 0.; 00873 Float_t zPos = 0.; 00874 Float_t linXPos = 0.; 00875 Float_t linYPos = 0.; 00876 00877 //loop over the strips in the track to find the sagitta 00878 fTrackInfo->sagitta = 0.; 00879 fTrackInfo->netDistFromLinearFit = 0.; 00880 fTrackInfo->meanDistFromLinearFit = 0.; 00881 fTrackInfo->rmsDistFromLinearFit = 0.; 00882 fTrackInfo->zeroCurvatureChi2 = 0.; 00883 00884 if(deltaZ == 0. || deltaX == 0. || deltaY == 0.){ 00885 fGoodTrack = false; 00886 return; 00887 } 00888 00889 TH1F dist("dist", "", 500, 0., 5.); 00890 00891 Float_t curDistFromLinearFit = 0.; 00892 00893 //loop over all strips in the track 00894 Int_t numStrips = (int)ntpTrack->nstrip; 00895 Int_t index = 0; 00896 double stripPosError = 0.041/TMath::Sqrt(12.); 00897 00898 for(Int_t ns = 0; ns < numStrips; ++ns){ 00899 00900 //get the index for this strip 00901 index = ntpTrack->stp[ns]; 00902 00903 //get the NtpSRStrip object 00904 ntpStrip = dynamic_cast<NtpSRStrip *>(stripArray->At(index)); 00905 00906 xPos = ntpTrack->stpx[ns]; 00907 yPos = ntpTrack->stpy[ns]; 00908 zPos = ntpTrack->stpz[ns]; 00909 00910 linXPos = ntpTrack->vtx.x + (dcosX/dcosZ) * (zPos - ntpTrack->vtx.z); 00911 linYPos = ntpTrack->vtx.y + (dcosY/dcosZ) * (zPos - ntpTrack->vtx.z); 00912 deltaX = xPos - linXPos; 00913 deltaY = yPos - linYPos; 00914 00915 fTrackInfo->zeroCurvatureChi2 += (deltaX*deltaX + deltaY*deltaY)/(stripPosError*stripPosError); 00916 00917 curDistFromLinearFit = TMath::Sqrt(deltaX*deltaX + deltaY*deltaY); 00918 dist.Fill(curDistFromLinearFit); 00919 00920 if(curDistFromLinearFit > fTrackInfo->sagitta) fTrackInfo->sagitta = curDistFromLinearFit; 00921 00922 } 00923 00924 if(numStrips>1) fTrackInfo->zeroCurvatureChi2 /= 1.*(numStrips-1); 00925 fTrackInfo->meanDistFromLinearFit = dist.GetMean(); 00926 fTrackInfo->rmsDistFromLinearFit = dist.GetRMS(); 00927 00928 //now to figure out the net pt kick to the muon 00929 00930 return; 00931 }
| void CondensedNtpModuleAtm::FillTrackTimingVariables | ( | NtpSRTrack * | ntpTrack, | |
| TClonesArray * | stripArray | |||
| ) | [private] |
Definition at line 710 of file CondensedNtpModuleAtm.cxx.
References ANtpHeaderInfo::dataType, fGoodTrack, fHeaderInfo, FitMinimizingResidual(), ANtpTrackInfoAtm::flatLineChi2, fTrackInfo, GetTimeWeight(), ANtpTrackInfoAtm::invBeta, ANtpTrackInfoAtm::invBetaChi2, Msg::kDebug, SimFlag::kMC, ANtpTrackInfo::length, MSG, NtpSRStrip::ndigit, Munits::ns, NtpSRTrack::nstrip, NtpSRStrip::ph0, NtpSRStrip::ph1, NtpSRStrip::plane, ANtpHeaderInfo::run, NtpSRPulseHeight::sigcor, ANtpHeaderInfo::snarl, NtpSRTrack::stp, NtpSRTrack::stpds, NtpSRTrack::stpt0, NtpSRTrack::stpt1, NtpSRTrack::stpx, NtpSRTrack::stpy, NtpSRTrack::stpz, NtpSRStrip::strip, ANtpHeaderInfo::subRun, ANtpTrackInfoAtm::timeSlopeY, and ANtpTrackInfoAtm::timeSlopeYChi2.
Referenced by Ana().
00712 { 00713 MSG("CondensedNtpModuleAtm", Msg::kDebug) << "in FillTrackTimingVariables" << endl; 00714 00715 NtpSRStrip *ntpStrip = 0; 00716 00717 Double_t time[2] = {0.}; 00718 Double_t xPos = 0., yPos = 0., zPos = 0.; 00719 Double_t yPosition[5000] = {0.}; 00720 Double_t yPosition1[5000] = {0.}; 00721 Double_t weight[5000] = {0.}; 00722 Double_t weight2[5000] = {0.}; 00723 Double_t pathLength[5000] = {0.}; 00724 Double_t weight1[5000] = {0.}; 00725 Double_t timeAv[5000] = {0.}; 00726 Double_t timeAv1[5000] = {0.}; 00727 Int_t arrayCtr = 0, betaCtr = 0, actr = 0, plane = 0, strip = 0, index = 0; 00728 Double_t minTime = 1.e20; 00729 00730 //loop over all strips in the track 00731 Int_t numStrips = (int)ntpTrack->nstrip; 00732 Int_t numDigits = 0; 00733 Double_t smear = 0.; 00734 00735 MSG("CondensedNtpModuleAtm", Msg::kDebug) << "fill track timing" << endl; 00736 00737 for(Int_t ns = 0; ns < numStrips; ++ns){ 00738 00739 if(fHeaderInfo->dataType == (int)SimFlag::kMC) 00740 smear = gRandom->Gaus(0., 0.75); 00741 else smear = 0.; 00742 00743 //get the index for this strip 00744 index = ntpTrack->stp[ns]; 00745 00746 //get the NtpSRStrip object 00747 ntpStrip = dynamic_cast<NtpSRStrip *>(stripArray->At(index)); 00748 plane = (int)ntpStrip->plane; 00749 strip = (int)ntpStrip->strip; 00750 numDigits = (int)ntpStrip->ndigit; 00751 00752 xPos = ntpTrack->stpx[ns]; 00753 yPos = ntpTrack->stpy[ns]; 00754 zPos = ntpTrack->stpz[ns]; 00755 00756 if(numDigits == 2){ 00757 00758 time[0] = ntpTrack->stpt0[ns]; 00759 time[1] = ntpTrack->stpt1[ns]; 00760 00761 if(time[0] < minTime) minTime = time[0]; 00762 if(time[1] < minTime) minTime = time[1]; 00763 00764 //have a point for each end of the strip 00765 timeAv[arrayCtr] = 1.e9*time[0] + smear; 00766 timeAv[arrayCtr+1] = 1.e9*time[1] + smear; 00767 pathLength[arrayCtr] = fTrackInfo->length - ntpTrack->stpds[ns]; 00768 pathLength[arrayCtr+1] = fTrackInfo->length - ntpTrack->stpds[ns]; 00769 00770 //stuff to get timing information 00771 yPosition[arrayCtr] = yPos; 00772 weight[arrayCtr] = GetTimeWeight(ntpStrip->ph0.sigcor); 00773 if (weight[arrayCtr]>0.4) weight[arrayCtr]=0.4; 00774 yPosition[arrayCtr+1] = yPos; 00775 weight[arrayCtr+1] = GetTimeWeight(ntpStrip->ph1.sigcor); 00776 if (weight[arrayCtr+1]>0.4) weight[arrayCtr+1]=0.4; 00777 MSG("CondensedNtpModuleAtm", Msg::kDebug) << weight[arrayCtr] << " " 00778 << weight[arrayCtr+1] << endl; 00779 arrayCtr += 2; 00780 00781 00782 } // end if the strip had double sided readout 00783 } // end loop over strips in the track 00784 00785 //subtract the minimum time and copy the weights into another array to do the fit to 00786 //a flat line later 00787 00788 for(Int_t ti = 0; ti< arrayCtr; ++ti){ 00789 timeAv[ti] -= 1.e9*minTime; 00790 weight2[ti] = weight[ti]; 00791 } 00792 00793 //cout << "got strip info" << endl; 00794 00795 //array for fit parameters 00796 Double_t param[] = {-10000., -10000.}; 00797 Double_t eparam[] = {-10000., -10000.}; 00798 00799 //find the inverse beta value 00800 betaCtr = arrayCtr; 00801 FitMinimizingResidual(betaCtr, pathLength, timeAv, 00802 weight, param, eparam); 00803 fTrackInfo->invBeta = param[1]*0.3; 00804 if(betaCtr > 2) fTrackInfo->invBetaChi2 = eparam[0] /(1.*(betaCtr-2)); 00805 00806 //find the chi^2 for a fit to a straight line so you 00807 //can see if the sloped fit is really better 00808 FitMinimizingResidual(betaCtr, pathLength, timeAv, 00809 weight2, param, eparam, false); 00810 if(betaCtr > 1) fTrackInfo->flatLineChi2 = eparam[0] / (1.*betaCtr-1); 00811 00812 //now loop over the points in the array and fill new arrays 00813 //containing only those points that have a nonzero weight 00814 actr = 0; 00815 MSG("CondensedNtpModuleAtm", Msg::kDebug) << fHeaderInfo->run << " " 00816 << fHeaderInfo->subRun 00817 << " " << fHeaderInfo->snarl 00818 << " " << arrayCtr << endl; 00819 for(Int_t ii = 0; ii < arrayCtr; ++ii){ 00820 MSG("CondensedNtpModuleAtm", Msg::kDebug) << " " << ii << "/" << arrayCtr << " " 00821 << yPosition[ii] << " " << weight[ii] 00822 << " " << timeAv[ii] << " " 00823 << actr << endl; 00824 if(weight[ii]>0.){ 00825 yPosition1[actr] = yPosition[ii]; 00826 weight1[actr] = weight[ii]; 00827 timeAv1[actr] = timeAv[ii]; 00828 ++actr; 00829 } 00830 }//end loop to select good points along the track 00831 00832 //get the fit for the time in the y direction 00833 FitMinimizingResidual(actr, yPosition1, timeAv1, weight1, param, eparam); 00834 fTrackInfo->timeSlopeY = param[1]; 00835 fTrackInfo->timeSlopeYChi2 = eparam[0]; 00836 00837 if(actr > 2) fTrackInfo->timeSlopeYChi2 /= 1.*(actr-2); 00838 00839 MSG("CondensedNtpModuleAtm", Msg::kDebug) << "actr = " << actr 00840 << "/" << arrayCtr 00841 << " goodtrack = " 00842 << (int)fGoodTrack << endl; 00843 00844 if(1.*actr < 0.5*arrayCtr) fGoodTrack = false; 00845 00846 return; 00847 }
| void CondensedNtpModuleAtm::FindChi2 | ( | Double_t * | x, | |
| Double_t * | y, | |||
| Double_t * | weight, | |||
| Int_t | nPoints, | |||
| Double_t & | a1, | |||
| Double_t & | a2, | |||
| Double_t & | chiSq | |||
| ) | [private] |
Definition at line 1277 of file CondensedNtpModuleAtm.cxx.
01281 { 01282 01283 double vtxX = x[0]; 01284 double vtxY = y[0]; 01285 double endX = x[nPoints-1]; 01286 double endY = y[nPoints-1]; 01287 01288 a2 = endY-vtxY; 01289 if(TMath::Abs(vtxX-endX)>0.) a2 /= (endX-vtxX); 01290 else{ 01291 cout << "cannot find slope for n = " << nPoints << " points " << endl; 01292 for(Int_t i = 0; i<nPoints; ++i) 01293 cout << x[i] << " " << y[i] << " " << weight[i] << endl; 01294 } 01295 01296 a1 = endY-a2*endX; 01297 01298 //find the chi^2 value for the fit 01299 chiSq = 0.; 01300 Double_t chi = 0.; 01301 for(Int_t j = 0; j<nPoints; j++){ 01302 01303 chi = (y[j] - (a1 + a2*x[j]))/weight[j]; 01304 01305 // cout << y[j] << " " << a1+a2*x[j] << " " << chiSq << endl; 01306 01307 chiSq += chi*chi; 01308 01309 } 01310 01311 // cout << "chi^2 = " << chiSq << " " << a1 << " " << a2 << endl; 01312 01313 //get the chiSq per degree of freedom 01314 01315 if(nPoints>2) chiSq /= 1.*(nPoints-2); 01316 01317 return; 01318 }
| void CondensedNtpModuleAtm::FindLinearFit | ( | Double_t * | x, | |
| Double_t * | y, | |||
| Double_t * | weight, | |||
| Int_t | nPoints, | |||
| Double_t & | a1, | |||
| Double_t & | a2, | |||
| Double_t & | chiSq | |||
| ) | [private] |
Definition at line 1179 of file CondensedNtpModuleAtm.cxx.
01183 { 01184 01185 //use method of determinants to do the least squares fit 01186 //this is from bevington chapter 6 01187 // 01188 //f_1(x) = 1 01189 //f_2(x) = x 01190 01191 //need to find the sums 01192 //SIGMA y_i*f_1(x_i)/weight_i 01193 //SIGMA y_i*f_2(x_i)/weight_i 01194 //SIGMA f_1(x_i)*f_1(x_i)/weight_i 01195 //SIGMA f_1(x_i)*f_2(x_i)/weight_i 01196 //SIGMA f_2(x_i)*f_1(x_i)/weight_i 01197 //SIGMA f_2(x_i)*f_2(x_i)/weight_i 01198 // 01199 //some of these are the same, ie f_2*f_1 = f_1*f_2 so only need 01200 //3 instead of 4 variables 01201 01202 Double_t yf1divWeight = 0.; 01203 Double_t yf2divWeight = 0.; 01204 Double_t f1f1divWeight = 0.; 01205 Double_t f1f2divWeight = 0.; 01206 Double_t f2f2divWeight = 0.; 01207 01208 for(Int_t i=0; i<nPoints; i++){ 01209 // cout << y[i] << " " << x[i] << " " << weight[i] << endl; 01210 01211 yf1divWeight += y[i]/(weight[i]*weight[i]); 01212 yf2divWeight += (y[i]*x[i])/(weight[i]*weight[i]); 01213 f1f1divWeight += 1./(weight[i]*weight[i]); 01214 f1f2divWeight += x[i]/(weight[i]*weight[i]); 01215 f2f2divWeight += (x[i]*x[i])/(weight[i]*weight[i]); 01216 } 01217 01218 //make a bunch of matricies to hold these values 01219 TMatrix deltaMatrix(2,2); 01220 TMatrix a1Matrix(2,2); 01221 TMatrix a2Matrix(2,2); 01222 01223 //fill the matricies 01224 deltaMatrix(0,0) = f1f1divWeight; 01225 deltaMatrix(0,1) = f1f2divWeight; 01226 deltaMatrix(1,0) = f1f2divWeight; 01227 deltaMatrix(1,1) = f2f2divWeight; 01228 01229 a1Matrix(0,0) = yf1divWeight; 01230 a1Matrix(0,1) = f1f2divWeight; 01231 a1Matrix(1,0) = yf2divWeight; 01232 a1Matrix(1,1) = f2f2divWeight; 01233 01234 a2Matrix(0,0) = f1f1divWeight; 01235 a2Matrix(0,1) = yf1divWeight; 01236 a2Matrix(1,0) = f1f2divWeight; 01237 a2Matrix(1,1) = yf2divWeight; 01238 01239 a1 = -10000.; 01240 a2 = -10000.; 01241 Double_t delta = deltaMatrix.Determinant(); 01242 01243 // cout << "delta = " << delta << endl; 01244 01245 if(delta != 0.){ 01246 a1 = a1Matrix.Determinant()/delta; 01247 a2 = a2Matrix.Determinant()/delta; 01248 } 01249 01250 //find the chi^2 value for the fit 01251 chiSq = 0.; 01252 Double_t chi = 0.; 01253 for(Int_t j = 0; j<nPoints; j++){ 01254 01255 chi = (y[j] - (a1 + a2*x[j]))/weight[j]; 01256 01257 // cout << chi*weight[j] << endl; 01258 01259 chiSq += chi*chi; 01260 01261 } 01262 01263 // cout << "chi^2 = " << chiSq << " " << a1 << " " << a2 << endl; 01264 01265 //get the chiSq per degree of freedom 01266 01267 chiSq /= (1.*nPoints - 2.); 01268 01269 return; 01270 }
| void CondensedNtpModuleAtm::FindMeanAndRMS | ( | double * | x, | |
| double * | weight, | |||
| int | nPoints, | |||
| double & | mean, | |||
| double & | rms | |||
| ) | [private] |
Definition at line 1366 of file CondensedNtpModuleAtm.cxx.
01369 { 01370 01371 // cout << "rms" << endl; 01372 01373 double sumXXW = 0.; 01374 double sumXW = 0.; 01375 double sumW = 0.; 01376 01377 for(Int_t i = 0; i<nPoints; ++i){ 01378 sumXXW += x[i]*x[i]*weight[i]; 01379 sumXW += x[i]*weight[i]; 01380 sumW += weight[i]; 01381 // cout << sumXW << " " << sumXXW << " " << sumW << " " 01382 // << x[i] << " " << weight[i] << endl; 01383 } 01384 01385 if(sumW > 0. && nPoints>1){ 01386 mean = sumXW/sumW; 01387 if(sumXXW >=0.) rms = TMath::Sqrt((sumXXW/sumW - mean*mean) 01388 *(nPoints*1./(1.*nPoints-1.))); 01389 } 01390 else{ 01391 rms = 0.; 01392 mean = x[0]; 01393 } 01394 return; 01395 }
| void CondensedNtpModuleAtm::FindStraightLineDeviation | ( | Double_t * | x, | |
| Double_t * | y, | |||
| Int_t | nPoints, | |||
| Double_t & | a1, | |||
| Double_t & | a2, | |||
| Double_t & | deviation | |||
| ) | [private] |
Definition at line 1325 of file CondensedNtpModuleAtm.cxx.
01331 { 01332 01333 double vtxX = x[0]; 01334 double vtxY = y[0]; 01335 double endX = x[nPoints-1]; 01336 double endY = y[nPoints-1]; 01337 01338 a2 = endY-vtxY; 01339 if(TMath::Abs(vtxX-endX)>0.) a2 /= (endX-vtxX); 01340 else{ 01341 cout << "cannot find slope for n = " << nPoints << " points " << endl; 01342 for(Int_t i = 0; i<nPoints; ++i) cout << x[i] << " " << y[i] << endl; 01343 } 01344 01345 a1 = endY-a2*endX; 01346 01347 //find the chi^2 value for the fit 01348 deviation = 0.; 01349 Double_t dev = 0.; 01350 for(Int_t j = 0; j<nPoints; j++){ 01351 01352 dev = (y[j] - (a1 + a2*x[j])); 01353 01354 // cout << y[j] << " " << a1+a2*x[j] << " " << chiSq << endl; 01355 01356 deviation += dev*dev; 01357 01358 } 01359 01360 return; 01361 }
| void CondensedNtpModuleAtm::FitMinimizingResidual | ( | Int_t & | n, | |
| Double_t * | x, | |||
| Double_t * | y, | |||
| Double_t * | w, | |||
| Double_t * | param, | |||
| Double_t * | eparam, | |||
| bool | findSlope = true | |||
| ) | [private] |
Definition at line 1118 of file CondensedNtpModuleAtm.cxx.
References LinearFit_Weighted(), and WeightedAverage().
Referenced by FillTrackTimingVariables().
01123 { 01124 // cout << "in LinearFit_MinimizeResid" << endl; 01125 01126 Int_t nremove = 0; 01127 Double_t maxresid = 11.; 01128 Double_t resid = 0.; 01129 Double_t parm[2] = {0.}; 01130 Double_t eparm[2] = {0.}; 01131 Double_t timefitchi2 = 0.; 01132 while(maxresid>10. 01133 &&n-nremove>4 01134 &&(1.*nremove)/(1.*n)<0.2){ 01135 01136 timefitchi2=0.; 01137 if(findSlope) LinearFit_Weighted(n,x,y,w,parm,eparm); 01138 else WeightedAverage(n,y,w,parm,eparm); 01139 maxresid = 0.; 01140 Int_t imaxresid=0; 01141 for (Int_t i=0; i<n; i++) { 01142 if (w[i]>0.) { 01143 if(findSlope) resid = (parm[1]*x[i]+parm[0]-y[i]); 01144 else resid = parm[0] - y[i]; 01145 timefitchi2 += w[i]*resid*resid; 01146 if (TMath::Abs(resid)>maxresid) { 01147 maxresid = TMath::Abs(resid); 01148 imaxresid = i; 01149 } 01150 } 01151 } 01152 if (maxresid>10.) { 01153 w[imaxresid] = 0.; 01154 nremove++; 01155 } 01156 }//end loop to minimize residuals 01157 01158 //get the chi^2 01159 timefitchi2 = 0.; 01160 for (Int_t i=0; i<n; i++) { 01161 if (w[i]>0.) { 01162 if(findSlope) resid = (parm[1]*x[i]+parm[0]-y[i]); 01163 else resid = parm[0] - y[i]; 01164 timefitchi2 += w[i]*resid*resid; 01165 } 01166 } 01167 01168 01169 n -= nremove; 01170 param[0] = parm[0]; 01171 param[1] = parm[1]; 01172 eparam[0] = timefitchi2; 01173 eparam[1] = eparm[1]; 01174 01175 return; 01176 }
| Double_t CondensedNtpModuleAtm::GetTimeWeight | ( | Float_t | ph | ) | [private] |
Definition at line 1000 of file CondensedNtpModuleAtm.cxx.
References ArrivalTime_Weight().
Referenced by FillTrackTimingVariables().
01001 { 01002 // cout << "in GetTimeWeight" << endl; 01003 01004 Double_t npe = ph/60.; 01005 if (npe<1.) npe=1.; 01006 01007 // cout << ArrivalTime_Weight(npe) << endl; 01008 01009 return ArrivalTime_Weight(npe); 01010 }
| void CondensedNtpModuleAtm::Help | ( | ) | [virtual] |
Implement to spew some useful help to cout
Reimplemented from JobCModule.
Definition at line 410 of file CondensedNtpModuleAtm.cxx.
References Msg::kInfo, and MSG.
00411 { 00412 MSG("JobC", Msg::kInfo) 00413 << "NearElectronicsCheck::Help\n" 00414 <<"CondensedNtpModuleAtm is a module which demultiplexes events " 00415 <<"in the far detector." 00416 << endl; 00417 }
| void CondensedNtpModuleAtm::LinearFit_Weighted | ( | const Int_t | n, | |
| const Double_t * | x, | |||
| const Double_t * | y, | |||
| const Double_t * | w, | |||
| Double_t * | parm, | |||
| Double_t * | eparm | |||
| ) | [private] |
Definition at line 1050 of file CondensedNtpModuleAtm.cxx.
Referenced by FitMinimizingResidual().
01055 { 01056 // cout << "in LinearFit_Weighted with errors" << endl; 01057 01058 Double_t sumx=0.; 01059 Double_t sumx2=0.; 01060 Double_t sumy=0.; 01061 Double_t sumy2=0.; 01062 Double_t sumxy=0.; 01063 Double_t sumw=0.; 01064 01065 for(Int_t i=0; i<n; ++i){ 01066 sumx += x[i]*w[i]; 01067 sumx2 += x[i]*x[i]*w[i]; 01068 sumy += y[i]*w[i]; 01069 sumy2 += y[i]*y[i]*w[i]; 01070 sumxy += x[i]*y[i]*w[i]; 01071 sumw += w[i]; 01072 01073 // cout << x[i] << " " << y[i] << " " << w[i] << endl; 01074 } 01075 01076 // cout << sumx2 << " " << sumw << " " << sumx << " " << sumxy << " " << sumy 01077 // << endl; 01078 01079 if(sumx2*sumw-sumx*sumx>0.&&sumx2<1.e5){ 01080 parm[0] = (sumy*sumx2-sumx*sumxy)/(sumx2*sumw-sumx*sumx); 01081 parm[1] = (sumxy*sumw-sumx*sumy)/(sumx2*sumw-sumx*sumx); 01082 01083 eparm[0] = TMath::Sqrt(sumx2*(sumx2*sumw-sumx*sumx))/(sumx2*sumw-sumx*sumx); 01084 eparm[1] = TMath::Sqrt(sumx2*sumw-sumx*sumx)/(sumx2*sumw-sumx*sumx); 01085 } 01086 01087 return; 01088 }
| void CondensedNtpModuleAtm::ResetTreeVariables | ( | ) | [private] |
Definition at line 989 of file CondensedNtpModuleAtm.cxx.
References fEventInfo, fGoodTrack, fTrackInfo, fTruthInfo, ANtpTruthInfoAtm::Reset(), ANtpTrackInfoAtm::Reset(), and ANtpEventInfo::Reset().
Referenced by Ana().
00990 { 00991 fEventInfo->Reset(); 00992 fTrackInfo->Reset(); 00993 fTruthInfo->Reset(); 00994 fGoodTrack = false; 00995 00996 return; 00997 }
| void CondensedNtpModuleAtm::WeightedAverage | ( | const Int_t | n, | |
| const Double_t * | x, | |||
| const Double_t * | w, | |||
| Double_t * | parm, | |||
| Double_t * | eparm | |||
| ) | [private] |
Definition at line 1091 of file CondensedNtpModuleAtm.cxx.
Referenced by FitMinimizingResidual().
01095 { 01096 // cout << "in LinearFit_Weighted with errors" << endl; 01097 01098 double sumWeights = 0.; 01099 double sumXWeights = 0.; 01100 //find the weighted average of the points 01101 for(int i = 0; i < n; ++i){ 01102 sumXWeights += x[i]*w[i]; 01103 sumWeights += w[i]; 01104 } 01105 01106 if(sumWeights > 0.) parm[0] = sumXWeights/sumWeights; 01107 else parm[0] = -9999.; 01108 parm[1] = 0.; 01109 01110 eparm[0] = 0.; 01111 eparm[1] = 0.; 01112 01113 01114 return; 01115 }
Double_t CondensedNtpModuleAtm::fBeamHPos [private] |
Definition at line 72 of file CondensedNtpModuleAtm.h.
ANtpBeamInfo* CondensedNtpModuleAtm::fBeamInfo [private] |
Definition at line 88 of file CondensedNtpModuleAtm.h.
Referenced by Ana(), BeginJob(), and ~CondensedNtpModuleAtm().
std::string CondensedNtpModuleAtm::fBeamPath [private] |
Double_t CondensedNtpModuleAtm::fBeamTimeSec [private] |
Definition at line 75 of file CondensedNtpModuleAtm.h.
TTree* CondensedNtpModuleAtm::fBeamTree [private] |
Definition at line 65 of file CondensedNtpModuleAtm.h.
std::string CondensedNtpModuleAtm::fBeamTreeName [private] |
Double_t CondensedNtpModuleAtm::fBeamVPos [private] |
Definition at line 73 of file CondensedNtpModuleAtm.h.
Double_t CondensedNtpModuleAtm::fColDenCosZen[1500] [private] |
Int_t CondensedNtpModuleAtm::fCtr [private] |
Double_t CondensedNtpModuleAtm::fCurrentPOT [private] |
Definition at line 71 of file CondensedNtpModuleAtm.h.
Int_t CondensedNtpModuleAtm::fDataType [private] |
Definition at line 68 of file CondensedNtpModuleAtm.h.
Referenced by Ana(), BeginJob(), and Config().
Double_t CondensedNtpModuleAtm::fDensity[1500] [private] |
Definition at line 83 of file CondensedNtpModuleAtm.h.
Referenced by BeginJob(), and FillTrackInformation().
Int_t CondensedNtpModuleAtm::fDetector [private] |
Definition at line 66 of file CondensedNtpModuleAtm.h.
Referenced by Ana(), Config(), and FillTrackInformation().
Int_t CondensedNtpModuleAtm::fEndPlane [private] |
Definition at line 69 of file CondensedNtpModuleAtm.h.
Referenced by Ana(), Config(), and FillTrackInformation().
ANtpEventInfo* CondensedNtpModuleAtm::fEventInfo [private] |
Definition at line 86 of file CondensedNtpModuleAtm.h.
Referenced by Ana(), BeginJob(), FillEventInformation(), FillTrackInformation(), ResetTreeVariables(), and ~CondensedNtpModuleAtm().
Int_t CondensedNtpModuleAtm::fFailDeMux [private] |
TTree* CondensedNtpModuleAtm::fFailTree [private] |
Definition at line 64 of file CondensedNtpModuleAtm.h.
Referenced by Ana(), BeginJob(), and EndJob().
std::string CondensedNtpModuleAtm::fFileName [private] |
Int_t CondensedNtpModuleAtm::fGoodTrack [private] |
Definition at line 67 of file CondensedNtpModuleAtm.h.
Referenced by Ana(), FillTrackInformation(), FillTrackMagneticBendingVariables(), FillTrackTimingVariables(), and ResetTreeVariables().
ANtpHeaderInfo* CondensedNtpModuleAtm::fHeaderInfo [private] |
Definition at line 87 of file CondensedNtpModuleAtm.h.
Referenced by Ana(), BeginJob(), FillTrackTimingVariables(), and ~CondensedNtpModuleAtm().
Double_t CondensedNtpModuleAtm::fHornCurrent [private] |
Definition at line 70 of file CondensedNtpModuleAtm.h.
Definition at line 92 of file CondensedNtpModuleAtm.h.
Referenced by Ana(), FillEventInformation(), FillMCInformation(), and FillTrackInformation().
Int_t CondensedNtpModuleAtm::fIsMultiple [private] |
Int_t CondensedNtpModuleAtm::fMajorRelease [private] |
bool CondensedNtpModuleAtm::fNewMC [private] |
TFile* CondensedNtpModuleAtm::fNtpFile [private] |
TTree* CondensedNtpModuleAtm::fNtuple [private] |
Definition at line 63 of file CondensedNtpModuleAtm.h.
Referenced by Ana(), BeginJob(), and EndJob().
Double_t CondensedNtpModuleAtm::fPathAzimuthDeg[1500] [private] |
Double_t CondensedNtpModuleAtm::fPathCosZenDeg[1500] [private] |
std::string CondensedNtpModuleAtm::fRockMapFileName [private] |
Double_t CondensedNtpModuleAtm::fTargetPos [private] |
Definition at line 74 of file CondensedNtpModuleAtm.h.
Definition at line 89 of file CondensedNtpModuleAtm.h.
Referenced by Ana(), BeginJob(), FillTrackInformation(), FillTrackMagneticBendingVariables(), FillTrackTimingVariables(), ResetTreeVariables(), and ~CondensedNtpModuleAtm().
std::string CondensedNtpModuleAtm::fTreeName [private] |
Definition at line 90 of file CondensedNtpModuleAtm.h.
Referenced by BeginJob(), FillMCInformation(), ResetTreeVariables(), and ~CondensedNtpModuleAtm().
Int_t CondensedNtpModuleAtm::fValidPlaneFail [private] |
1.4.7