Newer
Older
TB_Chris / TbUT / scripts / .svn / text-base / plotADCs.C.svn-base
  1. ///////////////////////////////////////////////////
  2. //
  3. // Macro: plotADCs.C
  4. // Author: S. Blusk
  5. // Date: Oct 19, 2015
  6. //
  7. // Purpose: To loop over events and look at the CMS ADC distribution vs channel number
  8. // Can also be run on selected events, liste in a file: "MissingDUTHits.dat", which can
  9. // be obtained by running ClusterAnaWithTrack
  10. // with the flag
  11. // writeEventsWithMissinhHitsToFile = true
  12. //
  13. // Note: The cursor needs to be on the event display window to advance to the next event.
  14. // Hit ENTER to go to next event.
  15. //
  16. ///////////////////////////////////////////////////
  17.  
  18. bool plotSelectedEvents = true;
  19. bool showEventDisplay = true;
  20.  
  21. bool fillStatDistribution = true; // Fill summary histogram
  22.  
  23. int lowCh = 127; // Low Channel number of ADC distribution
  24. int hiCh = 256; // High Channel number of ADC distribution
  25.  
  26. int evList[20000];
  27. double nomStrip[20000];
  28. double xtrk[20000];
  29. double ytrk[20000];
  30. // Various arrays of parameters common to this macro
  31. double tdcTime = 0;
  32. double adcs[512] = {0};
  33. double finalSignal[512] = {0};
  34.  
  35. int NCluster = 0;
  36. float ClusterCharge[100] = {0};
  37. int ClusterSize[100] ={0};
  38. float ClusterPosY[100] ={0};
  39. float seedChannel[100] ={0};
  40.  
  41. int iOffSet = 0; // Offset into tree, to get to specific events
  42. int printEventDisplays = 0; // Number of event displays to print to file
  43. float polarity = -1.0;
  44.  
  45.  
  46. /////////////////////////////////////////////////////////////////////////////////
  47.  
  48. void plotADCs(){
  49. gStyle->SetNdivisions(505,"Y");
  50. gStyle->SetNdivisions(505,"X");
  51. gStyle->SetOptStat(0);
  52. gStyle->SetOptFit(1);
  53. gStyle->SetPadTickX(1);
  54. gStyle->SetPadTickY(1);
  55.  
  56. //TString sigFile = getFileToOpen(runNum);
  57. TString sigFile = "/data2/sblusk/TB/July2015/AnaFiles/BoardA6/Run_Bias_Scan-B6-A-212-8358_Tuple.root";
  58.  
  59. // Open signal file
  60. TFile *fs = new TFile(sigFile);
  61. TTree* cmstree = (TTree*)fs->Get("TbUT/CMS");
  62. TTree* fChain = (TTree*)fs->Get("TbUT/Clusters");
  63. int nevents = cmstree->GetEntriesFast();
  64.  
  65. // Strip level data
  66. Double_t cmsData[512];
  67. // List of branches
  68. TBranch *b_cmsData; //!
  69. cmstree->SetBranchAddress("cmsData", cmsData, &b_cmsData);
  70.  
  71. // Clusters
  72. // Declaration of leaf types
  73. Int_t clusterNumberPerEvent;
  74. UInt_t clustersTDC;
  75. ULong64_t timestamps;
  76. Double_t clustersPosition[10];
  77. Int_t clustersSeedPosition[10];
  78. Double_t clustersCharge[10];
  79. Int_t clustersSize[10];
  80. Double_t clustersSeedCharge[10];
  81. Double_t clustersCharge2StripLeft[10];
  82. Double_t clustersCharge1StripLeft[10];
  83. Double_t clustersCharge1StripRight[10];
  84. Double_t clustersCharge2StripRight[10];
  85. TBranch *b_clusterNumberPerEvent; //!
  86. TBranch *b_clustersTDC; //!
  87. TBranch *b_timestamps; //!
  88. TBranch *b_clustersPosition; //!
  89. TBranch *b_clustersSeedPosition; //!
  90. TBranch *b_clustersCharge; //!
  91. TBranch *b_clustersSize; //!
  92. TBranch *b_clustersSeedCharge; //!
  93. TBranch *b_clustersCharge2StripLeft; //!
  94. TBranch *b_clustersCharge1StripLeft; //!
  95. TBranch *b_clustersCharge1StripRight; //!
  96. TBranch *b_clustersCharge2StripRight; //!
  97. fChain->SetBranchAddress("clusterNumberPerEvent", &clusterNumberPerEvent, &b_clusterNumberPerEvent);
  98. fChain->SetBranchAddress("clustersTDC", &clustersTDC, &b_clustersTDC);
  99. fChain->SetBranchAddress("timestamps", &timestamps, &b_timestamps);
  100. fChain->SetBranchAddress("clustersPosition", clustersPosition, &b_clustersPosition);
  101. fChain->SetBranchAddress("clustersSeedPosition", clustersSeedPosition, &b_clustersSeedPosition);
  102. fChain->SetBranchAddress("clustersCharge", clustersCharge, &b_clustersCharge);
  103. fChain->SetBranchAddress("clustersSize", clustersSize, &b_clustersSize);
  104. fChain->SetBranchAddress("clustersSeedCharge", clustersSeedCharge, &b_clustersSeedCharge);
  105. fChain->SetBranchAddress("clustersCharge2StripLeft", clustersCharge2StripLeft, &b_clustersCharge2StripLeft);
  106. fChain->SetBranchAddress("clustersCharge1StripLeft", clustersCharge1StripLeft, &b_clustersCharge1StripLeft);
  107. fChain->SetBranchAddress("clustersCharge1StripRight", clustersCharge1StripRight, &b_clustersCharge1StripRight);
  108. fChain->SetBranchAddress("clustersCharge2StripRight", clustersCharge2StripRight, &b_clustersCharge2StripRight);
  109. TH1F *h = new TH1F("h","ADC, CM + Step Subtracted",512,0.0,512);
  110. TH1F *hm1 = new TH1F("hm1","ADC of nomStrip",110,-100.0,1000);
  111. TH1F *hm = new TH1F("hm","ADC of nomStrip#pm1",110,-100.0,1000);
  112. h->SetTitle("");
  113. h->SetNdivisions(510,"X");
  114.  
  115. TCanvas *theCanvas;
  116. int nPrint = 0;
  117. char input;
  118. int nEvents;
  119. int jentry;
  120. double strip, xp, yp;
  121.  
  122. std::ifstream infile;
  123. if(plotSelectedEvents){
  124. // Read in event list
  125. TString filename = "MissingDUTHits.dat";
  126. cout << "------------------------------------------------------------" << endl;
  127. cout << " You've selected to view events from file: " << filename << endl;
  128. cout << "------------------------------------------------------------" << endl;
  129. infile.open(filename);
  130. if(!infile) { // file couldn't be opened
  131. cerr << "Error: file could not be opened" << endl;
  132. exit(1);
  133. }
  134. int il = 0;
  135. while ( !infile.eof() ) { // keep reading until end-of-file
  136. infile >> jentry >> strip >> xp >> yp;
  137. if(il>9999) break;
  138. evList[il] = jentry;
  139. nomStrip[il] = strip;
  140. xtrk[il] = xp;
  141. ytrk[il] = yp;
  142. il++;
  143. }
  144. infile.close();
  145. nEvents = nevents - 1000;
  146. }else{
  147. cout << "How many events do you want to look at ?" << endl;
  148. cin >> nEvents;
  149. cout << "--------------------------------------"<<endl;
  150. }
  151.  
  152. //------------------------
  153. // Begin Loop over events
  154. //------------------------
  155. il = 0;
  156. for(int i=iOffSet; i<iOffSet+nEvents; i++){
  157. //if(i%1==0) cout << "Processing Event " << i << endl;
  158. if(plotSelectedEvents){
  159. if(i != evList[il]) continue;
  160. //if(showEventDisplay){
  161. // cout << "-------------------------------------------------------------------------------"<<endl;
  162. // cout << "Found selected event: " << i << ", NomStrip of Missing Hit = " << nomStrip[il] << endl;
  163. // cout << "-------------------------------------------------------------------------------"<<endl;
  164. //}
  165. }
  166.  
  167. cmstree->GetEntry(i);
  168. fChain->GetEntry(i);
  169. int k=0;
  170.  
  171. //--------------
  172. // Get CMS data
  173. //--------------
  174. for(int j=0; j<512;j++){
  175. adcs[j] = cmsData[j];
  176. }
  177.  
  178. //------------------
  179. // Get Cluster data
  180. //------------------
  181. NCluster = clusterNumberPerEvent;
  182. if(NCluster>=10) NCluster=10;
  183. if(showEventDisplay) cout << "# Clusters = " << NCluster << endl;
  184. for(int j=0; j<NCluster; j++){
  185. ClusterSize[j] = clustersSize[j];
  186. ClusterPosY[j] = 1.0*clustersPosition[j];
  187. ClusterCharge[j] = clustersCharge[j];
  188. tdcTime = clustersTDC;
  189. }
  190.  
  191.  
  192. //-------------------------
  193. // Get signal distributions
  194. //-------------------------
  195. getSignal();
  196.  
  197. int nstrip = 0;
  198. double vstrip1 = 0, vstrip3 = 0;
  199. if(plotSelectedEvents){
  200. nstrip = nomStrip[il]+0.5;
  201. vstrip1 = finalSignal[nstrip];
  202. if(nstrip>1) vstrip3 = vstrip1 + finalSignal[nstrip-1];
  203. if(nstrip<511) vstrip3 = vstrip3 + finalSignal[nstrip+1];
  204. vstrip1 = -1.0*vstrip1;
  205. vstrip3 = -1.0*vstrip3;
  206. }
  207. //------------------------------------------------------------
  208. // Single Event displays of pedestal & noise subtracted data -
  209. //------------------------------------------------------------
  210. if(showEventDisplay){
  211. theCanvas = plotSingleEvent(i, h, nomStrip[il]);
  212. cout << "NomStrip = " << nstrip << " .... Charge in nomStrip and nomStrip +- 1: " << vstrip1 << " " << vstrip3 << endl;
  213. printf("\n <<< ENTER / RETURN to continue >>> \n");
  214. theCanvas->WaitPrimitive();
  215. if(nPrint < printEventDisplays){
  216. theCanvas->Print(Form("./Plots/EventDisplayCorr_%d.png",i));
  217. theCanvas->Print(Form("./Plots/EventDisplayCorr_%d.pdf",i));
  218. nPrint++;
  219. }
  220. }
  221. //----------------
  222. // Fill histograms
  223. //----------------
  224. if(fillStatDistribution && plotSelectedEvents){
  225. hm1->Fill(-1.0*vstrip1);
  226. hm->Fill(-1.0*vstrip3);
  227. }
  228. if(plotSelectedEvents) il++;
  229. }
  230.  
  231. // Plot statistical distributions
  232. TCanvas *c = new TCanvas("c","Plots",500,800);
  233. c->Divide(1,2);
  234. c->cd(1);
  235. hm1->Draw();
  236. c->cd(2);
  237. hm->Draw();
  238. }
  239.  
  240. void getSignal(){
  241. // -----------------
  242. // - Final result -
  243. // - Here, one could make any "corrections", or flip polarity, etc
  244. // -----------------
  245. for(int k=0;k<256;k++){finalSignal[k] = adcs[k];}
  246. return;
  247. }
  248.  
  249.  
  250. TCanvas* plotSingleEvent(int i, TH1F *h, double xstrip)
  251. {
  252. //-----------------
  253. // Final result
  254. //-----------------
  255. double low = 1000;
  256. double hi = -1000;
  257.  
  258. h->Reset();
  259. for(int k=0;k<512;k++) {
  260. h->Fill(k,finalSignal[k]);
  261. if(finalSignal[k]<low) low = finalSignal[k];
  262. if(finalSignal[k]>hi) hi = finalSignal[k];
  263. }
  264. h->SetMinimum(low-100);
  265. h->SetMaximum(hi+300);
  266. h->GetXaxis()->SetRangeUser(lowCh,hiCh);
  267.  
  268. TCanvas *c = new TCanvas("c","",500,0,800,800);
  269. h->Draw();
  270. TLine *l = new TLine(lowCh, 0, hiCh, 0);
  271. TLine *l2 = new TLine(xstrip, low-100, xstrip, hi+150);
  272. l2->SetLineColor(5); l2->SetLineWidth(4); l2->Draw();
  273. l->SetLineColor(3); l->SetLineStyle(2); l->SetLineWidth(2); l->Draw();
  274. h->Draw("same");
  275.  
  276. TLatex *myLatex = new TLatex();
  277. myLatex->SetTextFont(42); myLatex->SetTextColor(1);
  278. myLatex->SetTextAlign(12); myLatex->SetNDC(kTRUE); myLatex->SetTextSize(0.065);
  279. myLatex->DrawLatex(0.14,0.85,Form("Event %d",i));
  280. TString tdcs = Form("TDC time = %4.1f ns",tdcTime);
  281. myLatex->SetTextSize(0.045);myLatex->DrawLatex(0.18,0.79,tdcs);
  282.  
  283.  
  284. float yval = 0.79;
  285. myLatex->SetTextSize(0.035);
  286. yval = yval - 0.05;
  287. for(int k=0;k<NCluster;k++){
  288. TString clString = Form("Size=%d, Q=%3.0f, Pos=%4.1f",ClusterSize[k],ClusterCharge[k],ClusterPosY[k]);
  289. myLatex->DrawLatex(0.14,yval,clString);
  290. yval = yval - 0.05;
  291. }
  292.  
  293. return c;
  294. }
  295.  
  296.