#include "TROOT.h"#include <sstream>#include <string>#include <vector>#include "TTree.h"Go to the source code of this file.
Namespaces | |
| namespace | VHS |
Enumerations | |
| enum | evtType { vhsNC = 0, vhsCCe = 1, vhsCCmu = 2, vhsCCtau = 3, vhsUnknown = 4 } |
Functions | |
| bool | DrawEvent (int nEntry, const char *treeName, int nPlanes, int nStrips) |
| void | FillDiscriminants (NtpSREvent *evt, TClonesArray *stp, std::vector< double > avgNC, std::vector< double > avgCCe, std::vector< double > avgCCmu, std::vector< double > avgCCtau, std::vector< double > varNC, std::vector< double > varCCe, std::vector< double > varCCmu, std::vector< double > varCCtau, std::vector< double > pNC, std::vector< double > pCCe, std::vector< double > pCCmu, std::vector< double > pCCtau, const int nPlanes, const int nStrips, const bool bUnit, VHSevent *&vhsevt) |
| std::vector< double > | FindMedian (std::vector< std::vector< double > > pts, std::vector< double > avg, bool bUnit=true, bool bVerbose=true) |
| std::vector< double > | FullVec (std::vector< double > image, std::vector< int > index, int nPlanes, int nStrips) |
| double | GetDistance (std::vector< double > vec0, std::vector< double > vec1) |
| VHS::evtType | GetEvtType (int inu, int iaction) |
| void | GetImage (int *index, int nstp, TClonesArray *stp, int nPlanes, int nStrips, std::vector< double > &image, std::vector< int > &vecInd, std::vector< double > &theta) |
| double | GetLL (std::vector< double > fullImage, std::vector< double > pHit, std::vector< double > avg, std::vector< double > var) |
| int | GetPlane (int vecInd, int nPlanes) |
| int | GetStrip (int vecInd, int nPlanes) |
| void | GetThetaAxis (TClonesArray *stp, std::vector< int > index, double &theta, std::vector< int > ¢er) |
| std::vector< std::string > | ParseNmList (const char *cstr) |
| void | ReadFile (TFile *inFile, std::vector< double > &avgNC, std::vector< double > &avgCCe, std::vector< double > &avgCCmu, std::vector< double > &avgCCtau, std::vector< double > &varNC, std::vector< double > &varCCe, std::vector< double > &varCCmu, std::vector< double > &varCCtau, std::vector< double > &pNChit, std::vector< double > &pCCehit, std::vector< double > &pCCmuhit, std::vector< double > &pCCtauhit, const int nPlanes) |
| void | RotatePixel (int &plane, int &strip, const int avgPlane, const int avgStrip, const double theta) |
| void | SeparateViews (int *index, int nstp, TClonesArray *stp, std::vector< int > &ustp, std::vector< int > &vstp) |
| std::vector< VHSevent * > | Skim (NtpStRecord *ntpst, std::vector< double > avgNC, std::vector< double > avgCCe, std::vector< double > avgCCmu, std::vector< double > avgCCtau, std::vector< double > varNC, std::vector< double > varCCe, std::vector< double > varCCmu, std::vector< double > varCCtau, std::vector< double > pNC, std::vector< double > pCCe, std::vector< double > pCCmu, std::vector< double > pCCtau, const int nPlanes=20, const int nStrips=20, const bool bUnit=true) |
| std::vector< double > | SubtractVec (std::vector< double > image0, std::vector< double > image1) |
| template<class T> | |
| std::string | ToString (const T &thing, int w=0, int p=0) |
| void | ToTH2D (std::vector< double > vec, int nPlanes, int nStrips, std::string name, std::string title, TH2D *&Uview, TH2D *&Vview) |
| std::vector< double > | ToVector (TH2D *hist, int nPlanes) |
| int | Train (NtpStRecord *ntpst, std::vector< double > &avgNC, std::vector< double > &avgCCe, std::vector< double > &avgCCmu, std::vector< double > &avgCCtau, std::vector< double > &varNC, std::vector< double > &varCCe, std::vector< double > &varCCmu, std::vector< double > &varCCtau, std::vector< int > &numNC, std::vector< int > &numCCe, std::vector< int > &numCCmu, std::vector< int > &numCCtau, int &NCevts, int &eCCevts, int &muCCevts, int &tauCCevts, const int maxTrain=1000, const int nPlanes=20, const int nStrips=20, const int cutPlanes=40, const double eRecoMax=999., const double eTrueMax=999.) |
| void | TrainPost (TFile *outFile, std::vector< double > &avgNC, std::vector< double > &avgCCe, std::vector< double > &avgCCmu, std::vector< double > &avgCCtau, std::vector< double > &varNC, std::vector< double > &varCCe, std::vector< double > &varCCmu, std::vector< double > &varCCtau, std::vector< int > &numNC, std::vector< int > &numCCe, std::vector< int > &numCCmu, std::vector< int > &numCCtau, int NCevts, int eCCevts, int muCCevts, int tauCCevts, const int nPlanes=20, const int nStrips=20) |
| void | UnitVector (std::vector< double > &vec) |
| int | VecIndex (int ipln, int istp, int nPlanes) |
| void | WriteFile (TFile *outFile, std::vector< double > &avgNC, std::vector< double > &avgCCe, std::vector< double > &avgCCmu, std::vector< double > &avgCCtau, std::vector< double > &varNC, std::vector< double > &varCCe, std::vector< double > &varCCmu, std::vector< double > &varCCtau, std::vector< double > &pNChit, std::vector< double > &pCCehit, std::vector< double > &pCCmuhit, std::vector< double > &pCCtauhit, const int nPlanes, const int nStrips) |
|
|
Definition at line 28 of file VHS.h. 00028 { vhsNC=0, vhsCCe=1, vhsCCmu=2, vhsCCtau=3, vhsUnknown=4 };
|
|
||||||||||||||||||||
|
Definition at line 45 of file VHS.cxx. References VHSevent::dCCe, VHSevent::dCCmu, VHSevent::dCCtau, VHSevent::dNC, VHS::GetEvtType(), VHSevent::GeV, VHSevent::ind, VHSevent::MC, VHSevent::MCenu, VHSevent::MCiaction, VHSevent::MCinu, VHSevent::NumPln, VHSevent::pCCe, VHSevent::pCCmu, VHSevent::pCCtau, VHSevent::pNC, VHSevent::Sigcor, VHSevent::stp, VHS::ToString(), VHS::ToTH2D(), VHS::vhsCCmu, and VHS::vhsNC. 00047 {
00048 // Draw an event display for the given event from a TTree of VHSEvent
00049
00050 TTree* theTree = (TTree*)gROOT->FindObject(treeName);
00051
00052 if (!theTree) return false;
00053 const int kNumEntries = theTree->GetEntriesFast();
00054
00055 if ( nEntry > kNumEntries ) return false;
00056
00057 TCanvas* evtDisp = 0;
00058
00059 if ( !gROOT->FindObject("VHSEvtDisp") ) {
00060 evtDisp = new TCanvas("VHSEvtDisp","VHS Event Display",0,0,800,600);
00061 evtDisp->Divide(2,2);
00062 evtDisp->GetPad(1)->SetPad(0.00,0.5,0.33,1.0);
00063 evtDisp->GetPad(2)->SetPad(0.33,0.5,1.00,1.0);
00064 evtDisp->GetPad(3)->SetPad(0.00,0.0,0.33,0.5);
00065 evtDisp->GetPad(4)->SetPad(0.33,0.0,1.00,0.5);
00066
00067 } else {
00068 evtDisp = (TCanvas*)gROOT->FindObject("VHSEvtDisp");
00069 evtDisp->cd(1);
00070 }
00071 evtDisp->SetEditable(kTRUE);
00072
00073 std::string nextCom = (nEntry+1 < kNumEntries)?
00074 std::string("VHS::DrawEvent(") + VHS::ToString(nEntry+1) + std::string(",\"") +
00075 std::string(theTree->GetName()) + "\"," + VHS::ToString(nPlanes) + std::string(",") +
00076 VHS::ToString(nStrips) + std::string(");") :
00077 std::string("VHS::DrawEvent(") + VHS::ToString(kNumEntries-1) + std::string(",\"") +
00078 std::string(theTree->GetName()) + "\"," + VHS::ToString(nPlanes) + std::string(",") +
00079 VHS::ToString(nStrips) + std::string(");");
00080
00081 std::string prevCom = (nEntry > 0)?
00082 std::string("VHS::DrawEvent(") + VHS::ToString(nEntry-1) + std::string(",\"") +
00083 std::string(theTree->GetName()) + "\"," + VHS::ToString(nPlanes) + std::string(",") +
00084 VHS::ToString(nStrips) + std::string(");") :
00085 std::string("VHS::DrawEvent(0,\"") +
00086 std::string(theTree->GetName()) + "\"," + VHS::ToString(nPlanes) + std::string(",") +
00087 VHS::ToString(nStrips) + std::string(");");
00088
00089 // Draw navigation buttons
00090 evtDisp->cd(3);
00091 evtDisp->GetPad(3)->GetListOfPrimitives()->Delete();
00092 evtDisp->GetPad(3)->Range(0,0,1,1);
00093
00094 // Next button
00095 TButton* butNext =
00096 new TButton("Next Event",nextCom.c_str(),0.20,0.60,0.80,0.80);
00097 butNext->SetTextSize(0.3);
00098
00099 // Previous button
00100 TButton* butPrev =
00101 new TButton("Prev Event",prevCom.c_str(),0.20,0.40,0.80,0.60);
00102 butPrev->SetTextSize(0.3);
00103
00104 // Quit button
00105 TButton* butQuit = new TButton("Quit",".q",0.20,0.20,0.80,0.40);
00106 butQuit->SetTextSize(0.3);
00107
00108 butNext->Draw();
00109 butPrev->Draw();
00110 butQuit->Draw();
00111
00112 // Grab a copy of the event
00113 VHSevent* evt = new VHSevent();
00114 theTree->SetBranchAddress("VHSevent",&evt);
00115 theTree->GetEntry(nEntry);
00116
00117 // Draw UZ view
00118 TH2D* uzImage = (TH2D*)gROOT->FindObject("UZImage");
00119 TH2D* vzImage = (TH2D*)gROOT->FindObject("VZImage");
00120
00121 if ( uzImage ) { delete uzImage; uzImage = 0; }
00122 if ( vzImage ) { delete vzImage; vzImage = 0; }
00123
00124 VHS::ToTH2D( VHS::FullVec( evt->stp, evt->ind, nPlanes, nStrips ),
00125 nPlanes, nStrips, "ZImage", "Z",
00126 uzImage, vzImage);
00127
00128 evtDisp->cd(2);
00129 uzImage->SetXTitle("Plane");
00130 uzImage->SetYTitle("Strip");
00131 uzImage->SetStats(kFALSE);
00132 uzImage->Draw("COLZ");
00133
00134 // Draw VZ view
00135 evtDisp->cd(4);
00136 vzImage->SetXTitle("Plane");
00137 vzImage->SetYTitle("Strip");
00138 vzImage->SetStats(kFALSE);
00139 vzImage->Draw("COLZ");
00140
00141 // Draw Event information
00142 evtDisp->cd(1);
00143 evtDisp->GetPad(1)->GetListOfPrimitives()->Delete();
00144
00145 bool MC = evt->MC;
00146 double MCenu = evt->MCenu;
00147 int MCiaction = evt->MCiaction;
00148 int MCinu = TMath::Abs(evt->MCinu);
00149 //int MCrescode = evt->MCrescode;
00150
00151 double dNC = evt->dNC;
00152 double dCCe = evt->dCCe;
00153 double dCCmu = evt->dCCmu;
00154 double dCCtau = evt->dCCtau;
00155
00156 double pNC = evt->pNC;
00157 double pCCe = evt->pCCe;
00158 double pCCmu = evt->pCCmu;
00159 double pCCtau = evt->pCCtau;
00160
00161 double GeV = evt->GeV;
00162 double Sigcor = evt->Sigcor;
00163
00164 int NumPln = evt->NumPln;
00165
00166 VHS::evtType evtt = VHS::GetEvtType( MCinu, MCiaction );
00167
00168 // --> text boxes to dump info
00169 TPaveText* ptL = new TPaveText(0.02,0.02,0.5,0.45,"brNDC");
00170 ptL->SetLineWidth(2);
00171 ptL->SetBorderSize(1);
00172 ptL->SetTextSize(0.05);
00173 ptL->SetTextAlign(12);
00174
00175 std::string li1L = std::string("dCCe = ") + VHS::ToString(dCCe,7,5);
00176 std::string li2L = std::string("dCC#mu = ") + VHS::ToString(dCCmu,7,5);
00177 std::string li3L = std::string("dCC#tau = ") + VHS::ToString(dCCtau,7,5);
00178 std::string li4L = std::string("dNC = ") + VHS::ToString(dNC,7,5);
00179 ptL->AddText(li1L.c_str());
00180 ptL->AddText(li2L.c_str());
00181 ptL->AddText(li3L.c_str());
00182 ptL->AddText(li4L.c_str());
00183 ptL->Draw();
00184
00185 TPaveText* ptR = new TPaveText(0.5,0.02,0.98,0.45,"brNDC");
00186 ptR->SetLineWidth(2);
00187 ptR->SetBorderSize(1);
00188 ptR->SetTextSize(0.05);
00189 ptR->SetTextAlign(12);
00190
00191 std::string li1R = std::string("pCCe = ") + VHS::ToString(pCCe,7,5);
00192 std::string li2R = std::string("pCC#mu = ") + VHS::ToString(pCCmu,7,5);
00193 std::string li3R = std::string("pCC#tau = ") + VHS::ToString(pCCtau,7,5);
00194 std::string li4R = std::string("pNC = ") + VHS::ToString(pNC,7,5);
00195 ptR->AddText(li1R.c_str());
00196 ptR->AddText(li2R.c_str());
00197 ptR->AddText(li3R.c_str());
00198 ptR->AddText(li4R.c_str());
00199 ptR->Draw();
00200
00201 TPaveText* pt = new TPaveText(0.02,0.45,0.98,0.98,"brNDC");
00202 pt->SetLineWidth(2);
00203 pt->SetBorderSize(1);
00204 pt->SetTextSize(0.05);
00205 pt->SetTextAlign(12);
00206
00207 std::string li1 =
00208 std::string("Event #") + VHS::ToString(nEntry+1) +
00209 std::string(" of ") + VHS::ToString(kNumEntries);
00210 std::string li2 = std::string("");
00211 if (evtt == vhsNC ) li2 = "NC event ";
00212 if (evtt == vhsCCe ) li2 = "e CCevent ";
00213 if (evtt == vhsCCmu ) li2 = "#mu CC event ";
00214 if (evtt == vhsCCtau) li2 = "#tau CC event ";
00215 li2 += (MC)?
00216 std::string("(") + VHS::ToString(NumPln) + std::string(" total planes)") :
00217 VHS::ToString(NumPln) + std::string(" total planes") ;
00218 std::string li3 = std::string("E_{vis} = ") + VHS::ToString(GeV,7,5) +
00219 std::string(" GeV (") + VHS::ToString(Sigcor,7,7) + std::string(" pe)");
00220 std::string li4 = (MC)?
00221 std::string("E_{#nu} = ") +
00222 VHS::ToString(MCenu,7,5) +
00223 std::string(" GeV") :
00224 "";
00225 std::string li5 = "";
00226 if ( dNC < dCCmu )
00227 li5 = std::string("Determined: NC event");
00228 if ( dCCmu < dNC )
00229 li5 = std::string("Determined: #mu CC event");
00230 //if ( dNC < dCCe && dNC < dCCmu && dNC < dCCtau )
00231 //li5 = std::string("Determined: NC event");
00232 //if ( dCCe < dNC && dCCe < dCCmu && dCCe < dCCtau )
00233 //li5 = std::string("Determined: e CC event");
00234 //if ( dCCmu < dNC && dCCmu < dCCe && dCCmu < dCCtau )
00235 //li5 = std::string("Determined: #mu CC event");
00236 //if ( dCCtau < dNC && dCCtau < dCCe && dCCtau < dCCmu )
00237 //li5 = std::string("Determined: #tau CC event");
00238 li5 += std::string(" (distance) ");
00239 pt->AddText(li1.c_str());
00240 pt->AddText(li2.c_str());
00241 pt->AddText(li3.c_str());
00242 if (MC) pt->AddText(li4.c_str());
00243 if ( (MC) ){
00244 if ((evtt == vhsNC &&
00245 dNC > dCCmu ) ||
00246 (evtt == vhsCCmu &&
00247 dCCmu > dNC ) )
00248 /*
00249 ( dNC > dCCe || dNC > dCCmu || dNC > dCCtau )) ||
00250 (evtt == vhsCCe &&
00251 ( dCCe > dNC || dCCe > dCCmu || dCCe > dCCtau )) ||
00252 (evtt == vhsCCmu &&
00253 ( dCCmu > dNC || dCCmu > dCCe || dCCmu > dCCtau )) ||
00254 (evtt == vhsCCtau &&
00255 ( dCCtau > dNC || dCCtau > dCCe || dCCtau > dCCmu )) )
00256 */
00257 pt->AddText(li5.c_str())->SetTextColor(2);
00258 else
00259 pt->AddText(li5.c_str());
00260 } else {
00261 pt->AddText(li5.c_str());
00262 }
00263
00264 std::string li6 = "";
00265 if ( pNC > pCCmu )
00266 li6 = std::string("Determined: NC event");
00267 if ( pCCmu > pNC )
00268 li6 = std::string("Determined: #mu CC event");
00269 /*
00270 if ( pNC > pCCe && pNC > pCCmu && pNC > pCCtau )
00271 li6 = std::string("Determined: NC event");
00272 if ( pCCe > pNC && pCCe > pCCmu && pCCe > pCCtau )
00273 li6 = std::string("Determined: e CC event");
00274 if ( pCCmu > pNC && pCCmu > pCCe && pCCmu > pCCtau )
00275 li6 = std::string("Determined: #mu CC event");
00276 if ( pCCtau > pNC && pCCtau > pCCe && pCCtau > pCCmu )
00277 li6 = std::string("Determined: #tau CC event");
00278 */
00279 li6 += std::string(" (loglikelihood) ");
00280 if ( (MC) ) {
00281 if ((evtt == vhsNC &&
00282 pNC < pCCmu ) ||
00283 (evtt == vhsCCmu &&
00284 pCCmu < pNC ) )
00285 /*
00286 ( pNC < pCCe || pNC < pCCmu || pNC < pCCtau )) ||
00287 (evtt == vhsCCe &&
00288 ( pCCe < pNC || pCCe < pCCmu || pCCe < pCCtau )) ||
00289 (evtt == vhsCCmu &&
00290 ( pCCmu < pNC || pCCmu < pCCe || pCCmu < pCCtau )) ||
00291 (evtt == vhsCCtau &&
00292 ( pCCtau < pNC || pCCtau < pCCe || pCCtau < pCCmu )) )
00293 */
00294 pt->AddText(li6.c_str())->SetTextColor(2);
00295 else
00296 pt->AddText(li6.c_str());
00297 } else {
00298 pt->AddText(li6.c_str());
00299 }
00300
00301 pt->Draw();
00302
00303 delete evt; evt = 0;
00304
00305 evtDisp->cd(3);
00306
00307 evtDisp->SetEditable(kFALSE);
00308
00309 return true;
00310
00311 }
|
|
||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
|
Definition at line 315 of file VHS.cxx. References NtpSRPlane::begu, NtpSRPlane::begv, VHSevent::dCCe, VHSevent::dCCmu, VHSevent::dCCtau, VHSevent::dNC, VHS::FullVec(), VHS::GetDistance(), VHS::GetImage(), VHS::GetLL(), VHSevent::ind, VHSevent::MinPlnU, VHSevent::MinPlnV, VHSevent::MOIThetaU, VHSevent::MOIThetaV, NtpSREvent::nstrip, VHSevent::pCCe, VHSevent::pCCmu, VHSevent::pCCtau, NtpSREvent::plane, VHSevent::pNC, VHSevent::stp, NtpSREvent::stp, and VHS::UnitVector(). Referenced by ANtpInfoObjectFillerNC::FillVHSEvtInfo(), and VHS::Skim(). 00333 {
00334
00335 // Gather our image and index sparse vectors
00336 std::vector<double> image;
00337 std::vector<int> index;
00338 std::vector<double> theta;
00339 VHS::GetImage(evt->stp, evt->nstrip,
00340 stp, nPlanes, nStrips, image, index, theta);
00341
00342 vhsevt->stp = image;
00343 vhsevt->ind = index;
00344
00345 vhsevt->MOIThetaU = theta[0];
00346 vhsevt->MOIThetaV = theta[1];
00347
00348 // Get minimum plane in U & V views
00349 vhsevt->MinPlnU = evt->plane.begu;
00350 vhsevt->MinPlnV = evt->plane.begv;
00351
00352 // Fill relative loglikelihoods
00353 std::vector<double> fullImg = FullVec(image,index,nPlanes,nStrips);
00354
00355 vhsevt->pNC =
00356 VHS::GetLL(fullImg, pNC , avgNC , varNC );
00357 vhsevt->pCCe =
00358 VHS::GetLL(fullImg, pCCe , avgCCe , varCCe );
00359 vhsevt->pCCmu =
00360 VHS::GetLL(fullImg, pCCmu , avgCCmu , varCCmu );
00361 vhsevt->pCCtau =
00362 VHS::GetLL(fullImg, pCCtau, avgCCtau, varCCtau);
00363
00364 // Fill mean distances
00365 if (bUnit) {
00366 VHS::UnitVector(fullImg );
00367 VHS::UnitVector(avgNC );
00368 VHS::UnitVector(avgCCe );
00369 VHS::UnitVector(avgCCmu );
00370 VHS::UnitVector(avgCCtau);
00371 }
00372
00373 vhsevt->dNC = VHS::GetDistance( fullImg, avgNC );
00374 vhsevt->dCCe = VHS::GetDistance( fullImg, avgCCe );
00375 vhsevt->dCCmu = VHS::GetDistance( fullImg, avgCCmu );
00376 vhsevt->dCCtau = VHS::GetDistance( fullImg, avgCCtau );
00377
00378 }
|
|
||||||||||||||||||||
|
Definition at line 382 of file VHS.cxx. References VHS::GetDistance(), and VHS::UnitVector(). 00386 {
00387 const unsigned int nDim = avg.size(); // number of dimensions each point lives in
00388 const unsigned int nPts = pt.size(); // number of points for median-finding
00389 const double kEpsilon = 1.e-6; // precision to zero in calculations
00390 const double kDelta = 1.e-4; // cutoff for iterative median displacement
00391
00392 if (bUnit) {
00393 VHS::UnitVector(avg);
00394 for (unsigned int i = 0; i < nPts; i++)
00395 VHS::UnitVector(pt[i]);
00396 }
00397
00398 // Find the median via Weiszfeld iteration
00399 vector<double> median = avg;
00400 vector<double> lastMedian(nDim, 0.);
00401 double dDist = 999.; // iterative difference in avg distance to iterative median
00402 double lastAvgD = 0.; // average distance to the last iterative median
00403 double avgD = 0.; // average distance to the current iterative median
00404
00405 // Find the average distance to our initial guess point (average)
00406 for (unsigned int i = 0; i < nPts; i++) {
00407 double dPt = VHS::GetDistance(median, pt[i]);
00408 avgD += dPt/nPts;
00409 }
00410 if (bVerbose)
00411 cout << "-- Iteration 0 : dItr = N/A : avgDist = " << avgD << " del(N/A)" << endl;
00412
00413 // Iteratively move closer to the true median point
00414 double dItr = 999.; // distance between last iterative median and new one
00415 int nItr = 0; // number of iterations
00416 while( dItr > kDelta ) {
00417 lastMedian = median;
00418 for (unsigned int j = 0; j < nDim; j++) median[j] = 0.;
00419 lastAvgD = avgD;
00420 avgD = 0.;
00421
00422 // Calculate our new median guess
00423 double norm = 0.;
00424 for (unsigned int i = 0; i < nPts; i++) {
00425 double dPt = VHS::GetDistance(lastMedian, pt[i]);
00426
00427 if (!bUnit) norm += (dPt > kEpsilon)? 1./dPt : 1./kEpsilon;
00428
00429 for (unsigned int j = 0; j < nDim; j++)
00430 median[j] += (dPt > kEpsilon)?
00431 pt[i][j]/dPt : pt[i][j]/TMath::Abs(kEpsilon) ;
00432 } // end for loop
00433 // Normalize our iterative median properly
00434 if (bUnit) VHS::UnitVector(median);
00435 else for (unsigned int j = 0; j < nDim; j++) median[j] *= 1./norm;
00436
00437 dItr = VHS::GetDistance( lastMedian, median );
00438 // Find the average distance to our new guess point
00439 for (unsigned int i = 0; i < nPts; i++) {
00440 double dPt = VHS::GetDistance(median, pt[i]);
00441 avgD += dPt/nPts;
00442 }
00443
00444 dDist = lastAvgD-avgD;
00445 nItr++;
00446 if (bVerbose)
00447 cout << "-- Iteration " << nItr << " : dItr = " << dItr
00448 << " : avgDist = " << avgD << " del(" << dDist
00449 << ")" << endl;
00450 } // end while
00451 if (bVerbose) {
00452 double dMedAvg = VHS::GetDistance( median, avg );
00453 cout << "|median-average| = " << dMedAvg << endl;
00454 }
00455
00456 return median;
00457
00458 }
|
|
||||||||||||||||||||
|
Definition at line 462 of file VHS.cxx. Referenced by VHS::FillDiscriminants(). 00466 {
00467 std::vector<double> full( 2*nPlanes*nStrips, 0.);
00468
00469 for (unsigned int i = 0; i < index.size(); i++)
00470 full[ index[i] ] = image[i];
00471
00472 return full;
00473
00474 }
|
|
||||||||||||
|
Definition at line 478 of file VHS.cxx. Referenced by VHS::FillDiscriminants(), and VHS::FindMedian(). 00479 {
00480 // Return the Euclidean distance between the points described
00481 // by the input vectors
00482
00483 if (vec0.size() != vec1.size()) return -1.;
00484
00485 double dist = 0.;
00486
00487 for (unsigned int i = 0; i < vec0.size(); i++)
00488 dist += TMath::Power( vec0[i]-vec1[i], 2. );
00489
00490 return TMath::Sqrt( dist );
00491
00492 }
|
|
||||||||||||
|
Definition at line 496 of file VHS.cxx. Referenced by VHS::DrawEvent(), VHS::Skim(), and VHS::Train(). 00497 {
00498 if ( iaction == 0 ) return vhsNC;
00499
00500 if ( iaction == 1 && TMath::Abs(inu) == 12 ) return vhsCCe;
00501
00502 if ( iaction == 1 && TMath::Abs(inu) == 14 ) return vhsCCmu;
00503
00504 if ( iaction == 1 && TMath::Abs(inu) == 16 ) return vhsCCtau;
00505
00506 return vhsUnknown;
00507
00508 }
|
|
||||||||||||||||||||||||||||||||||||
|
Definition at line 512 of file VHS.cxx. References VHS::GetThetaAxis(), NtpSRStrip::ph0, NtpSRStrip::ph1, NtpSRStrip::plane, VHS::RotatePixel(), VHS::SeparateViews(), NtpSRPulseHeight::sigcor, NtpSRStrip::strip, and VHS::VecIndex(). Referenced by VHS::FillDiscriminants(), and VHS::Train(). 00520 {
00521 // Return a standard vector of hit values representing
00522 // the event image for both U & V views
00523
00524 // The images will be centered on strip ~nStrips/2 & plane ~nPlanes/4
00525
00526 image.clear(); // empty the return vector of values
00527 vecInd.clear(); // empty the return vector of indices
00528 theta.clear(); // empty the return vector of event angle
00529 theta = std::vector<double>(2,-999.);
00530
00531 // Separate the input index vector & stp array into u & v
00532 std::vector< std::vector<int> > uvindex(2); //uvindex index 0=u, 1=v
00533 VHS::SeparateViews(index,nstp,stp,uvindex[0],uvindex[1]);
00534
00535 const int kOffsetInd = nPlanes*nStrips;
00536
00537 for (unsigned int uv = 0; uv <= 1 ; uv++) {
00538 std::vector<int> ind = uvindex[uv];
00539
00540 // Derive image rotation info
00541 std::vector<int> centerPix(2,0);
00542 VHS::GetThetaAxis(stp, uvindex[uv], theta[uv], centerPix);
00543
00544 // Find the minimum plane in the view
00545 int minPln = 999;
00546 for (unsigned int j = 0; j < ind.size(); j++) {
00547 NtpSRStrip* stpEntry = dynamic_cast<NtpSRStrip*>(stp->At(ind[j]));
00548 if (stpEntry) {
00549 int plane = stpEntry->plane;
00550 int strip = stpEntry->strip;
00551 VHS::RotatePixel(plane,strip,centerPix[0],centerPix[1],theta[uv]);
00552 if (plane < minPln) minPln = plane;
00553 } // end if stp ok
00554 } // end loop over strips
00555
00556 // Find the end plane in our window
00557 int endPln = minPln + 2*nPlanes;
00558
00559 // Find the weighted average plane in the view
00560 double stpXph = 0.;
00561 double totSigcor = 0.;
00562 for (unsigned int j = 0; j < ind.size(); j++) {
00563 NtpSRStrip* stpEntry = dynamic_cast<NtpSRStrip*>(stp->At(ind[j]));
00564 if (stpEntry) {
00565 int plane = stpEntry->plane;
00566 int strip = stpEntry->strip;
00567 VHS::RotatePixel(plane,strip,centerPix[0],centerPix[1],theta[uv]);
00568 if (plane >= minPln && plane < endPln) {
00569 double sigcor = stpEntry->ph0.sigcor + stpEntry->ph1.sigcor;
00570
00571 stpXph += plane*sigcor;
00572 totSigcor += sigcor;
00573 } // end if strip is in our plane window
00574 } // end if stp ok
00575 } // end loop over strips
00576
00577 // Find the weighted mean plane# in our window (~minos convention)
00578 // --> adjust the min & end plane numbers accordingly
00579 const int avgPln = TMath::FloorNint(stpXph/totSigcor);
00580 minPln = avgPln - TMath::FloorNint(0.25*2.*nPlanes);
00581 endPln = minPln + 2*nPlanes;
00582
00583 // Find the weighted average strip in the view
00584 stpXph = 0.;
00585 totSigcor = 0.;
00586 for (unsigned int j = 0; j < ind.size(); j++) {
00587 NtpSRStrip* stpEntry = dynamic_cast<NtpSRStrip*>(stp->At(ind[j]));
00588 if (stpEntry) {
00589 int plane = stpEntry->plane;
00590 int strip = stpEntry->strip;
00591 VHS::RotatePixel(plane,strip,centerPix[0],centerPix[1],theta[uv]);
00592 if (plane >= minPln && plane < endPln) {
00593 double sigcor = stpEntry->ph0.sigcor + stpEntry->ph1.sigcor;
00594
00595 stpXph += strip*sigcor;
00596 totSigcor += sigcor;
00597 } // end if strip is in our plane window
00598 } // end if stp ok
00599 } // end loop over strips
00600
00601 // Find the weighted mean strip in our window
00602 const int avgStp = TMath::FloorNint(stpXph/totSigcor);
00603 const int minStp = avgStp - nStrips/2;
00604
00605 // Loop over strips
00606 for (unsigned int i = 0; i < ind.size(); i++) {
00607 NtpSRStrip* stpEntry = dynamic_cast<NtpSRStrip*>(stp->At(ind[i]));
00608 if (stpEntry) {
00609 int plane = stpEntry->plane;
00610 int strip = stpEntry->strip;
00611 VHS::RotatePixel(plane,strip,centerPix[0],centerPix[1],theta[uv]);
00612 double sigcor = stpEntry->ph0.sigcor + stpEntry->ph1.sigcor;
00613
00614 int iPln = (plane - minPln)/2; // img coords
00615 int iStp = (strip - minStp); // img coords
00616 int iVec = VecIndex(iPln,iStp,nPlanes); // vec img coords
00617
00618 if (iPln >= 0 && iPln < nPlanes &&
00619 iStp >= 0 && iStp < nStrips ) {
00620 image.push_back(sigcor);
00621 if (uv == 0)
00622 vecInd.push_back(iVec);
00623 else
00624 vecInd.push_back(iVec+kOffsetInd);
00625 } // end if hit placement OK
00626 } // end if strip info OK
00627 } // end loop over strips
00628 } // end loop over u & v views
00629
00630 return;
00631
00632 }
|
|
||||||||||||||||||||
|
Definition at line 636 of file VHS.cxx. Referenced by VHS::FillDiscriminants(). 00640 {
00641 // Given a set of averages and variances,
00642 // calculate the Gaussian probability of the collection of
00643 // particular hit values convoluted with the probability of the
00644 // particular hit pattern
00645
00646 double ll = 0.;
00647 const double kCutoff = 1.e-10;
00648 const double kLnCutoff = TMath::Log10(kCutoff);
00649
00650 for (unsigned int i = 0; i < fullImage.size(); i++) {
00651 double x = fullImage[i];
00652 double mean = avg[i];
00653 double sigma = TMath::Sqrt( var[i] );
00654
00655 if ( x < kCutoff || pHit[i] < kCutoff )
00656 ll += TMath::Log10(1.-pHit[i]);
00657 else {
00658 ll += TMath::Log10(pHit[i]);
00659 ll += (TMath::Gaus(x,mean,sigma,true) > kCutoff)?
00660 TMath::Log10(TMath::Gaus(x,mean,sigma,true)) : kLnCutoff;
00661 ll -= (TMath::Gaus(mean,mean,sigma,true) > kCutoff)?
00662 TMath::Log10(TMath::Gaus(mean,mean,sigma,true)) : kLnCutoff;
00663 } // end if values within ranges
00664 } // end loop over image entries
00665
00666 return ll;
00667
00668 }
|
|
||||||||||||
|
Definition at line 672 of file VHS.cxx. 00673 {
00674 // Returns relative plane # (*not* MINOS global plane convention)
00675 return (nPlanes != 0)? vecInd % nPlanes : -1;
00676 }
|
|
||||||||||||
|
Definition at line 680 of file VHS.cxx. 00681 {
00682 // Returns relative strip # (*not* MINOS global strip convention)
00683 return (nPlanes != 0)? vecInd / nPlanes : -1;
00684 }
|
|
||||||||||||||||||||
|
Definition at line 688 of file VHS.cxx. References NtpSRStrip::ph0, NtpSRStrip::ph1, NtpSRStrip::plane, NtpSRPulseHeight::sigcor, and NtpSRStrip::strip. Referenced by VHS::GetImage(). 00692 {
00693 const unsigned int kNhits = index.size();
00694
00695 // Find the weighted average plane & strip in the view
00696 // And the minimum plane/strip
00697 double stpXph = 0.;
00698 double plnXph = 0.;
00699 double totSigcor = 0.;
00700 int minPln = 999;
00701 int minPlnStp = 999;
00702 for (unsigned int j = 0; j < kNhits; j++) {
00703 NtpSRStrip* stpEntry = dynamic_cast<NtpSRStrip*>(stp->At(index[j]));
00704 if (stpEntry) {
00705 int plane = stpEntry->plane;
00706 int strip = stpEntry->strip;
00707 double sigcor = stpEntry->ph0.sigcor + stpEntry->ph1.sigcor;
00708
00709 plnXph += plane*sigcor;
00710 stpXph += strip*sigcor;
00711 totSigcor += sigcor;
00712
00713 if (plane < minPln) {
00714 minPln = plane;
00715 minPlnStp = strip;
00716 } // end if minimum plane number
00717 } // end if stp ok
00718 } // end loop over strips
00719
00720 center[0] = (int)TMath::Floor(plnXph/totSigcor);
00721 center[1] = (int)TMath::Floor(stpXph/totSigcor);
00722
00723 double relPln = plnXph/totSigcor - 1.*minPln;
00724 double relStp = stpXph/totSigcor - 1.*minPlnStp;
00725
00726 theta = (TMath::Abs(relPln) > 1.e-10)?
00727 TMath::ATan( relStp/relPln ) : 0.;
00728
00729 return;
00730 }
|
|
|
Definition at line 734 of file VHS.cxx. 00735 {
00736 // Turn the input space delimited cstring into a vector of strings
00737
00738 std::string strList = cstr;
00739
00740 std::vector<std::string> strVec;
00741
00742 int len = strList.length();
00743
00744 unsigned int lastLoc = 0;
00745 unsigned int loc = strList.find(" ",lastLoc);
00746 while (lastLoc != std::string::npos) {
00747 int subLen = (loc != std::string::npos)? loc-lastLoc : len-lastLoc;
00748 std::string sub = strList.substr(lastLoc, subLen);
00749 if (sub.length() > 0) strVec.push_back(sub);
00750 lastLoc = (loc != std::string::npos)? loc+1 : loc;
00751 loc = strList.find(" ",lastLoc);
00752 }
00753
00754 return strVec;
00755 }
|
|
||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
|
Definition at line 759 of file VHS.cxx. References VHS::ToVector(). Referenced by ANtpInfoObjectFillerNC::InitializeVHSTraining(). 00773 {
00774 //
00775 // inFile must be open and valid
00776 //
00777 // Translate averages & completed variances from TH2D and write
00778 // to vectors for use
00779 //
00780
00781 // Prepare to write our training info to output file
00782 inFile->cd();
00783
00784 // NC hit prob hist
00785 const std::string UpNCName = "UpNCHit";
00786 const std::string VpNCName = "VpNCHit";
00787 TH2D* UpNCHist = (TH2D*)inFile->Get(UpNCName.c_str());
00788 TH2D* VpNCHist = (TH2D*)inFile->Get(VpNCName.c_str());
00789 std::vector<double> UpNChit = VHS::ToVector( UpNCHist, nPlanes);
00790 std::vector<double> VpNChit = VHS::ToVector( VpNCHist, nPlanes);
00791 pNChit.clear();
00792 for (unsigned int i = 0; i < UpNChit.size(); i++)
00793 pNChit.push_back(UpNChit[i]);
00794 for (unsigned int i = 0; i < VpNChit.size(); i++)
00795 pNChit.push_back(VpNChit[i]);
00796
00797 // CCe hit prob hist
00798 const std::string UpCCeName = "UpCCeHit";
00799 const std::string VpCCeName = "VpCCeHit";
00800 TH2D* UpCCeHist = (TH2D*)inFile->Get(UpCCeName.c_str());
00801 TH2D* VpCCeHist = (TH2D*)inFile->Get(VpCCeName.c_str());
00802 std::vector<double> UpCCehit = VHS::ToVector( UpCCeHist, nPlanes);
00803 std::vector<double> VpCCehit = VHS::ToVector( VpCCeHist, nPlanes);
00804 pCCehit.clear();
00805 for (unsigned int i = 0; i < UpCCehit.size(); i++)
00806 pCCehit.push_back(UpCCehit[i]);
00807 for (unsigned int i = 0; i < VpCCehit.size(); i++)
00808 pCCehit.push_back(VpCCehit[i]);
00809
00810 // CCmu hit prob hist
00811 const std::string UpCCmuName = "UpCCmuHit";
00812 const std::string VpCCmuName = "VpCCmuHit";
00813 TH2D* UpCCmuHist = (TH2D*)inFile->Get(UpCCmuName.c_str());
00814 TH2D* VpCCmuHist = (TH2D*)inFile->Get(VpCCmuName.c_str());
00815 std::vector<double> UpCCmuhit = VHS::ToVector( UpCCmuHist, nPlanes);
00816 std::vector<double> VpCCmuhit = VHS::ToVector( VpCCmuHist, nPlanes);
00817 pCCmuhit.clear();
00818 for (unsigned int i = 0; i < UpCCmuhit.size(); i++)
00819 pCCmuhit.push_back(UpCCmuhit[i]);
00820 for (unsigned int i = 0; i < VpCCmuhit.size(); i++)
00821 pCCmuhit.push_back(VpCCmuhit[i]);
00822
00823 // CCtau hit prob hist
00824 const std::string UpCCtauName = "UpCCtauHit";
00825 const std::string VpCCtauName = "VpCCtauHit";
00826 TH2D* UpCCtauHist = (TH2D*)inFile->Get(UpCCtauName.c_str());
00827 TH2D* VpCCtauHist = (TH2D*)inFile->Get(VpCCtauName.c_str());
00828 std::vector<double> UpCCtauhit = VHS::ToVector( UpCCtauHist, nPlanes);
00829 std::vector<double> VpCCtauhit = VHS::ToVector( VpCCtauHist, nPlanes);
00830 pCCtauhit.clear();
00831 for (unsigned int i = 0; i < UpCCtauhit.size(); i++)
00832 pCCtauhit.push_back(UpCCtauhit[i]);
00833 for (unsigned int i = 0; i < VpCCtauhit.size(); i++)
00834 pCCtauhit.push_back(VpCCtauhit[i]);
00835
00836 // NC Avg hist
00837 const std::string UavgNCName = "UAvgNCImg";
00838 const std::string VavgNCName = "VAvgNCImg";
00839 TH2D* UavgNCHist = (TH2D*)inFile->Get(UavgNCName.c_str());
00840 TH2D* VavgNCHist = (TH2D*)inFile->Get(VavgNCName.c_str());
00841 std::vector<double> UavgNC = VHS::ToVector( UavgNCHist, nPlanes);
00842 std::vector<double> VavgNC = VHS::ToVector( VavgNCHist, nPlanes);
00843 avgNC.clear();
00844 for (unsigned int i = 0; i < UavgNC.size(); i++)
00845 avgNC.push_back(UavgNC[i]);
00846 for (unsigned int i = 0; i < VavgNC.size(); i++)
00847 avgNC.push_back(VavgNC[i]);
00848
00849 // CCe Avg hist
00850 const std::string UavgCCeName = "UAvgCCeImg";
00851 const std::string VavgCCeName = "VAvgCCeImg";
00852 TH2D* UavgCCeHist = (TH2D*)inFile->Get(UavgCCeName.c_str());
00853 TH2D* VavgCCeHist = (TH2D*)inFile->Get(VavgCCeName.c_str());
00854 std::vector<double> UavgCCe = VHS::ToVector( UavgCCeHist, nPlanes);
00855 std::vector<double> VavgCCe = VHS::ToVector( VavgCCeHist, nPlanes);
00856 avgCCe.clear();
00857 for (unsigned int i = 0; i < UavgCCe.size(); i++)
00858 avgCCe.push_back(UavgCCe[i]);
00859 for (unsigned int i = 0; i < VavgCCe.size(); i++)
00860 avgCCe.push_back(VavgCCe[i]);
00861
00862 // CCmu Avg hist
00863 const std::string UavgCCmuName = "UAvgCCmuImg";
00864 const std::string VavgCCmuName = "VAvgCCmuImg";
00865 TH2D* UavgCCmuHist = (TH2D*)inFile->Get(UavgCCmuName.c_str());
00866 TH2D* VavgCCmuHist = (TH2D*)inFile->Get(VavgCCmuName.c_str());
00867 std::vector<double> UavgCCmu = VHS::ToVector( UavgCCmuHist, nPlanes);
00868 std::vector<double> VavgCCmu = VHS::ToVector( VavgCCmuHist, nPlanes);
00869 avgCCmu.clear();
00870 for (unsigned int i = 0; i < UavgCCmu.size(); i++)
00871 avgCCmu.push_back(UavgCCmu[i]);
00872 for (unsigned int i = 0; i < VavgCCmu.size(); i++)
00873 avgCCmu.push_back(VavgCCmu[i]);
00874
00875 // CCtau Avg hist
00876 const std::string UavgCCtauName = "UAvgCCtauImg";
00877 const std::string VavgCCtauName = "VAvgCCtauImg";
00878 TH2D* UavgCCtauHist = (TH2D*)inFile->Get(UavgCCtauName.c_str());
00879 TH2D* VavgCCtauHist = (TH2D*)inFile->Get(VavgCCtauName.c_str());
00880 std::vector<double> UavgCCtau = VHS::ToVector( UavgCCtauHist, nPlanes);
00881 std::vector<double> VavgCCtau = VHS::ToVector( VavgCCtauHist, nPlanes);
00882 avgCCtau.clear();
00883 for (unsigned int i = 0; i < UavgCCtau.size(); i++)
00884 avgCCtau.push_back(UavgCCtau[i]);
00885 for (unsigned int i = 0; i < VavgCCtau.size(); i++)
00886 avgCCtau.push_back(VavgCCtau[i]);
00887
00888 // NC variance hist
00889 const std::string UNCVarName = "UNCVar";
00890 const std::string VNCVarName = "VNCVar";
00891 TH2D* UNCVarHist = (TH2D*)inFile->Get(UNCVarName.c_str());
00892 TH2D* VNCVarHist = (TH2D*)inFile->Get(VNCVarName.c_str());
00893 std::vector<double> UNCVar = VHS::ToVector( UNCVarHist, nPlanes);
00894 std::vector<double> VNCVar = VHS::ToVector( VNCVarHist, nPlanes);
00895 varNC.clear();
00896 for (unsigned int i = 0; i < UNCVar.size(); i++)
00897 varNC.push_back(UNCVar[i]);
00898 for (unsigned int i = 0; i < VpNChit.size(); i++)
00899 varNC.push_back(VNCVar[i]);
00900
00901 // CCe variance hist
00902 const std::string UCCeVarName = "UCCeVar";
00903 const std::string VCCeVarName = "VCCeVar";
00904 TH2D* UCCeVarHist = (TH2D*)inFile->Get(UCCeVarName.c_str());
00905 TH2D* VCCeVarHist = (TH2D*)inFile->Get(VCCeVarName.c_str());
00906 std::vector<double> UCCeVar = VHS::ToVector( UCCeVarHist, nPlanes);
00907 std::vector<double> VCCeVar = VHS::ToVector( VCCeVarHist, nPlanes);
00908 varCCe.clear();
00909 for (unsigned int i = 0; i < UCCeVar.size(); i++)
00910 varCCe.push_back(UCCeVar[i]);
00911 for (unsigned int i = 0; i < VpCCehit.size(); i++)
00912 varCCe.push_back(VCCeVar[i]);
00913
00914 // CCmu variance hist
00915 const std::string UCCmuVarName = "UCCmuVar";
00916 const std::string VCCmuVarName = "VCCmuVar";
00917 TH2D* UCCmuVarHist = (TH2D*)inFile->Get(UCCmuVarName.c_str());
00918 TH2D* VCCmuVarHist = (TH2D*)inFile->Get(VCCmuVarName.c_str());
00919 std::vector<double> UCCmuVar = VHS::ToVector( UCCmuVarHist, nPlanes);
00920 std::vector<double> VCCmuVar = VHS::ToVector( VCCmuVarHist, nPlanes);
00921 varCCmu.clear();
00922 for (unsigned int i = 0; i < UCCmuVar.size(); i++)
00923 varCCmu.push_back(UCCmuVar[i]);
00924 for (unsigned int i = 0; i < VpCCmuhit.size(); i++)
00925 varCCmu.push_back(VCCmuVar[i]);
00926
00927 // CCtau variance hist
00928 const std::string UCCtauVarName = "UCCtauVar";
00929 const std::string VCCtauVarName = "VCCtauVar";
00930 TH2D* UCCtauVarHist = (TH2D*)inFile->Get(UCCtauVarName.c_str());
00931 TH2D* VCCtauVarHist = (TH2D*)inFile->Get(VCCtauVarName.c_str());
00932 std::vector<double> UCCtauVar = VHS::ToVector( UCCtauVarHist, nPlanes);
00933 std::vector<double> VCCtauVar = VHS::ToVector( VCCtauVarHist, nPlanes);
00934 varCCtau.clear();
00935 for (unsigned int i = 0; i < UCCtauVar.size(); i++)
00936 varCCtau.push_back(UCCtauVar[i]);
00937 for (unsigned int i = 0; i < VpCCtauhit.size(); i++)
00938 varCCtau.push_back(VCCtauVar[i]);
00939
00940 inFile->Close();
00941
00942 return;
00943
00944 }
|
|
||||||||||||||||||||||||
|
Definition at line 948 of file VHS.cxx. Referenced by VHS::GetImage(). 00951 {
00952 // rotate around pixel center to new coords (pln, stp);
00953 // --> round floats to nearest integer
00954
00955 if ( TMath::Abs(theta) < 1.e-10 ) return; // ~zero rotation
00956
00957 const double cosTheta = TMath::Cos(-theta); // rotate opposite dir as given
00958 const double sinTheta = TMath::Sin(-theta); // rotate opposite dir as given
00959
00960 double newPlnD = cosTheta*(plane-avgPlane) + sinTheta*(strip-avgStrip);
00961 double newStpD = cosTheta*(strip-avgStrip) - sinTheta*(plane-avgPlane);
00962
00963 if (newStpD-TMath::Floor(newStpD) < 0.5) newStpD = TMath::Floor(newStpD);
00964 else newStpD = TMath::Ceil(newStpD);
00965
00966 if (newPlnD-TMath::Floor(newPlnD) < 0.5) newPlnD = TMath::Floor(newPlnD);
00967 else newPlnD = TMath::Ceil(newPlnD);
00968
00969 strip = avgStrip + (int)newStpD;
00970 plane = avgPlane + (int)newPlnD;
00971
00972 return;
00973
00974 }
|
|
||||||||||||||||||||||||
|
Definition at line 978 of file VHS.cxx. References NtpSRStrip::planeview. Referenced by VHS::GetImage(). 00980 {
00981 // Sort hit strips into separate views
00982
00983 for (int i = 0; i < nstp; i++) {
00984 int ind = index[i];
00985 NtpSRStrip* strip = dynamic_cast<NtpSRStrip*>(stp->At(ind));
00986 if (strip) {
00987 unsigned int pv = strip->planeview;
00988 if (pv == PlaneView::kU) ustp.push_back(ind);
00989 if (pv == PlaneView::kV) vstp.push_back(ind);
00990 }
00991 }
00992 return;
00993 }
|
|
||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
|
Definition at line 997 of file VHS.cxx. References NtpSRPlane::beg, NtpSRPlane::end, NtpStRecord::evt, NtpStRecord::evthdr, VHS::FillDiscriminants(), VHS::GetEvtType(), RecRecordImp< T >::GetHeader(), NtpSRStripPulseHeight::gev, VHSevent::GeV, NtpMCTruth::iaction, NtpMCTruth::inu, NtpMCTruth::iresonance, VHSevent::MC, NtpStRecord::mc, VHSevent::MCenu, VHSevent::MCiaction, VHSevent::MCinu, VHSevent::MCQsq, VHSevent::MCrescode, VHSevent::MCWsq, VHSevent::MCx, VHSevent::MCy, NtpTHEvent::neumc, NtpSREventSummary::nevent, NtpSREvent::nstrip, VHSevent::NumPln, VHSevent::NumStp, NtpMCTruth::p4neu, NtpSREvent::ph, NtpSREvent::plane, NtpMCTruth::q2, NtpSRPulseHeight::sigcor, VHSevent::Sigcor, NtpStRecord::stp, NtpStRecord::thevt, NtpMCTruth::w2, NtpMCTruth::x, and NtpMCTruth::y. 01013 {
01014 // Inputs: NtpStRecord, average images, variances,
01015 // vector of raw hit probabilities, window size (pln x stp),
01016 // flag signalling whether to compute distances on unit or
01017 // absolute image vectors
01018 // Output: Vector of pointers to complete VHSevent records
01019 //
01020
01021 // Define the vector we will return
01022 vector< VHSevent* > vhslist;
01023
01024 // If we have no record, do nothing
01025 if ( !ntpst ) return vhslist;
01026
01027 const NtpSREventSummary& evthdr = ntpst->evthdr;
01028 TClonesArray* evt = ntpst->evt;
01029 TClonesArray* stp = ntpst->stp;
01030 TClonesArray* mc = ntpst->mc;
01031 TClonesArray* thevt = ntpst->thevt;
01032
01033 // Event quantities
01034 if (evt && stp) {
01035 const unsigned int nEvt = evthdr.nevent; // evt->GetEntries();
01036
01037 for (unsigned int i = 0; i < nEvt; i++){
01038 NtpSREvent* evtEntry = dynamic_cast<NtpSREvent*>( evt->At(i) );
01039 if ( !evtEntry ) continue;
01040
01041 // Define our VHSevent to fill
01042 VHSevent* vhsevt = new VHSevent( ntpst->GetHeader() ) ;
01043 vhsevt->SetName("VHSevent");
01044 vhsevt->SetTitle("VHS Event");
01045
01046 // Fill NtpSREvent info
01047 vhsevt->GeV = evtEntry->ph.gev;
01048 vhsevt->NumPln = evtEntry->plane.end - evtEntry->plane.beg;
01049 vhsevt->NumStp = evtEntry->nstrip;
01050 vhsevt->Sigcor = evtEntry->ph.sigcor;
01051
01052 VHS::evtType evtt = vhsUnknown;
01053
01054 // Fill MC information
01055 if ( mc && thevt ) {
01056 NtpTHEvent* thevtEntry = dynamic_cast<NtpTHEvent*>(thevt->At(i));
01057 NtpMCTruth* mcEntry = (thevtEntry)?
01058 dynamic_cast<NtpMCTruth*>( mc->At( (thevtEntry->neumc) ) ) : 0;
01059 if (mcEntry) {
01060 vhsevt->MCenu = mcEntry->p4neu[3];
01061 vhsevt->MCinu = mcEntry->inu;
01062 vhsevt->MCiaction = mcEntry->iaction;
01063 vhsevt->MCrescode = mcEntry->iresonance;
01064 vhsevt->MCQsq = mcEntry->q2;
01065 vhsevt->MCWsq = mcEntry->w2;
01066 vhsevt->MCx = mcEntry->x;
01067 vhsevt->MCy = mcEntry->y;
01068
01069 evtt = VHS::GetEvtType(mcEntry->inu, mcEntry->iaction);
01070 } // end if mcEntry
01071 } // end if mc && thevt
01072 vhsevt->MC = (evtt != vhsUnknown);
01073
01074 VHS::FillDiscriminants(evtEntry, stp,
01075 avgNC, avgCCe, avgCCmu, avgCCtau,
01076 varNC, varCCe, varCCmu, varCCtau,
01077 pNC, pCCe, pCCmu, pCCtau,
01078 nPlanes, nStrips, bUnit, vhsevt );
01079
01080 vhslist.push_back( vhsevt ); // add this event to the list
01081
01082 } // end loop over event entries
01083 } // end if evt and stp arrays exist
01084
01085 return vhslist;
01086 }
|
|
||||||||||||
|
Definition at line 1090 of file VHS.cxx. 01092 {
01093 // 'nuff said
01094
01095 std::vector<double> subt(image0.size(), 0.);
01096
01097 if (image0.size() != image1.size() ) return subt;
01098
01099 for (unsigned int i = 0; i < image0.size(); i++)
01100 subt[i] = image0[i] - image1[i];
01101
01102 return subt;
01103
01104 }
|
|
||||||||||||||||||||
|
Definition at line 1109 of file VHS.cxx. Referenced by VHS::DrawEvent(). 01110 {
01111
01112 std::ostringstream os;
01113 os << std::setw(w) << std::setprecision(p) << thing;
01114 return os.str();
01115
01116 }
|
|
||||||||||||||||||||||||||||||||
|
Definition at line 1120 of file VHS.cxx. Referenced by VHS::DrawEvent(), and VHS::WriteFile(). 01124 {
01125 // 'nuff said
01126
01127 std::string Uname = "U" + name;
01128 std::string Utitle = "U" + title;
01129
01130 std::string Vname = "V" + name;
01131 std::string Vtitle = "V" + title;
01132
01133 //if (Uview) delete Uview;
01134 //if (Vview) delete Vview;
01135
01136 Uview = new TH2D(Uname.c_str(), Utitle.c_str(),
01137 nPlanes, 0, nPlanes,
01138 nStrips, 0, nStrips);
01139 Vview = new TH2D(Vname.c_str(), Vtitle.c_str(),
01140 nPlanes, 0, nPlanes,
01141 nStrips, 0, nStrips);
01142 const unsigned int kOffset = nPlanes*nStrips;
01143 for (unsigned int i = 0; i < vec.size(); i++){
01144 if ( i < kOffset )
01145 Uview->SetBinContent( VHS::GetPlane(i,nPlanes)+1,
01146 VHS::GetStrip(i,nPlanes)+1,
01147 vec[i] );
01148 else
01149 Vview->SetBinContent( VHS::GetPlane(i-kOffset,nPlanes)+1,
01150 VHS::GetStrip(i-kOffset,nPlanes)+1,
01151 vec[i] );
01152 } // end loop over vec entries
01153
01154 return;
01155
01156 }
|
|
||||||||||||
|
Definition at line 1160 of file VHS.cxx. References VHS::VecIndex(). Referenced by VHS::ReadFile(). 01161 {
01162 // 'nuff said
01163
01164 const int nPln = hist->GetNbinsX();
01165 const int nStp = hist->GetNbinsY();
01166 const int nBins = nPln*nStp;
01167
01168 std::vector<double> vec(nBins,0.);
01169
01170 if (nPln != nPlanes) return vec;
01171
01172 for (int ipln = 0; ipln < nPln; ipln++)
01173 for (int istp = 0; istp < nStp; istp++)
01174 vec[VecIndex(ipln,istp,nPlanes)] =
01175 hist->GetBinContent(ipln+1, istp+1);
01176
01177 return vec;
01178
01179 }
|
|
||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
|
Definition at line 1183 of file VHS.cxx. References NtpSRPlane::beg, NtpSRPlane::end, NtpStRecord::evt, NtpStRecord::evthdr, VHS::GetEvtType(), VHS::GetImage(), NtpSRStripPulseHeight::gev, NtpMCTruth::iaction, NtpMCTruth::inu, NtpStRecord::mc, NtpTHEvent::neumc, NtpSREventSummary::nevent, NtpSREvent::nstrip, NtpMCTruth::p4neu, NtpSREvent::ph, NtpSREvent::plane, NtpSREvent::stp, NtpStRecord::stp, NtpStRecord::thevt, VHS::vhsCCe, VHS::vhsCCmu, VHS::vhsCCtau, VHS::vhsNC, and VHS::vhsUnknown. 01206 {
01207 // Input : NtpStRecord, average images to modify, variances to
01208 // modify, vectors of number of hits to each pixel to
01209 // modify, #planes in each view window, #strips in each view
01210 // window
01211 // Output: Modified average images & modified covariance matrices
01212 //
01213
01214 int nEvtUsed = 0;
01215
01216 // If we have no record, do nothing
01217 if ( !ntpst ) return nEvtUsed;
01218
01219 const NtpSREventSummary& evthdr = ntpst->evthdr;
01220 TClonesArray* evt = ntpst->evt;
01221 TClonesArray* stp = ntpst->stp;
01222 TClonesArray* mc = ntpst->mc;
01223 TClonesArray* thevt = ntpst->thevt;
01224
01225 // Event quantities
01226 if (evt && stp && mc && thevt) {
01227 const unsigned int nEvt = evthdr.nevent; // evt->GetEntries();
01228
01229 for (unsigned int i = 0; i < nEvt; i++){
01230 NtpSREvent* evtEntry = dynamic_cast<NtpSREvent*>( evt->At(i) );
01231 if ( !evtEntry ) continue;
01232
01233 NtpTHEvent* thevtEntry = dynamic_cast<NtpTHEvent*>(thevt->At(i));
01234 NtpMCTruth* mcEntry = (thevtEntry)?
01235 dynamic_cast<NtpMCTruth*>( mc->At( (thevtEntry->neumc) ) ) : 0;
01236 if ( !mcEntry ) continue;
01237
01238 const double eReco = evtEntry->ph.gev;
01239 const double eTrue = mcEntry->p4neu[3];
01240
01241 if (eReco < eRecoMax || eTrue < eTrueMax ) continue;
01242
01243 VHS::evtType evtt = VHS::GetEvtType(mcEntry->inu, mcEntry->iaction);
01244
01245 if ( evtt == vhsUnknown ) continue;
01246
01247 if (evtEntry->plane.end - evtEntry->plane.beg > cutPlanes) continue;
01248
01249 if ( ( evtt == vhsNC && NCevts >= maxTrain ) ||
01250 ( evtt == vhsCCe && eCCevts >= maxTrain ) ||
01251 ( evtt == vhsCCmu && muCCevts >= maxTrain ) ||
01252 ( evtt == vhsCCtau && tauCCevts >= maxTrain ) )
01253 continue;
01254
01255 // Gather our image and index sparse vectors
01256 std::vector<double> image;
01257 std::vector<int> index;
01258 std::vector<double> theta;
01259 VHS::GetImage(evtEntry->stp, evtEntry->nstrip,
01260 stp, nPlanes, nStrips, image, index, theta);
01261
01262 // Break up into cases by event type
01263 switch (evtt) {
01264 case vhsNC :
01265 for (unsigned int ind = 0; ind < image.size(); ind++) {
01266 int ii = index[ind];
01267 double xi = image[ind];
01268 int ni = numNC[ii];
01269
01270 avgNC[ii] = (avgNC[ii]*ni + xi )/(ni+1);
01271 varNC[ii] = (varNC[ii]*ni + xi*xi)/(ni+1);
01272 numNC[ii]++;
01273 } // end outer for loop
01274 NCevts++;
01275
01276 break;
01277 case vhsCCe :
01278 for (unsigned int ind = 0; ind < image.size(); ind++) {
01279 int ii = index[ind];
01280 double xi = image[ind];
01281 int ni = numCCe[ii];
01282
01283 avgCCe[ii] = (avgCCe[ii]*ni + xi )/(ni+1);
01284 varCCe[ii] = (varCCe[ii]*ni + xi*xi)/(ni + 1);
01285 numCCe[ii]++;
01286 } // end outer for loop
01287 eCCevts++;
01288
01289 break;
01290 case vhsCCmu :
01291 for (unsigned int ind = 0; ind < image.size(); ind++) {
01292 int ii = index[ind];
01293 double xi = image[ind];
01294 int ni = numCCmu[ii];
01295
01296 avgCCmu[ii] = (avgCCmu[ii]*ni + xi )/(ni+1);
01297 varCCmu[ii] = (varCCmu[ii]*ni + xi*xi)/(ni+1);
01298 numCCmu[ii]++;
01299 } // end outer for loop
01300 muCCevts++;
01301
01302 break;
01303 case vhsCCtau :
01304 for (unsigned int ind = 0; ind < image.size(); ind++) {
01305 int ii = index[ind];
01306 double xi = image[ind];
01307 int ni = numCCtau[ii];
01308
01309 avgCCtau[ii] = (avgCCtau[ii]*ni + xi )/(ni+1);
01310 varCCtau[ii] = (varCCtau[ii]*ni + xi*xi)/(ni+1);
01311 numCCtau[ii]++;
01312 } // end outer for loop
01313 tauCCevts++;
01314
01315 break;
01316 case vhsUnknown :
01317 continue;
01318 } // end switch (evtType)
01319
01320 nEvtUsed++;
01321
01322 } // end for loop over events
01323 } // end if the ntpst arrays are valid
01324
01325 return nEvtUsed;
01326
01327 }
|
|
||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
|
Definition at line 1331 of file VHS.cxx. References VHS::WriteFile(). 01350 {
01351 //
01352 // Post-processing of the variances, assuming input "variances"
01353 // are entries simply of <x^2>
01354 //
01355 // Translate averages & completed variances into TH2D and write
01356 // to file, outFile
01357 //
01358 // Call this only once at the end of the total training job
01359 //
01360
01361 const unsigned int kImgSz = avgNC.size();
01362
01363 // Calculate the simple raw hit probabilities
01364 std::vector<double> pNChit (kImgSz, 0.0);
01365 std::vector<double> pCCehit (kImgSz, 0.0);
01366 std::vector<double> pCCmuhit (kImgSz, 0.0);
01367 std::vector<double> pCCtauhit(kImgSz, 0.0);
01368
01369 for (unsigned int pix = 0; pix < kImgSz; pix++) {
01370 pNChit[pix] = (NCevts > 0)?
01371 (double)numNC[pix] / (double)NCevts : 0.;
01372 pCCehit[pix] = (eCCevts > 0)?
01373 (double)numCCe[pix] / (double)eCCevts : 0.;
01374 pCCmuhit[pix] = (muCCevts > 0)?
01375 (double)numCCmu[pix] / (double)muCCevts : 0.;
01376 pCCtauhit[pix] = (tauCCevts > 0)?
01377 (double)numCCtau[pix] / (double)tauCCevts : 0.;
01378 }
01379
01380 // Complete the variance calculation: var = <x^2> - <x>^2
01381 for (unsigned int ii = 0; ii < kImgSz; ii++) {
01382 varNC[ii] -= avgNC[ii] * avgNC[ii];
01383 varCCe[ii] -= avgCCe[ii] * avgCCe[ii];
01384 varCCmu[ii] -= avgCCmu[ii] * avgCCmu[ii];
01385 varCCtau[ii] -= avgCCtau[ii] * avgCCtau[ii];
01386
01387 } // end strip loop
01388
01389 VHS::WriteFile( outFile,
01390 avgNC, avgCCe, avgCCmu, avgCCtau,
01391 varNC, varCCe, varCCmu, varCCtau,
01392 pNChit, pCCehit, pCCmuhit, pCCtauhit,
01393 nPlanes, nStrips);
01394
01395 return;
01396
01397 }
|
|
|
Definition at line 1401 of file VHS.cxx. References len. Referenced by VHS::FillDiscriminants(), and VHS::FindMedian(). 01402 {
01403 double len = 0.;
01404
01405 for (unsigned int i = 0; i < vec.size(); i++)
01406 len += TMath::Power( vec[i], 2. );
01407
01408 len = TMath::Sqrt( len );
01409
01410 if (len <= 0.) return;
01411
01412 for (unsigned int i = 0; i < vec.size(); i++)
01413 vec[i] = vec[i] / len ;
01414
01415 return;
01416
01417 }
|
|
||||||||||||||||
|
Definition at line 1421 of file VHS.cxx. Referenced by VHS::GetImage(), and VHS::ToVector(). 01422 {
01423 // A particular choice of convention
01424
01425 return (ipln + nPlanes * istp);
01426
01427 }
|
|
||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
|
Definition at line 1431 of file VHS.cxx. References VHS::ToTH2D(). Referenced by VHS::TrainPost(). 01446 {
01447 //
01448 // outFile must be open and valid
01449 //
01450 // Translate averages & completed variances into TH2D and write
01451 // to file, outFile
01452 //
01453 // Call this only once at the end of the total training job
01454 //
01455
01456 // Prepare to write our training info to output file
01457 outFile->cd();
01458
01459 const std::string imgXTitle = "Plane";
01460 const std::string imgYTitle = "Strip";
01461
01462 // NC hit prob hist
01463 const std::string pNCTitle = "NC Hit Probabilities";
01464 const std::string pNCName = "pNCHit";
01465 TH2D* UpNCHist; TH2D* VpNCHist;
01466 VHS::ToTH2D( pNChit,
01467 nPlanes, nStrips,
01468 pNCName.c_str(),
01469 pNCTitle.c_str(),
01470 UpNCHist, VpNCHist);
01471 UpNCHist->SetXTitle( imgXTitle.c_str() );
01472 UpNCHist->SetYTitle( imgYTitle.c_str() );
01473 UpNCHist->Write();
01474 VpNCHist->SetXTitle( imgXTitle.c_str() );
01475 VpNCHist->SetYTitle( imgYTitle.c_str() );
01476 VpNCHist->Write();
01477
01478 // CCe hit prob hist
01479 const std::string pCCeTitle = "CCe Hit Probabilities";
01480 const std::string pCCeName = "pCCeHit";
01481 TH2D* UpCCeHist; TH2D* VpCCeHist;
01482 VHS::ToTH2D( pCCehit,
01483 nPlanes, nStrips,
01484 pCCeName.c_str(),
01485 pCCeTitle.c_str(),
01486 UpCCeHist, VpCCeHist);
01487 UpCCeHist->SetXTitle( imgXTitle.c_str() );
01488 UpCCeHist->SetYTitle( imgYTitle.c_str() );
01489 UpCCeHist->Write();
01490 VpCCeHist->SetXTitle( imgXTitle.c_str() );
01491 VpCCeHist->SetYTitle( imgYTitle.c_str() );
01492 VpCCeHist->Write();
01493
01494 // CCmu hit prob hist
01495 const std::string pCCmuTitle = "CC#mu Hit Probabilities";
01496 const std::string pCCmuName = "pCCmuHit";
01497 TH2D* UpCCmuHist; TH2D* VpCCmuHist;
01498 VHS::ToTH2D( pCCmuhit,
01499 nPlanes, nStrips,
01500 pCCmuName.c_str(),
01501 pCCmuTitle.c_str(),
01502 UpCCmuHist, VpCCmuHist);
01503 UpCCmuHist->SetXTitle( imgXTitle.c_str() );
01504 UpCCmuHist->SetYTitle( imgYTitle.c_str() );
01505 UpCCmuHist->Write();
01506 VpCCmuHist->SetXTitle( imgXTitle.c_str() );
01507 VpCCmuHist->SetYTitle( imgYTitle.c_str() );
01508 VpCCmuHist->Write();
01509
01510 // CCtau hit prob hist
01511 const std::string pCCtauTitle = "CC#tau Hit Probabilities";
01512 const std::string pCCtauName = "pCCtauHit";
01513 TH2D* UpCCtauHist; TH2D* VpCCtauHist;
01514 VHS::ToTH2D( pCCtauhit,
01515 nPlanes, nStrips,
01516 pCCtauName.c_str(),
01517 pCCtauTitle.c_str(),
01518 UpCCtauHist, VpCCtauHist);
01519 UpCCtauHist->SetXTitle( imgXTitle.c_str() );
01520 UpCCtauHist->SetYTitle( imgYTitle.c_str() );
01521 UpCCtauHist->Write();
01522 VpCCtauHist->SetXTitle( imgXTitle.c_str() );
01523 VpCCtauHist->SetYTitle( imgYTitle.c_str() );
01524 VpCCtauHist->Write();
01525
01526 // NC Avg hist
01527 const std::string avgNCTitle = "Average NC Image";
01528 const std::string avgNCName = "AvgNCImg";
01529 TH2D* UavgNCHist; TH2D* VavgNCHist;
01530 VHS::ToTH2D( avgNC,
01531 nPlanes, nStrips,
01532 avgNCName.c_str(),
01533 avgNCTitle.c_str(),
01534 UavgNCHist, VavgNCHist);
01535 UavgNCHist->SetXTitle( imgXTitle.c_str() );
01536 UavgNCHist->SetYTitle( imgYTitle.c_str() );
01537 UavgNCHist->Write();
01538 VavgNCHist->SetXTitle( imgXTitle.c_str() );
01539 VavgNCHist->SetYTitle( imgYTitle.c_str() );
01540 VavgNCHist->Write();
01541
01542 // CCe Avg hist
01543 const std::string avgCCeTitle = "Average CCe Image";
01544 const std::string avgCCeName = "AvgCCeImg";
01545 TH2D* UavgCCeHist; TH2D* VavgCCeHist;
01546 VHS::ToTH2D( avgCCe,
01547 nPlanes, nStrips,
01548 avgCCeName.c_str(),
01549 avgCCeTitle.c_str(),
01550 UavgCCeHist, VavgCCeHist);
01551 UavgCCeHist->SetXTitle( imgXTitle.c_str() );
01552 UavgCCeHist->SetYTitle( imgYTitle.c_str() );
01553 UavgCCeHist->Write();
01554 VavgCCeHist->SetXTitle( imgXTitle.c_str() );
01555 VavgCCeHist->SetYTitle( imgYTitle.c_str() );
01556 VavgCCeHist->Write();
01557
01558 // CCmu Avg hist
01559 const std::string avgCCmuTitle = "Average CC#mu Image";
01560 const std::string avgCCmuName = "AvgCCmuImg";
01561 TH2D* UavgCCmuHist; TH2D* VavgCCmuHist;
01562 VHS::ToTH2D( avgCCmu,
01563 nPlanes, nStrips,
01564 avgCCmuName.c_str(),
01565 avgCCmuTitle.c_str(),
01566 UavgCCmuHist, VavgCCmuHist);
01567 UavgCCmuHist->SetXTitle( imgXTitle.c_str() );
01568 UavgCCmuHist->SetYTitle( imgYTitle.c_str() );
01569 UavgCCmuHist->Write();
01570 VavgCCmuHist->SetXTitle( imgXTitle.c_str() );
01571 VavgCCmuHist->SetYTitle( imgYTitle.c_str() );
01572 VavgCCmuHist->Write();
01573
01574 // CCtau Avg hist
01575 const std::string avgCCtauTitle = "Average CC#tau Image";
01576 const std::string avgCCtauName = "AvgCCtauImg";
01577 TH2D* UavgCCtauHist; TH2D* VavgCCtauHist;
01578 VHS::ToTH2D( avgCCtau,
01579 nPlanes, nStrips,
01580 avgCCtauName.c_str(),
01581 avgCCtauTitle.c_str(),
01582 UavgCCtauHist, VavgCCtauHist);
01583 UavgCCtauHist->SetXTitle( imgXTitle.c_str() );
01584 UavgCCtauHist->SetYTitle( imgYTitle.c_str() );
01585 UavgCCtauHist->Write();
01586 VavgCCtauHist->SetXTitle( imgXTitle.c_str() );
01587 VavgCCtauHist->SetYTitle( imgYTitle.c_str() );
01588 VavgCCtauHist->Write();
01589
01590 // NC variance hist
01591 const std::string varNCTitle = "NC Variance";
01592 const std::string varNCName = "NCVar";
01593 TH2D* UvarNCHist; TH2D* VvarNCHist;
01594 VHS::ToTH2D( varNC,
01595 nPlanes, nStrips,
01596 varNCName.c_str(),
01597 varNCTitle.c_str(),
01598 UvarNCHist, VvarNCHist);
01599 UvarNCHist->SetXTitle( imgXTitle.c_str() );
01600 UvarNCHist->SetYTitle( imgYTitle.c_str() );
01601 UvarNCHist->Write();
01602 VvarNCHist->SetXTitle( imgXTitle.c_str() );
01603 VvarNCHist->SetYTitle( imgYTitle.c_str() );
01604 VvarNCHist->Write();
01605
01606 // CCe variance hist
01607 const std::string varCCeTitle = "CCe Variance";
01608 const std::string varCCeName = "CCeVar";
01609 TH2D* UvarCCeHist; TH2D* VvarCCeHist;
01610 VHS::ToTH2D( varCCe,
01611 nPlanes, nStrips,
01612 varCCeName.c_str(),
01613 varCCeTitle.c_str(),
01614 UvarCCeHist, VvarCCeHist);
01615 UvarCCeHist->SetXTitle( imgXTitle.c_str() );
01616 UvarCCeHist->SetYTitle( imgYTitle.c_str() );
01617 UvarCCeHist->Write();
01618 VvarCCeHist->SetXTitle( imgXTitle.c_str() );
01619 VvarCCeHist->SetYTitle( imgYTitle.c_str() );
01620 VvarCCeHist->Write();
01621
01622 // CCmu variance hist
01623 const std::string varCCmuTitle = "CCmu Variance";
01624 const std::string varCCmuName = "CCmuVar";
01625 TH2D* UvarCCmuHist; TH2D* VvarCCmuHist;
01626 VHS::ToTH2D( varCCmu,
01627 nPlanes, nStrips,
01628 varCCmuName.c_str(),
01629 varCCmuTitle.c_str(),
01630 UvarCCmuHist, VvarCCmuHist);
01631 UvarCCmuHist->SetXTitle( imgXTitle.c_str() );
01632 UvarCCmuHist->SetYTitle( imgYTitle.c_str() );
01633 UvarCCmuHist->Write();
01634 VvarCCmuHist->SetXTitle( imgXTitle.c_str() );
01635 VvarCCmuHist->SetYTitle( imgYTitle.c_str() );
01636 VvarCCmuHist->Write();
01637
01638 // CCtau variance hist
01639 const std::string varCCtauTitle = "CCtau Variance";
01640 const std::string varCCtauName = "CCtauVar";
01641 TH2D* UvarCCtauHist; TH2D* VvarCCtauHist;
01642 VHS::ToTH2D( varCCtau,
01643 nPlanes, nStrips,
01644 varCCtauName.c_str(),
01645 varCCtauTitle.c_str(),
01646 UvarCCtauHist, VvarCCtauHist);
01647 UvarCCtauHist->SetXTitle( imgXTitle.c_str() );
01648 UvarCCtauHist->SetYTitle( imgYTitle.c_str() );
01649 UvarCCtauHist->Write();
01650 VvarCCtauHist->SetXTitle( imgXTitle.c_str() );
01651 VvarCCtauHist->SetYTitle( imgYTitle.c_str() );
01652 VvarCCtauHist->Write();
01653
01654 return;
01655
01656 }
|
1.3.9.1