source: trunk/MagicSoft/Mars/mraw/MRawEvtData.cc@ 4642

Last change on this file since 4642 was 4576, checked in by tbretz, 20 years ago
*** empty log message ***
File size: 16.3 KB
Line 
1/* ======================================================================== *\
2!
3! *
4! * This file is part of MARS, the MAGIC Analysis and Reconstruction
5! * Software. It is distributed to you in the hope that it can be a useful
6! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
7! * It is distributed WITHOUT ANY WARRANTY.
8! *
9! * Permission to use, copy, modify and distribute this software and its
10! * documentation for any purpose is hereby granted without fee,
11! * provided that the above copyright notice appear in all copies and
12! * that both that copyright notice and this permission notice appear
13! * in supporting documentation. It is provided "as is" without express
14! * or implied warranty.
15! *
16!
17!
18! Author(s): Thomas Bretz, 12/2000 <mailto:tbretz@astro.uni-wuerzburg.de>
19!
20! Copyright: MAGIC Software Development, 2000-2004
21!
22!
23\* ======================================================================== */
24
25/////////////////////////////////////////////////////////////////////////////
26//
27// MRawEvtData
28//
29// Storage container to store the raw FADC values.
30//
31// MArrayS fHiGainPixId
32// ---------------------
33// Array of Pixel Numbers for their high voltage channel in the order the
34// FADC values are stored in fHiGainFadcSamples
35//
36// MArrayB fHiGainFadcSaples
37// -------------------------
38// FADC samples (hi gain) of all pixels
39//
40// MArrayS fLoGainPixId
41// --------------------
42// see fHiGainPixId
43//
44// MArrayB fLoGainFadcSamples
45// --------------------------
46// see fHiGainFadcSaples
47//
48//
49// Version 3:
50// ----------
51// - Added fABFlags
52//
53// Version 2:
54// ----------
55// - Derives from MCamEvent now
56//
57// Version 1:
58// ----------
59// - First implementation
60//
61/////////////////////////////////////////////////////////////////////////////
62
63#include "MRawEvtData.h"
64
65#include <fstream>
66
67#include <TH1.h>
68#include <TGraph.h>
69#include <TArrayC.h>
70#include <TVirtualPad.h>
71
72#include "MLog.h"
73#include "MLogManip.h"
74
75#include "MArrayS.h"
76#include "MArrayB.h"
77#include "MGeomCam.h"
78
79#include "MRawCrateArray.h"
80#include "MRawCrateData.h"
81
82#include "MRawRunHeader.h"
83#include "MRawEvtPixelIter.h"
84
85ClassImp(MRawEvtData);
86
87using namespace std;
88
89// --------------------------------------------------------------------------
90//
91// Default constructor. It initializes all arrays with zero size.
92//
93MRawEvtData::MRawEvtData(const char *name, const char *title)
94{
95 fName = name ? name : "MRawEvtData";
96 fTitle = title ? title : "Raw Event Data Information";
97
98 InitArrays();
99}
100
101// --------------------------------------------------------------------------
102//
103// Destructor. Deletes all the arrays.
104//
105MRawEvtData::~MRawEvtData()
106{
107 DeleteArrays();
108}
109
110// --------------------------------------------------------------------------
111//
112// reset all arrays
113//
114void MRawEvtData::Clear(Option_t *)
115{
116 /*
117 FIXME:
118 Is Reset (set all entries to zero) what you want to do
119 or Set(0) (delete the array)
120 */
121 fHiGainPixId->Reset();
122 fLoGainPixId->Reset();
123 fHiGainFadcSamples->Reset();
124 fLoGainFadcSamples->Reset();
125 fABFlags->Reset();
126}
127
128// --------------------------------------------------------------------------
129//
130// return the number of hi gain samples per pixel
131//
132Byte_t MRawEvtData::GetNumHiGainSamples() const
133{
134 return fHiGainPixId->GetSize() ? fHiGainFadcSamples->GetSize()/fHiGainPixId->GetSize() : 0;
135}
136
137// --------------------------------------------------------------------------
138//
139// return the number of lo gain samples per pixel
140//
141Byte_t MRawEvtData::GetNumLoGainSamples() const
142{
143 return fLoGainPixId->GetSize() ? fLoGainFadcSamples->GetSize()/fLoGainPixId->GetSize() : 0;
144}
145
146// --------------------------------------------------------------------------
147//
148// return the number of stored pixel
149//
150UShort_t MRawEvtData::GetNumPixels() const
151{
152 return fHiGainPixId->GetSize();
153}
154
155
156// --------------------------------------------------------------------------
157//
158// Print out the onformation to *fLog.
159// Options:
160// "hex" Prints the time slices hexadecimal (default)
161// "dec" Prints the time slices decimal
162//
163void MRawEvtData::Print(Option_t *opt) const
164{
165 //
166 // print fadc inforation to screen
167 // Possible Options:
168 // - DEC: Print values decimal instead of hexadecimal (default)
169 //
170 const Byte_t nHiSamp = GetNumHiGainSamples();
171 const Byte_t nLoSamp = GetNumLoGainSamples();
172
173 const UShort_t nHiPix = fHiGainPixId->GetSize();;
174 const UShort_t nLoPix = fLoGainPixId->GetSize();;
175
176 fLog->unsetf(ios::showbase);
177
178 *fLog << dec << all;
179 *fLog << "HiGain: " << nHiPix << " Pixels with " << (Int_t)nHiSamp << " Samples" << endl;
180 *fLog << "LoGain: " << nLoPix << " Pixels with " << (Int_t)nLoSamp << " Samples";;
181
182 TString str(opt);
183 Int_t manip = str.Contains("dec", TString::kIgnoreCase);
184
185 Int_t l=0;
186 for (int i=0; i<nHiPix; i++)
187 {
188 const UShort_t idx = (*fHiGainPixId)[i];
189
190 *fLog << endl;
191 *fLog << " " << setfill(' ') << setw(3) << dec << i << " - " << setw(3);
192 *fLog << dec << idx << " ";
193
194 if (fABFlags->GetSize())
195 {
196 *fLog << "<" << hex << setfill('0') << setw(2);
197 *fLog << ((Int_t)(*fABFlags)[idx/8]&0xff) << "> ";
198 }
199
200 *fLog << (manip?dec:hex) << (manip ? setfill(' ') : setfill('0'));
201
202 for (int j=0; j<nHiSamp; j++)
203 {
204 *fLog << setw(manip?3:2);
205 *fLog << ((UShort_t)(*fHiGainFadcSamples)[j+i*nHiSamp]&0xff);
206 if (manip)
207 *fLog << ' ';
208 }
209
210 if (!(l<nLoPix && (*fLoGainPixId)[l]==idx))
211 continue;
212
213 if (manip)
214 *fLog << "/ ";
215
216 for (int j=0; j<nLoSamp; j++)
217 {
218 *fLog << setw(manip?3:2);
219 *fLog << ((UShort_t)(*fLoGainFadcSamples)[j+i*nLoSamp]&0xff);
220 if (manip)
221 *fLog << ' ';
222 }
223 l++;
224 }
225 *fLog << endl;
226}
227
228// --------------------------------------------------------------------------
229//
230// Draw a pixel. A Histogram or Graph is created and it's draw function is
231// called.
232// Options:
233// "GRAPH" A graph is drawn
234// "HIST" A histogram is drawn
235// <index> The pixel with the given index is drawn
236//
237void MRawEvtData::Draw(Option_t *opt)
238{
239 if (GetNumPixels()==0)
240 {
241 *fLog << warn << "Sorry, no pixel to draw!" << endl;
242 return;
243 }
244
245 TString str(opt);
246 str.ToLower();
247
248 UInt_t id = 0;
249
250 if (str.BeginsWith("graph"))
251 if (str.Length()>5)
252 sscanf(&str[5], "%d", &id);
253 if (str.BeginsWith("hist"))
254 if (str.Length()>4)
255 sscanf(&str[4], "%d", &id);
256
257 MRawEvtPixelIter pix(this);
258 if (!pix.Jump(id))
259 {
260 *fLog << warn << dec << "Pixel Idx #" << id << " doesn't exist!" << endl;
261 return;
262 }
263
264 const Byte_t *higains = pix.GetHiGainSamples();
265 const Byte_t *logains = pix.GetLoGainSamples();
266
267 const Int_t nh = GetNumHiGainSamples();
268 const Int_t nl = GetNumLoGainSamples();
269
270 TString name = "Pixel Idx.";
271 name += pix.GetPixelId();
272
273 Bool_t same = str.Contains("same");
274
275 if (str.BeginsWith("graph"))
276 {
277 *fLog << inf << "Drawing Graph: Pixel Idx #" << dec << pix.GetPixelId();
278 *fLog << " of " << (int)GetNumPixels() << "Pixels" << endl;
279
280 TGraph *graphhi = new TGraph;
281
282 for (int i=0; i<nh; i++)
283 graphhi->SetPoint(graphhi->GetN(), i, higains[i]);
284
285 graphhi->SetMaximum(256);
286 graphhi->SetMinimum(0);
287
288 graphhi->SetBit(kCanDelete);
289 graphhi->Draw(same ? "C*" : "AC*");
290
291 TH1F *histhi = graphhi->GetHistogram();
292 histhi->SetMinimum(0);
293 histhi->SetMaximum(255);
294
295 histhi->SetXTitle("Time/FADC Slices");
296 histhi->SetYTitle("Signal/FADC Units");
297
298 if (nl>0)
299 {
300 TGraph *graphlo = new TGraph;
301
302 for (int i=0; i<nl; i++)
303 graphlo->SetPoint(graphlo->GetN(), i, logains[i]);
304
305 graphlo->SetMaximum(256);
306 graphlo->SetMinimum(0);
307 graphlo->SetLineColor(kBlue);
308
309 graphlo->SetBit(kCanDelete);
310 graphlo->Draw("C*");
311
312 TH1F *histlo = graphlo->GetHistogram();
313
314 histlo->SetXTitle("Time/FADC Slices");
315 histlo->SetYTitle("Signal/FADC Units");
316 }
317
318 return;
319 }
320
321 if (str.BeginsWith("hist"))
322 {
323 // FIXME: Add Legend
324 *fLog << inf << "Drawing Histogram of Pixel with Idx #" << dec << pix.GetPixelId() << " to " << gPad << endl;
325
326 TH1F *histh = new TH1F(name, "FADC Samples", nh, -0.5, nh-.5);
327 histh->SetMinimum(0);
328 histh->SetMaximum(255);
329 histh->SetXTitle("Time [FADC Slices]");
330 histh->SetYTitle("Signal [FADC Units]");
331 histh->SetDirectory(NULL);
332 for (int i=0; i<nh; i++)
333 histh->Fill(i, higains[i]);
334 histh->SetBit(kCanDelete);
335 histh->Draw(same ? "same" : "");
336
337 if (nl>0)
338 {
339 TH1F *histl = new TH1F(name+";2", "FADC Samples", nl, -0.5, nl-.5);
340 histl->SetLineColor(kBlue);
341 histl->SetDirectory(NULL);
342 for (int i=0; i<nl; i++)
343 histl->Fill(i, logains[i]);
344 histl->SetBit(kCanDelete);
345 histl->Draw("same");
346 }
347 return;
348 }
349
350 *fLog << warn << dbginf << "Warning - You must specify either 'GRAPH' or 'HIST'" << endl;
351}
352
353// --------------------------------------------------------------------------
354//
355// Deletes all arrays describing the pixel Id and Samples in pixels.
356// The flag is for future usage.
357//
358void MRawEvtData::DeletePixels(Bool_t flag)
359{
360 if (fRunHeader && flag)
361 {
362 //const int npix = fRunHeader->GetNumCrates()*fRunHeader->GetNumPixInCrate();
363 const int npix = fRunHeader->GetNumConnectedPixels();
364
365 if (fArraySize == npix)
366 {
367 fPosInArray = 0;
368 fConnectedPixels = 0;
369 return;
370 }
371 }
372
373 DeleteArrays();
374 InitArrays(flag);
375}
376
377// --------------------------------------------------------------------------
378//
379// Deletes all the arrays
380//
381void MRawEvtData::DeleteArrays()
382{
383 delete fHiGainPixId;
384 delete fLoGainPixId;
385 delete fHiGainFadcSamples;
386 delete fLoGainFadcSamples;
387 delete fABFlags;
388}
389
390// --------------------------------------------------------------------------
391//
392// Deletes all the arrays
393// The flag is for future usage.
394//
395void MRawEvtData::InitArrays(Bool_t flag)
396{
397 if (flag && fRunHeader)
398 {
399 //const int npix = fRunHeader->GetNumCrates()*fRunHeader->GetNumPixInCrate();
400 const int npix = fRunHeader->GetNumConnectedPixels();
401
402 fHiGainPixId = new MArrayS(npix);
403 fLoGainPixId = new MArrayS(npix);
404 fHiGainFadcSamples = new MArrayB(npix*fRunHeader->GetNumSamplesHiGain());
405 fLoGainFadcSamples = new MArrayB(npix*fRunHeader->GetNumSamplesLoGain());
406 fABFlags = new TArrayC(fRunHeader->GetMaxPixId());
407
408 fArraySize = npix;
409 }
410 else
411 {
412 fHiGainPixId = new MArrayS(0);
413 fLoGainPixId = new MArrayS(0);
414 fHiGainFadcSamples = new MArrayB(0);
415 fLoGainFadcSamples = new MArrayB(0);
416 fABFlags = new TArrayC(0);
417
418 fArraySize = 0;
419 }
420
421 fPosInArray = 0;
422 fConnectedPixels = 0;
423}
424
425// --------------------------------------------------------------------------
426//
427// This is to fill the data of one pixel to the MRawEvtHeader Class.
428// The parameters are the pixelnumber and the FADC_SLICES values of ADCs
429// Add to lo gains if lflag = 1
430//
431void MRawEvtData::AddPixel(UShort_t nOfPixel, TArrayC *data, Bool_t lflag)
432{
433 MArrayS *arrpix = lflag ? fLoGainPixId : fHiGainPixId;
434 MArrayB *arrsam = lflag ? fLoGainFadcSamples : fHiGainFadcSamples;
435
436 //
437 // check whether we got the right number of new samples
438 // if there are no samples already stored: this is the new number of samples
439 //
440 const Byte_t ns = data->GetSize();
441 const Byte_t nSamp = lflag ? GetNumLoGainSamples() : GetNumHiGainSamples();
442 if (nSamp && ns!=nSamp)
443 {
444 *fLog << err << "RawEvtData::AddPixel: Error, number of samples in ";
445 *fLog << "TArrayC " << ns << " doesn't match current number " << nSamp << endl;
446 return;
447 }
448
449 //
450 // enhance pixel array by one
451 //
452 arrpix->Set(arrpix->GetSize()+1);
453
454 //
455 // add the number of the new pixel to the array as last entry
456 //
457 arrpix->AddAt(nOfPixel, arrpix->GetSize()-1);
458
459 //
460 // enhance the array by the number of new samples
461 //
462 arrsam->Set(arrsam->GetSize()+ns);
463
464 //
465 // add the new slices as last entries to array
466 //
467 arrsam->AddAt((Byte_t*)data->GetArray(), arrsam->GetSize()-ns, ns);
468}
469
470// --------------------------------------------------------------------------
471//
472// Fills members with information from a magic binary file.
473// WARNING: you have to use Init() before you can do this
474//
475void MRawEvtData::ReadEvt(istream &fin)
476{
477 const UShort_t nlo = fRunHeader->GetNumSamplesLoGain();
478 const UShort_t nhi = fRunHeader->GetNumSamplesHiGain();
479
480 const UShort_t npic = fRunHeader->GetNumPixInCrate();
481
482 const UShort_t npos = npic*fPosInArray;
483
484 const Byte_t ab = fCrateArray->GetEntry(fPosInArray)->GetABFlags();
485
486 Byte_t *higainsam = fHiGainFadcSamples->GetArray()+nhi*fConnectedPixels;
487 Byte_t *logainsam = fLoGainFadcSamples->GetArray()+nlo*fConnectedPixels;
488
489 for (int i=0; i<npic; i++)
490 {
491 fin.read((char*)higainsam, nhi);
492 fin.read((char*)logainsam, nlo);
493
494 // calc the spiral hardware pixel number
495 const UShort_t ipos = npos+i;
496
497 // -1 converts the hardware pixel Id into the software pixel index
498 const Int_t npix = (Int_t)fRunHeader->GetPixAssignment(ipos)-1;
499
500 // Check whether the pixel is connected or not
501 if (npix<0)
502 continue;
503
504 //
505 // This is to fill the data of one pixel to the MRawEvtHeader Class.
506 // The parameters are the pixelnumber and the FADC_SLICES values of ADCs
507 // Add to lo gains if lflag = 1
508 //
509 fHiGainPixId->AddAt(npix, fConnectedPixels);
510 higainsam += nhi;
511
512 // FIXME: Not implemented in the raw files yet
513 //if (IsLoGainOn(i, j))
514 //{
515 fLoGainPixId->AddAt(npix, fConnectedPixels);
516 logainsam += nlo;
517 //}
518
519 if (TESTBIT(ab, i))
520 SETBIT((*fABFlags)[npix/8], npix%8);
521 else
522 CLRBIT((*fABFlags)[npix/8], npix%8);
523
524 fConnectedPixels++;
525 }
526
527 fPosInArray++;
528}
529
530// --------------------------------------------------------------------------
531//
532// Return the size in bytes of one event data block
533//
534Int_t MRawEvtData::GetNumBytes() const
535{
536 const UShort_t nlo = fRunHeader->GetNumSamplesLoGain();
537 const UShort_t nhi = fRunHeader->GetNumSamplesHiGain();
538 const UShort_t npic = fRunHeader->GetNumPixInCrate();
539
540 return (nhi+nlo)*npic;
541}
542
543// --------------------------------------------------------------------------
544//
545// Make sure, that you skip the whole event. This function only skips a part
546// of the event - see MRawRead::SkipEvent
547//
548void MRawEvtData::SkipEvt(istream &fin)
549{
550 fin.seekg(GetNumBytes(), ios::cur);
551}
552
553Bool_t MRawEvtData::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const
554{
555 MRawEvtPixelIter Next(const_cast<MRawEvtData*>(this));
556 if (!Next.Jump(idx))
557 return kFALSE;
558
559 switch (type)
560 {
561 case 0:
562 val = Next.GetSumHiGainSamples()-(float)GetNumHiGainSamples()*fHiGainFadcSamples->GetArray()[0];
563 val*= cam.GetPixRatio(idx);
564 break;
565 case 1:
566 val = Next.GetMaxHiGainSample();
567 break;
568 case 2:
569 val = Next.GetMaxLoGainSample();
570 break;
571 case 3:
572 val = Next.GetIdxMaxHiGainSample();
573 break;
574 case 4:
575 val = Next.GetIdxMaxLoGainSample();
576 break;
577 case 5:
578 val = Next.GetIdxMaxHiLoGainSample();
579 return val >= 0;
580 }
581
582 return kTRUE;
583}
584
585void MRawEvtData::Copy(TObject &named)
586#if ROOT_VERSION_CODE > ROOT_VERSION(3,04,01)
587const
588#endif
589{
590 MRawEvtData &evt = (MRawEvtData &)named;
591
592 *evt.fHiGainPixId = *fHiGainPixId;
593 *evt.fLoGainPixId = *fLoGainPixId;
594
595 *evt.fHiGainFadcSamples = *fHiGainFadcSamples;
596 *evt.fLoGainFadcSamples = *fLoGainFadcSamples;
597
598 *evt.fABFlags = *fABFlags;
599
600 evt.fPosInArray = fPosInArray;
601 evt.fConnectedPixels = fConnectedPixels;
602 evt.fArraySize = fArraySize;
603}
Note: See TracBrowser for help on using the repository browser.