HARPO  5.1.1
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
HarpoAnalyseRunNoZS.cxx
Go to the documentation of this file.
1 //
2 // File HarpoAnalyseRunNoZS.cxx
3 //
12 #include "HarpoAnalyseRunNoZS.h"
13 #include "HarpoConfig.h"
14 #include "HarpoDetSet.h"
15 #include "HarpoDebug.h"
16 #include "HarpoDccEvent.h"
17 #include "Pmm2Event.h"
18 #include "HarpoEvent.h"
19 #include "HarpoRecoEvent.h"
20 #include "MakeNiceHisto.h"
21 
22 #include "TFile.h"
23 #include "TStyle.h"
24 #include "TCanvas.h"
25 #include "TLatex.h"
26 #include "TGraph.h"
27 #include "TF1.h"
28 #include "TMath.h"
29 #include "TSystem.h"
30 
31 #include <cstdlib>
32 #include <cstring>
33 #include <cassert>
34 #include <fstream>
35 #include <iostream>
36 
37 ClassImp(HarpoAnalyseRunNoZS);
38 
40  {
41  unsigned int nd; // number of detectors
43 
44  assert(hdr != NULL);
45  hdr->print();
46 
47  for (nd = 0; nd < gkNDetectors; nd++) {
48  // if (fEvt->isdataExist(nd)) {
49  HarpoDccMap *plane = fEvt->GetDccMap(nd);
50  if (plane != NULL )plane->print();
51  }
52  }
53 
55 {
56  // Bool_t drawEvent = kFALSE;
57  nEvents++;
58 
59  Info("process","Event %ld",nEvents);
60 
61  Int_t startSatTmp = -1;
62 
63 
64  Double_t mean[2][NCHAN];
65  Double_t sigma[2][NCHAN];
66  Double_t qtot[2][NCHAN];
67  // Double_t widt[2][NCHAN];
68 
69 
70  for(Int_t ndet=0;ndet<2;ndet++) { // Detectors 0 and 1 = X and Y planes
71  for(Int_t j=0;j<NCHAN;j++){ // Channels
72  mean[ndet][j] = 0;
73  sigma[ndet][j] = 0;
74  qtot[ndet][j] = 0;
75  // widt[ndet][j] = 0;
76  }
77  }
78 
79 #if 0
80  // Process RAW data
81  for(Int_t ndet=0;ndet<2;ndet++) { // Detectors 0 and 1 = X and Y planes
82  if (!gHDetSet->isExist(ndet)) continue;
83  HarpoDccMap* m = fEvt->GetDccMap(ndet);
84  if ( m == NULL ) continue;
85 
86  for(Int_t j=0;j<NCHAN;j++){ // Channels
87  Double_t meanQ = 0, sigmaQ = 0;
88  TArrayD* arr = new TArrayD(NADC);
89  Int_t nPositive = 0;
90  for(Int_t i=0;i<NADC;i++){ //Time bins
91  Double_t q = m->GetData(j,i);
92  //if(q>40) continue;
93  // meanQ += q;
94  // sigmaQ += q*q;
95  if(q<40)
96  arr->AddAt(q,i);
97  else
98  arr->AddAt(-10000,i);
99  if(q>0) nPositive++;
100  }
101  hNpositive[ndet]->Fill(nPositive);
102  for(Int_t i=0;i<NADC;i++){ //Time bins
103  Double_t q = m->GetData(j,i);
104  if(nPositive>320){
105  hTest[ndet]->Fill(i,q);
106  continue;
107  }
108  if(nPositive<180){
109  hTest2[ndet]->Fill(i,q);
110  continue;
111  }
112  hTest3[ndet]->Fill(i,q);
113  }
114 
115 
116  TruncSigma(arr,sigmaQ,meanQ,0.2,0.8);
117  sigma[ndet][j] = sigmaQ;
118  mean[ndet][j] = meanQ;
119 
120  hMeanCh[ndet]->Fill(meanQ);
121  hSigmaCh[ndet]->Fill(sigmaQ);
122  hSigmaMeanCh[ndet]->Fill(meanQ,sigmaQ);
123  hSigmaChVsCh[ndet]->Fill(j,sigmaQ);
124  hMeanChVsCh[ndet]->Fill(j,meanQ);
125 
126  }
127  }
128 
129  HarpoDccMap* m0 = fEvt->GetDccMap(0);
130  if ( m0 == NULL ) return;
131  HarpoDccMap* m1 = fEvt->GetDccMap(1);
132  if ( m1 == NULL ) return;
133  for(Int_t j0=0;j0<NCHAN;j0++){ // Channels
134  for(Int_t j1=0;j1<NCHAN;j1++){ // Channels
135 
136  // meanQ /= NADC;
137  // sigmaQ = TMath::Sqrt((sigmaQ/NADC - meanQ*meanQ));
138  if(sigma[0][j0] == 0) continue;
139  if(sigma[1][j1] == 0) continue;
140 
141  Int_t startSC = -1, startSat = -1, startSat2 = -1;
142  Int_t width = 0, width1 = 0, qTot1 = 0;
143 
144  for(Int_t i=0;i<NADC;i++){ //Time bins
145  Double_t q0 = m0->GetData(j0,i);
146  Double_t q1 = m1->GetData(j1,i);
147  if(q0>mean[0][j0]+5*sigma[0][j0] && q1>mean[1][j1]+5*sigma[1][j1]){
148  width++;
149  width1++;
150  qTot1 += q0 + q1;
151  startSC = i;
152  }else{
153  width1 = 0;
154  qTot1 = 0;
155  }
156  if(q0>3700 && q1>3700)
157  startSat = i;
158  if(q0>100 && q1>100)
159  startSat2 = i;
160  }
161  if(startSat2 != -1)
162  startSat2 = fStartSatOld;
163  if(startSat!=-1)
164  startSatTmp = startSat;
165  Double_t qTot = 0;
166  for(Int_t i=0;i<NADC;i++){ //Time bins
167  Double_t q0 = m0->GetData(j0,i);
168  Double_t q1 = m1->GetData(j1,i);
169  Double_t w = 1;
170  Double_t q = (q0+q1)/2;
171  if(q<50) w = 10;
172  if(startSat!=-1){
173  hSaturation[0]->Fill(i-startSat,q0,w);
174  hSaturation[1]->Fill(i-startSat,q1,w);
175  }
176  if(startSat2!=-1){
177  hSaturation2[0]->Fill(i-startSat2,q0,w);
178  hSaturation2[1]->Fill(i-startSat2,q1,w);
179  }
180  if(mean[0][j0]<5 && sigma[0][j0]>5 && mean[1][j1]<5 && sigma[1][j1]>5){
181  hSpaceChargeRef[0]->Fill(i,q0);
182  hSpaceCharge2Ref[0]->Fill(i-startSC,q0,w);
183  hSpaceChargeRef[1]->Fill(i,q1);
184  hSpaceCharge2Ref[1]->Fill(i-startSC,q1,w);
185  }else{
186  hSpaceCharge[0]->Fill(i,q0);
187  hSpaceCharge2[0]->Fill(i-startSC,q0,w);
188  hSpaceCharge[1]->Fill(i,q1);
189  hSpaceCharge2[1]->Fill(i-startSC,q1,w);
190  }
191  if(i<startSC)
192  qTot += q;
193  }
194  if(mean[0][j0]<5 && sigma[0][j0]>5 && mean[1][j1]<5 && sigma[1][j1]>5){
195  hQtotRef[0]->Fill(qTot);
196  hQtotWidthRef[0]->Fill(width,qTot);
197  hWidthRef[0]->Fill(width);
198  // hQtotRef[ndet]->Fill(qTot);
199  // hQtotWidthRef[ndet]->Fill(width,qTot);
200  // hWidthRef[ndet]->Fill(width);
201  }else{
202  hQtot[0]->Fill(qTot);
203  hQtotWidth[0]->Fill(width,qTot);
204  hWidth[0]->Fill(width);
205  }
206  }
207  }
208 
209 
210 
211 
212 #else
213 
214 
215  // Process RAW data
216  for(Int_t ndet=0;ndet<2;ndet++) { // Detectors 0 and 1 = X and Y planes
217  if (!gHDetSet->isExist(ndet)) continue;
218  HarpoDccMap* m = fEvt->GetDccMap(ndet);
219  if ( m == NULL ) continue;
220 
221  // for(Int_t i=0;i<NADC;i++){ //Time bins
222  // for(Int_t j=0;j<NCHAN;j++){ // Channels
223  // Double_t q = m->GetData(j,i);
224  // if(q == -1000) continue;
225  // if(q <= 0){
226  // hNegative[ndet]->Fill(i,j);
227  // hQnegative[ndet]->Fill(i,j,q);
228  // }
229  // hQ[ndet]->Fill(i,j,q);
230  // hQraw[ndet]->Fill(q);
231  // }
232  // }
233  for(Int_t j=0;j<NCHAN;j++){ // Channels
234  Double_t meanQ = 0, sigmaQ = 0;
235  TArrayD* arr = new TArrayD(NADC);
236  Int_t nPositive = 0;
237  for(Int_t i=0;i<NADC;i++){ //Time bins
238  Double_t q = m->GetData(j,i);
239  //if(q>40) continue;
240  // meanQ += q;
241  // sigmaQ += q*q;
242  if(q<40)
243  arr->AddAt(q,i);
244  else
245  arr->AddAt(-10000,i);
246  if(q>0) nPositive++;
247  }
248  hNpositive[ndet]->Fill(nPositive);
249  for(Int_t i=0;i<NADC;i++){ //Time bins
250  Double_t q = m->GetData(j,i);
251  if(nPositive>320){
252  hTest[ndet]->Fill(i,q);
253  continue;
254  }
255  if(nPositive<180){
256  hTest2[ndet]->Fill(i,q);
257  continue;
258  }
259  hTest3[ndet]->Fill(i,q);
260  }
261 
262 
263  TruncSigma(arr,sigmaQ,meanQ,0.2,0.8);
264  sigma[ndet][j] = sigmaQ;
265  mean[ndet][j] = meanQ;
266 
267  // meanQ /= NADC;
268  // sigmaQ = TMath::Sqrt((sigmaQ/NADC - meanQ*meanQ));
269  hMeanCh[ndet]->Fill(meanQ);
270  hSigmaCh[ndet]->Fill(sigmaQ);
271  hSigmaMeanCh[ndet]->Fill(meanQ,sigmaQ);
272  hSigmaChVsCh[ndet]->Fill(j,sigmaQ);
273  hMeanChVsCh[ndet]->Fill(j,meanQ);
274 
275 
276  Int_t startSC = -1, startSat = -1, startSat2 = -1;
277  Int_t width = 0, width1 = 0, qTot1 = 0;
278  for(Int_t i=0;i<NADC;i++){ //Time bins
279  Double_t q = m->GetData(j,i);
280  Double_t qraw = m->GetRawData(j,i);
281  if(q>meanQ+5*sigmaQ){
282  width++;
283  width1++;
284  qTot1 += q;
285  startSC = i;
286  }else{
287  // if(meanQ<5 && sigmaQ>5){
288  // hQtotWidthRef[ndet]->Fill(width,qTot1);
289  // hWidthRef[ndet]->Fill(width);
290  // }else{
291  // hQtotWidth[ndet]->Fill(width,qTot1);
292  // hWidth[ndet]->Fill(width);
293  // }
294  width1 = 0;
295  qTot1 = 0;
296  }
297  if(qraw>4090)
298  startSat = i;
299  if(q>100)
300  startSat2 = i;
301  }
302  if(startSat2 != -1)
303  startSat2 = fStartSatOld;
304  if(startSat!=-1)
305  startSatTmp = startSat;
306  Double_t qTot = 0;
307  for(Int_t i=0;i<NADC;i++){ //Time bins
308  Double_t q = m->GetData(j,i);
309  Double_t w = 1;
310  if(q<50) w = 10;
311  if(startSat!=-1)
312  hSaturation[ndet]->Fill(i-startSat,q,w);
313  if(startSat2!=-1)
314  hSaturation2[ndet]->Fill(i-startSat2,q,w);
315  if(meanQ<5 && sigmaQ>5){
316  hSpaceChargeRef[ndet]->Fill(i,q);
317  hSpaceCharge2Ref[ndet]->Fill(i-startSC,q,w);
318  }else{
319  hSpaceCharge[ndet]->Fill(i,q);
320  hSpaceCharge2[ndet]->Fill(i-startSC,q,w);
321  }
322  if(i<startSC)
323  qTot += q;
324  }
325  qtot[ndet][j] = qTot;
326  // widt[ndet][j] = width;
327  if(meanQ<5 && sigmaQ>5){
328  // hQtotRef[ndet]->Fill(qTot);
329  hQtotWidthRef[ndet]->Fill(width,qTot);
330  hWidthRef[ndet]->Fill(width);
331  }else{
332  // hQtot[ndet]->Fill(qTot);
333  hQtotWidth[ndet]->Fill(width,qTot);
334  hWidth[ndet]->Fill(width);
335  }
336  }
337  }
338  for(Int_t j0=0;j0<NCHAN;j0++){ // Channels
339  for(Int_t j1=0;j1<NCHAN;j1++){ // Channels
340 
341  // meanQ /= NADC;
342  // sigmaQ = TMath::Sqrt((sigmaQ/NADC - meanQ*meanQ));
343  if(sigma[0][j0] == 0) continue;
344  if(sigma[1][j1] == 0) continue;
345 
346  if(mean[0][j0]<5 && sigma[0][j0]>5 && mean[1][j1]<5 && sigma[1][j1]>5)
347  hQtotRef[0]->Fill((qtot[0][j0] + qtot[0][j0])/2);
348  else
349  hQtot[0]->Fill((qtot[0][j0] + qtot[0][j0])/2);
350  }
351  }
352 #endif
353 
354  // if(startSatTmp!=-1)
355  fStartSatOld = startSatTmp;
356 
357 }
358 
360  {
361 
362  fStartSatOld = -1;
363 
364  const Int_t nbinsQ1 = 100;
365  const Int_t nbinsQ2 = 395;
366  Double_t xminQ1 = -50;
367  Double_t xminQ2 = xminQ1 + nbinsQ1;
368  // Double_t xmaxQ = 4000;
369  Double_t xbinsQ[nbinsQ1+nbinsQ2+1];
370  // xbinsQ[0] = xminQ;
371  // xbinsQ[nbinsQ1+nbinsQ2] = xmaxQ;
372  for (Int_t i=0;i<=nbinsQ1;i++)
373  xbinsQ[i] = xminQ1 + i;
374  for (Int_t i=0;i<=nbinsQ2;i++)
375  xbinsQ[nbinsQ1+i] = xminQ2 + 10*i;
376 
377 
378  hSpaceCharge[0] = new TH2F("hSpaceChargeX","map",512,0,512,500,-100,400);
379  hSpaceCharge[1] = new TH2F("hSpaceChargeY","map",512,0,512,500,-100,400);
380  hSpaceCharge2[0] = new TH2F("hSpaceCharge2X","map",1024,-512,512,nbinsQ1+nbinsQ2,xbinsQ);
381  hSpaceCharge2[1] = new TH2F("hSpaceCharge2Y","map",1024,-512,512,nbinsQ1+nbinsQ2,xbinsQ);
382  hSpaceChargeRef[0] = new TH2F("hSpaceChargeRefX","map",512,0,512,500,-100,400);
383  hSpaceChargeRef[1] = new TH2F("hSpaceChargeRefY","map",512,0,512,500,-100,400);
384  hSpaceCharge2Ref[0] = new TH2F("hSpaceCharge2RefX","map",1024,-512,512,nbinsQ1+nbinsQ2,xbinsQ);
385  hSpaceCharge2Ref[1] = new TH2F("hSpaceCharge2RefY","map",1024,-512,512,nbinsQ1+nbinsQ2,xbinsQ);
386  hSaturation[0] = new TH2F("hSaturationX","map",1024,-512,512,nbinsQ1+nbinsQ2,xbinsQ);
387  hSaturation[1] = new TH2F("hSaturationY","map",1024,-512,512,nbinsQ1+nbinsQ2,xbinsQ);
388  hSaturation2[0] = new TH2F("hSaturation2X","map",1024,-512,512,nbinsQ1+nbinsQ2,xbinsQ);
389  hSaturation2[1] = new TH2F("hSaturation2Y","map",1024,-512,512,nbinsQ1+nbinsQ2,xbinsQ);
390 
391 
392  // Initialise histograms here
393  hWidth[0] = new TH1F("hWidthX","map",512,0,512);
394  hWidth[1] = new TH1F("hWidthY","map",512,0,512);
395  hWidthRef[0] = new TH1F("hWidthRefX","map",512,0,512);
396  hWidthRef[1] = new TH1F("hWidthRefY","map",512,0,512);
397  hNpositive[0] = new TH1F("hNpositiveX","map",512,0,512);
398  hNpositive[1] = new TH1F("hNpositiveY","map",512,0,512);
399  hTest[0] = new TH2F("hTestX","map",512,0,512,500,-100,400);
400  hTest[1] = new TH2F("hTestY","map",512,0,512,500,-100,400);
401  hTest2[0] = new TH2F("hTest2X","map",512,0,512,500,-100,400);
402  hTest2[1] = new TH2F("hTest2Y","map",512,0,512,500,-100,400);
403  hTest3[0] = new TH2F("hTest3X","map",512,0,512,500,-100,400);
404  hTest3[1] = new TH2F("hTest3Y","map",512,0,512,500,-100,400);
405  // hNegative[0] = new TH2F("hNegativeX","map",512,0,512,288,0,288);
406  // hNegative[1] = new TH2F("hNegativeY","map",512,0,512,288,0,288);
407  // hQnegative[0] = new TH2F("hQnegativeX","map",512,0,512,288,0,288);
408  // hQnegative[1] = new TH2F("hQnegativeY","map",512,0,512,288,0,288);
409  // hQ[0] = new TH3F("hQX","",512,0,512,288,0,288,1000,-300,4700);
410  // hQ[1] = new TH3F("hQY","",512,0,512,288,0,288,1000,-300,4700);
411  // hQraw[0] = new TH1F("hQrawX",";Q_{raw}",4096,0,4096);
412  // hQraw[1] = new TH1F("hQrawY",";Q_{raw}",4096,0,4096);
413 
414  hSigmaCh[0] = new TH1F("hSigmaChX",";#sigma_{Q}",1000,0,20);
415  hSigmaCh[1] = new TH1F("hSigmaChY",";#sigma_{Q}",1000,0,20);
416  hMeanCh[0] = new TH1F("hMeanChX",";Q_{raw}",1000,-50,50);
417  hMeanCh[1] = new TH1F("hMeanChY",";Q_{raw}",1000,-50,50);
418  hSigmaMeanCh[0] = new TH2F("hSigmaMeanChX",";#sigma_{Q};Q_{raw}",1000,-50,50,500,0,20);
419  hSigmaMeanCh[1] = new TH2F("hSigmaMeanChY",";#sigma_{Q};Q_{raw}",1000,-50,50,500,0,20);
420  hSigmaChVsCh[0] = new TH2F("hSigmaChVsChX",";#sigma_{Q};Q_{raw}",288,0,288,500,0,20);
421  hSigmaChVsCh[1] = new TH2F("hSigmaChVsChY",";#sigma_{Q};Q_{raw}",288,0,288,500,0,20);
422  hMeanChVsCh[0] = new TH2F("hMeanChVsChX",";#sigma_{Q};Q_{raw}",288,0,288,1000,-50,50);
423  hMeanChVsCh[1] = new TH2F("hMeanChVsChY",";#sigma_{Q};Q_{raw}",288,0,288,1000,-50,50);
424 
425 
426  const Int_t nbinsQcl = 500;
427  Double_t xminQcl = 10;
428  Double_t xmaxQcl = 1e6;
429  Double_t logxminQcl = TMath::Log(xminQcl);
430  Double_t logxmaxQcl = TMath::Log(xmaxQcl);
431  Double_t binwidthQcl = (logxmaxQcl-logxminQcl)/nbinsQcl;
432  Double_t xbinsQcl[nbinsQcl+1];
433  xbinsQcl[0] = xminQcl;
434  xbinsQcl[nbinsQcl] = xmaxQcl;
435  for (Int_t i=1;i<nbinsQcl;i++)
436  xbinsQcl[i] = TMath::Exp(logxminQcl+i*binwidthQcl);
437 
438  hQtot[0] = new TH1F("hQtotX",";Q_{tot}",nbinsQcl,xbinsQcl);
439  hQtot[1] = new TH1F("hQtotY",";Q_{tot}",nbinsQcl,xbinsQcl);
440  hQtotRef[0] = new TH1F("hQtotRefX",";Q_{tot}",nbinsQcl,xbinsQcl);
441  hQtotRef[1] = new TH1F("hQtotRefY",";Q_{tot}",nbinsQcl,xbinsQcl);
442  hQtotWidth[0] = new TH2F("hQtotWidthX",";Q_{tot}",512,0,512,nbinsQcl,xbinsQcl);
443  hQtotWidth[1] = new TH2F("hQtotWidthY",";Q_{tot}",512,0,512,nbinsQcl,xbinsQcl);
444  hQtotWidthRef[0] = new TH2F("hQtotWidthRefX",";Q_{tot}",512,0,512,nbinsQcl,xbinsQcl);
445  hQtotWidthRef[1] = new TH2F("hQtotWidthRefY",";Q_{tot}",512,0,512,nbinsQcl,xbinsQcl);
446 
447  }
448 
449 void HarpoAnalyseRunNoZS::Save(char * /* mode */)
450  {
451 
452 
453  OpenHistFile("runnozs");
454 
455  // Save histograms here
456  hSpaceCharge[0]->Write();
457  hSpaceCharge[1]->Write();
458  hSpaceCharge2[0]->Write();
459  hSpaceCharge2[1]->Write();
460  hSpaceChargeRef[0]->Write();
461  hSpaceChargeRef[1]->Write();
462  hSpaceCharge2Ref[0]->Write();
463  hSpaceCharge2Ref[1]->Write();
464  hQtot[0]->Write();
465  hQtot[1]->Write();
466  hQtotRef[0]->Write();
467  hQtotRef[1]->Write();
468  hQtotWidth[0]->Write();
469  hQtotWidth[1]->Write();
470  hQtotWidthRef[0]->Write();
471  hQtotWidthRef[1]->Write();
472  hWidth[0]->Write();
473  hWidth[1]->Write();
474  hWidthRef[0]->Write();
475  hWidthRef[1]->Write();
476  hQtotWidth[0]->Write();
477  hQtotWidth[1]->Write();
478  hQtotWidthRef[0]->Write();
479  hQtotWidthRef[1]->Write();
480 
481  hSaturation[0]->Write();
482  hSaturation[1]->Write();
483  hSaturation2[0]->Write();
484  hSaturation2[1]->Write();
485  hTest[0]->Write();
486  hTest[1]->Write();
487  hTest2[0]->Write();
488  hTest2[1]->Write();
489  hTest3[0]->Write();
490  hTest3[1]->Write();
491  hNpositive[0]->Write();
492  hNpositive[1]->Write();
493  // hNegative[0]->Write();
494  // hNegative[1]->Write();
495  // hQnegative[0]->Write();
496  // hQnegative[1]->Write();
497  // hQraw[0]->Write();
498  // hQraw[1]->Write();
499  // hQ[0]->Write();
500  // hQ[1]->Write();
501  hMeanCh[0]->Write();
502  hMeanCh[1]->Write();
503  hSigmaCh[0]->Write();
504  hSigmaCh[1]->Write();
505  hSigmaMeanCh[0]->Write();
506  hSigmaMeanCh[1]->Write();
507  hSigmaChVsCh[0]->Write();
508  hSigmaChVsCh[1]->Write();
509  hMeanChVsCh[0]->Write();
510  hMeanChVsCh[1]->Write();
511 
512  }
513 
514 
515 
516 Double_t HarpoAnalyseRunNoZS::TruncSigma(TArrayD* vect, Double_t &truncS, Double_t &truncM, Double_t tl, Double_t th)
517 {
518  //
519  // Calculates the truncated mean mean of the non zero value in vect
520  // The truncation is done on the 100*tl% lowest and 100*(1-th) highest value
521  //
522  // debug(2,"truncated mean");
523 
524 
525  Int_t size = vect->GetSize();
526  // Info("TruncSigma","size = %d",size);
527  Double_t truncMean = 0;
528  Double_t truncSigma = 0;
529  Int_t* index = new Int_t[size];
530  TMath::Sort(size,vect->GetArray(),index,kFALSE);
531  Int_t t = 0, tLow, tHigh;
532 
533 
534  while(vect->At(index[t])==-10000&&t<size+1) t++;
535  // while(t<size+1) t++;
536  tLow = TMath::FloorNint(t + (size - t)*tl);
537  tHigh = TMath::FloorNint(t + (size - t)*th);
538  //tLow = TMath::FloorNint(size*tl);
539  //tHigh = TMath::FloorNint(size*th);
540 
541  for(Int_t tind = tLow; tind<tHigh; tind++){
542  truncMean += vect->At(index[tind]);
543  truncSigma += vect->At(index[tind])*vect->At(index[tind]);
544  }
545 
546  delete[] index;
547 
548  if(tHigh-tLow == 0) return 0;
549  //if(tHigh-tLow < 15) return 0;
550 
551  truncMean /= (tHigh-tLow);
552  truncSigma /= (tHigh-tLow);
553  truncSigma -= truncMean*truncMean;
554 
555  truncM = truncMean;
556  truncS = TMath::Sqrt(truncSigma);
557 
558  return truncSigma;
559 }
#define NCHAN
Definition: HarpoDccMap.h:16
Int_t fStartSatOld
Redefine empty default.
Double_t GetRawData(Int_t i, Int_t j)
Definition: HarpoDccMap.cxx:75
Double_t GetData(Int_t i, Int_t j)
Set/Get Data Cell.
Definition: HarpoDccMap.cxx:84
Bool_t isExist(ULong_t det)
Detecror date exist //! Number of Real Detectors.
Definition: HarpoDetSet.h:33
Double_t TruncSigma(TArrayD *vect, Double_t &truncS, Double_t &truncM, Double_t tl, Double_t th)
virtual void print()
FullEvent Header not scecific to the detectors The class is ....
void print()
Ovreloaded medod whic do all job.
virtual void print()
TFile * OpenHistFile(const char *ananame)
unpacked dcc data The class contains the data map for DCC or Feminos The data is stored as a 2D TMatr...
Definition: HarpoDccMap.h:29
Dummy analysis to run as test and example. Give basic histograms of the data.
#define NADC
Definition: HarpoDccMap.h:18
HarpoEventHeader * GetHeader()
Definition: HarpoEvent.cxx:80
HarpoEvent * fEvt
Definition: HarpoAnalyse.h:70
ULong_t nEvents
Definition: HarpoAnalyse.h:75
const ULong_t gkNDetectors
Definition: HarpoDet.h:14
HarpoDccMap * GetDccMap(Long_t plane=XDCC)
Definition: HarpoEvent.cxx:108
R__EXTERN HarpoDetSet * gHDetSet
Definition: HarpoDetSet.h:71
void Save(char *mode=NULL)