FMD/scripts/MakeSpectra.C

   1 #ifdef SPECTRA_BUILD
   2 #include <TAxis.h>
   3 #include <TNamed.h>
   4 #include <TH1.h>
   5 #include <TH2.h>
   6 #include <TCanvas.h>
   7 #include <TObjArray.h>
   8 #include <TDirectory.h>
   9 #include <TBrowser.h>
  10 #include <TFile.h>
  11
  12 namespace Spectra
  13 {
  14   //__________________________________________________________________
  15   /**
  16    * An Element
  17    *
  18    */
  19   struct Element : public TNamed
  20   {
  21     /**
  22      * Destructor
  23      */
  24     virtual ~Element()
  25     {
  26       if (fFull) delete fFull;
  27       if (fFull) delete fCut;
  28       fChildren.Delete();
  29     }
  30     /**
  31      * Draw this element.  Draws the full histogram and the selected
  32      * range on top.
  33      *
  34      * @param option Drawing option
  35      * @param l      Low cut
  36      * @param h      High cut
  37      */
  38     virtual void DrawIt(Option_t* option="", Double_t l=-1, Double_t h=-1) //*MENU*
  39     {
  40       if (!fFull) return;
  41       if (fFull->GetMinimum() < 0) gPad->SetLogy(kFALSE);
  42       fFull->Draw(option);
  43       if (!fCut || l == h) return;
  44       Double_t lx = fFull->GetXaxis()->GetXmin();
  45       Double_t hx = fFull->GetXaxis()->GetXmax();
  46       Double_t rr = (hx-lx);
  47       Double_t ll = rr * l + lx;
  48       Double_t hh = rr * h + lx;
  49       for (Int_t i = 1; i <= fFull->GetNbinsX(); i++) {
  50         if (fFull->GetBinCenter(i) <= ll ||
  51             fFull->GetBinCenter(i) >  hh) {
  52           fCut->SetBinContent(i, 0);
  53           continue;
  54         }
  55         fCut->SetBinContent(i, fFull->GetBinContent(i));
  56       }
  57       fCut->Draw(Form("%s same", option));
  58       gPad->Modified();
  59       gPad->Update();
  60       gPad->cd();
  61     }//*MENU*
  62     /**
  63      * Draw
  64      *
  65      * @param option
  66      */
  67     virtual void Draw(Option_t* option="")
  68     {
  69       if (!fFull) return;
  70       if (fFull->GetMinimum() < 0) gPad->SetLogy(kFALSE);
  71       fFull->Draw(option);
  72       gPad->Modified();
  73       gPad->Update();
  74       gPad->cd();
  75     } //*MENU*
  76     /**
  77      * Get the top-level node
  78      *
  79      * @return Top-level node
  80      */
  81     virtual Element& GetParent() { return *fParent; }
  82     /**
  83      * Make the histograms
  84      *
  85      * @param axis Axis to use
  86      */
  87     virtual void MakeHistograms(const TAxis& axis)
  88     {
  89       if (fFull) return;
  90       if (axis.GetNbins() <= 1) return;
  91       if (axis.IsVariableBinSize())
  92         fFull = new TH1F(Form("f_%s", GetName()), Form("%s spectra", GetName()),
  93                          axis.GetXbins()->fN, axis.GetXbins()->fArray);
  94       else
  95         fFull = new TH1F(Form("f_%s", GetName()), Form("%s spectra", GetName()),
  96                          axis.GetNbins(), axis.GetXmin(), axis.GetXmax());
  97       fCut = static_cast<TH1F*>(fFull->Clone(Form("c_%s", GetName())));
  98       fCut->SetTitle(Form("%s restricted spectra", GetName()));
  99       fFull->SetDirectory(0);
 100       fFull->SetFillColor(kRed);
 101       fFull->SetFillStyle(3001);
 102       fCut->SetDirectory(0);
 103       fCut->SetFillColor(kBlue);
 104       fCut->SetFillStyle(3001);
 105     }
 106     /**
 107      * Get Children
 108      *
 109      * @return Array of children
 110      */
 111     TObjArray& Children() { return fChildren; }
 112     /**
 113      * Get Children
 114      *
 115      * @return Array of children
 116      */
 117     const TObjArray& Children() const { return fChildren; }
 118     /**
 119      * Write to a directory
 120      *
 121      * @param d Directory to write to
 122      */
 123     void WriteOut(TDirectory* d)
 124     {
 125       TDirectory* dd = d->mkdir(GetName());
 126       dd->cd();
 127       if (fFull) fFull->Write();
 128       if (fCut)  fCut->Write();
 129
 130       TIter next(&fChildren);
 131       Element* e = 0;
 132       while ((e = static_cast<Element*>(next())))
 133         e->WriteOut(dd);
 134       d->cd();
 135     }
 136     /**
 137      * This is a folder
 138      *
 139      * @return Always true
 140      */
 141     virtual Bool_t IsFolder() const { return kTRUE; }
 142     /**
 143      * Browse this element
 144      *
 145      * @param b Brower to use
 146      */
 147     virtual void Browse(TBrowser* b)
 148     {
 149       if (fFull) b->Add(fFull);
 150       if (fCut)  b->Add(fCut);
 151
 152       TIter next(&fChildren);
 153       TObject* o = 0;
 154       while ((o = next())) b->Add(o);
 155       // if (fChildren.GetEntriesFast() > 0) b->Add(&fChildren);
 156     }
 157   protected:
 158     Element() : TNamed(), fFull(0), fCut(0), fParent(0), fChildren(0) {}
 159     /**
 160      * Constructor
 161      *
 162      * @param name   Name
 163      * @param title  Title
 164      */
 165     Element(const char* name, const char* title,  Element& parent)
 166       : TNamed(name, title), fFull(0), fCut(0),
 167         fParent(&parent), fChildren(0)
 168     {
 169       fChildren.SetOwner();
 170       fChildren.SetName("children");
 171     }
 172     /**
 173      * Fill value
 174      *
 175      * @param v Value to fill
 176      */
 177     void DoFill(Double_t v) { if (fFull) fFull->Fill(v); }
 178     /**
 179      * Get a child or null
 180      *
 181      * @param id Id of child
 182      *
 183      * @return Pointer to child, or null
 184      */
 185     Element* GetChild(UShort_t id) const
 186     {
 187       if (id >= fChildren.GetEntriesFast()) return 0;
 188       return static_cast<Element*>(fChildren.At(id));
 189     }
 190     TH1*      fFull; // Full histogram
 191     TH1*      fCut;  // Selected data
 192     Element*  fParent;
 193     TObjArray fChildren;
 194
 195     ClassDef(Element,1);
 196   };
 197
 198   //__________________________________________________________________
 199   struct Detector;
 200   /**
 201    * Top-level node
 202    *
 203    */
 204   struct Top : public Element
 205   {
 206     /**
 207      * Constructor
 208      */
 209     Top() : Element("top", "Top", *this), fAxis(1,0,1) {}
 210     /**
 211      * Get or add a detector
 212      *
 213      * @param d Detector (1-3)
 214      *
 215      * @return Reference to detector
 216      */
 217     Detector& GetOrAdd(UShort_t d);
 218     /**
 219      * Fill in values
 220      *
 221      * @param d  Detector
 222      * @param r  Ring
 223      * @param s  Sector
 224      * @param t  Strip
 225      * @param v  Value
 226      */
 227     void Fill(UShort_t d, Char_t r, UShort_t s, UShort_t t, Double_t v);
 228     /**
 229      * Get the axis to use
 230      *
 231      * @return Reference to axis
 232      */
 233     const TAxis& GetAxis() const { return fAxis; }
 234     /**
 235      * Set the axis to use
 236      *
 237      * @param a Axis used
 238      */
 239     void SetAxis(const TAxis& a)
 240     {
 241       fAxis.Set(a.GetNbins(),a.GetXmin(),a.GetXmax());
 242     }
 243     void SetAxis(Int_t n, Double_t l, Double_t h) { fAxis.Set(n, l, h); }
 244   protected:
 245     TAxis fAxis; // The axis
 246
 247     ClassDef(Top,1);
 248   };
 249
 250   //__________________________________________________________________
 251   struct Ring;
 252   /**
 253    * Detector
 254    *
 255    */
 256   struct Detector : public Element
 257   {
 258     Detector() : Element(), fId(0) {}
 259     /**
 260      * Constrictor
 261      *
 262      * @param d    Id
 263      * @param top  PArent node
 264      * @param a    Axis
 265      */
 266     Detector(UShort_t d, Top& top)
 267       : Element(Form("FMD%d", d), Form("FMD%d",d), top), fId(d)
 268     {}
 269     /**
 270      * Get ID
 271      *
 272      * @return The ID
 273      */
 274     UShort_t Id() const { return fId; }
 275     /**
 276      * Get or add a sub element
 277      *
 278      * @param id Id of sub-element
 279      *
 280      * @return Sub element
 281      */
 282     Ring& GetOrAdd(Char_t id);
 283     /**
 284      * Fill in value
 285      *
 286      * @param r  Ring
 287      * @param s  Sector
 288      * @param t  Strip
 289      * @param v  Value
 290      */
 291     void Fill(Char_t r, UShort_t s, UShort_t t, Double_t v);
 292
 293     /**
 294      * Get Top
 295      *
 296      * @return top
 297      */
 298     Top& M() { return static_cast<Top&>(*fParent); }
 299     /**
 300      * Get Top
 301      *
 302      * @return top
 303      */
 304     const Top& M() const { return static_cast<Top&>(*fParent); }
 305   protected:
 306     UShort_t R2Id(Char_t r) { return (r == 'I' || r == 'i') ? 0 : 1; }
 307     UShort_t fId;
 308
 309     ClassDef(Detector,1);
 310   };
 311
 312   //__________________________________________________________________
 313   struct Sector;
 314   /**
 315    * A ring
 316    *
 317    */
 318   struct Ring : public Element
 319   {
 320     Ring() : Element(), fId('\0') {}
 321     Ring(UShort_t r, Detector& d)
 322       : Element(Form("%s%c", d.GetName(), r),
 323                 Form("%s%c", d.GetName(), r), d), fId(r)
 324     {}
 325     /**
 326      * Get ID
 327      *
 328      * @return The ID
 329      */
 330     Char_t Id() const { return fId; }
 331     /**
 332      * Get or add a sub element
 333      *
 334      * @param id Id of sub-element
 335      *
 336      * @return Sub element
 337      */
 338     Sector& GetOrAdd(UShort_t id);
 339     /**
 340      * Fill in value
 341      *
 342      * @param s  Sector
 343      * @param t  Strip
 344      * @param v  Value
 345      */
 346     void Fill(UShort_t s, UShort_t t, Double_t v);
 347
 348     /**
 349      * Get parent detector
 350      *
 351      * @return Parent detector
 352      */
 353     Detector& D() { return static_cast<Detector&>(*fParent); }
 354     /**
 355      * Get parent detector
 356      *
 357      * @return Parent detector
 358      */
 359     const Detector& D() const { return static_cast<Detector&>(*fParent); }
 360     /**
 361      * Get Top
 362      *
 363      * @return top
 364      */
 365     Top& M() { return D().M(); }
 366     /**
 367      * Get Top
 368      *
 369      * @return top
 370      */
 371     const Top& M() const { return D().M(); }
 372
 373     /**
 374      * Draw
 375      *
 376      * @param option
 377      */
 378     virtual void Draw(Option_t* option="lego2z")
 379     {
 380       if (!fFull) return;
 381       gPad->SetLogy(fFull->GetMinimum() > 0);
 382       fFull->Draw(option);
 383       gPad->cd();
 384     } //*MENU*
 385     void MakeHistograms(const TAxis& axis)
 386     {
 387       if (fFull) return;
 388       if (axis.GetNbins() <= 1) return;
 389       Int_t nSec = (fId == 'I' || fId == 'i' ? 20 : 40);
 390       if (axis.IsVariableBinSize())
 391         fFull = new TH2F(Form("f_%s", GetName()), Form("%s spectra", GetName()),
 392                          nSec, -.5, nSec-.5,
 393                          axis.GetXbins()->fN, axis.GetXbins()->fArray);
 394       else
 395         fFull = new TH2F(Form("f_%s", GetName()), Form("%s spectra", GetName()),
 396                          nSec, -.5, nSec-.5,
 397                          axis.GetNbins(), axis.GetXmin(), axis.GetXmax());
 398     }
 399   protected:
 400     Char_t fId;
 401
 402     ClassDef(Ring,1);
 403   };
 404
 405
 406   //__________________________________________________________________
 407   struct Strip;
 408   /**
 409    * A sector
 410    *
 411    */
 412   struct Sector : public Element
 413   {
 414     Sector() : Element(), fId(1024) {}
 415     Sector(UShort_t s, Ring& r)
 416       : Element(Form("%s[%02d]", r.GetName(), s),
 417                 Form("%s[%02d]", r.GetName(), s), r), fId(s)
 418     {
 419     }
 420     /**
 421      * Get ID
 422      *
 423      * @return The ID
 424      */
 425     UShort_t Id() const { return fId; }
 426     /**
 427      * Get or add a sub element
 428      *
 429      * @param id Id of sub-element
 430      *
 431      * @return Sub element
 432      */
 433     Strip& GetOrAdd(UShort_t id);
 434     /**
 435      * Fill in value
 436      *
 437      * @param t  Strip
 438      * @param v  Value
 439      */
 440     void Fill(UShort_t t, Double_t v);
 441
 442     /**
 443      * Draw
 444      *
 445      * @param option
 446      */
 447     virtual void Draw(Option_t* option="lego2z")
 448     {
 449       if (!fFull) return;
 450       gPad->SetLogy(fFull->GetMinimum() > 0);
 451       fFull->Draw(option);
 452       gPad->cd();
 453     } //*MENU*
 454     void MakeHistograms(const TAxis& axis)
 455     {
 456       if (fFull) return;
 457       if (axis.GetNbins() <= 1) return;
 458       Int_t nStr = (R().Id() == 'I' || fId == 'i' ? 512 : 256);
 459       if (axis.IsVariableBinSize())
 460         fFull = new TH2F(Form("f_%s", GetName()), Form("%s spectra", GetName()),
 461                          nStr, -.5, nStr-.5,
 462                          axis.GetXbins()->fN, axis.GetXbins()->fArray);
 463       else
 464         fFull = new TH2F(Form("f_%s", GetName()), Form("%s spectra", GetName()),
 465                          nStr, -.5, nStr-.5,
 466                          axis.GetNbins(), axis.GetXmin(), axis.GetXmax());
 467     }
 468
 469     /**
 470      * Get parent ring
 471      *
 472      * @return Parent ring
 473      */
 474     Ring& R() { return static_cast<Ring&>(*fParent); }
 475     /**
 476      * Get parent ring
 477      *
 478      * @return Parent ring
 479      */
 480     const Ring& R() const { return static_cast<Ring&>(*fParent); }
 481     /**
 482      * Get parent detector
 483      *
 484      * @return Parent detector
 485      */
 486     Detector& D() { return R().D(); }
 487     /**
 488      * Get parent detector
 489      *
 490      * @return Parent detector
 491      */
 492     const Detector& D() const { return R().D(); }
 493     /**
 494      * Get Top
 495      *
 496      * @return top
 497      */
 498     Top& M() { return D().M(); }
 499     /**
 500      * Get Top
 501      *
 502      * @return top
 503      */
 504     const Top& M() const { return D().M(); }
 505   protected:
 506     UShort_t fId;
 507
 508     ClassDef(Sector,1);
 509   };
 510
 511   //__________________________________________________________________
 512   /**
 513    * A stripo
 514    *
 515    */
 516   struct Strip : public Element
 517   {
 518     Strip() : Element(), fId(1024) {}
 519     Strip(UShort_t t, Sector& s)
 520       : Element(Form("%s[%03d]", s.GetName(), t),
 521                 Form("%s[%03d]", s.GetName(), t), s), fId(t) {}
 522     /**
 523      * Get ID
 524      *
 525      * @return The ID
 526      */
 527     UShort_t Id() const { return fId; }
 528     /**
 529      * Fill in value
 530      *
 531      * @param v  Value
 532      */
 533     void Fill(Double_t v)
 534     {
 535       // Info("Fill", "%s: Filling with %f", GetName(), v);
 536       DoFill(v);
 537     }
 538     Bool_t IsFolder() const { return 0; }
 539     void Browse(TBrowser* b) { Draw(b->GetDrawOption()); }
 540     /**
 541      * Get parent sector
 542      *
 543      * @return Parent sector
 544      */
 545     Sector& S() { return static_cast<Sector&>(*fParent); }
 546     /**
 547      * Get parent sector
 548      *
 549      * @return Parent sector
 550      */
 551     const Sector& S() const { return static_cast<Sector&>(*fParent); }
 552     /**
 553      * Get parent ring
 554      *
 555      * @return Parent ring
 556      */
 557     Ring& R() { return S().R(); }
 558     /**
 559      * Get parent ring
 560      *
 561      * @return Parent ring
 562      */
 563     const Ring& R() const { return S().R(); }
 564     /**
 565      * Get parent detector
 566      *
 567      * @return Parent detector
 568      */
 569     Detector& D() { return R().D(); }
 570     /**
 571      * Get parent detector
 572      *
 573      * @return Parent detector
 574      */
 575     const Detector& D() const { return R().D(); }
 576     /**
 577      * Get Top
 578      *
 579      * @return top
 580      */
 581     Top& M() { return D().M(); }
 582     /**
 583      * Get Top
 584      *
 585      * @return top
 586      */
 587     const Top& M() const { return D().M(); }
 588   protected:
 589     UShort_t fId;
 590
 591     ClassDef(Strip,1);
 592   };
 593
 594   //==================================================================
 595   Strip& Sector::GetOrAdd(UShort_t id)
 596   {
 597     Strip* t = static_cast<Strip*>(GetChild(id));
 598     if (!t) {
 599       t = new Strip(id, *this);
 600       t->MakeHistograms(M().GetAxis());
 601       // Info("GetOrAdd", "%s: Adding strip @ %d", GetName(), id);
 602       fChildren.AddAtAndExpand(t, id);
 603     }
 604     return *t;
 605   }
 606   //__________________________________________________________________
 607   Sector& Ring::GetOrAdd(UShort_t id)
 608   {
 609     Sector* s = static_cast<Sector*>(GetChild(id));
 610     if (!s) {
 611       s = new Sector(id, *this);
 612       s->MakeHistograms(M().GetAxis());
 613       // Info("GetOrAdd", "%s: Adding sector @ %d", GetName(), id);
 614       fChildren.AddAtAndExpand(s, id);
 615     }
 616     return *s;
 617   }
 618   //__________________________________________________________________
 619   Ring& Detector::GetOrAdd(Char_t id)
 620   {
 621     UShort_t idd = R2Id(id);
 622     Ring*    r   = static_cast<Ring*>(GetChild(idd));
 623     if (!r) {
 624       r = new Ring(id, *this);
 625       r->MakeHistograms(M().GetAxis());
 626       // Info("GetOrAdd", "%s: Adding ring @ %d", GetName(), idd);
 627       fChildren.AddAtAndExpand(r, idd);
 628     }
 629     return *r;
 630   }
 631   //__________________________________________________________________
 632   Detector& Top::GetOrAdd(UShort_t id)
 633   {
 634     Int_t     idx = id - 1;
 635     Detector* d   = static_cast<Detector*>(fChildren.At(idx));
 636     if (!d) {
 637       d = new Detector(id, *this);
 638       d->MakeHistograms(fAxis);
 639       // Info("GetOrAdd", "%s: Adding detector @ %d", GetName(), id);
 640       fChildren.AddAtAndExpand(d, idx);
 641     }
 642     return *d;
 643   }
 644   //__________________________________________________________________
 645   void Top::Fill(UShort_t d, Char_t r, UShort_t s, UShort_t t, Double_t v)
 646   {
 647     DoFill(v);
 648     Detector& sub = GetOrAdd(d);
 649     //Info("Fill", "%s: Filling %d,%c,%d,%d with %f", GetName(), d, r, s, t, v);
 650     sub.Fill(r, s, t, v);
 651   }
 652   //__________________________________________________________________
 653   void Detector::Fill(Char_t r, UShort_t s, UShort_t t, Double_t v)
 654   {
 655     DoFill(v);
 656     Ring& sub = GetOrAdd(r);
 657     //Info("Fill", "%s: Filling %c,%d,%d with %f", GetName(), r, s, t, v);
 658     sub.Fill(s, t, v);
 659   }
 660   //__________________________________________________________________
 661   void Ring::Fill(UShort_t s, UShort_t t, Double_t v)
 662   {
 663     fFull->Fill(s, v);
 664     Sector& sub = GetOrAdd(s);
 665     //Info("Fill", "%s: Filling %d,%d with %f", GetName(), s, t, v);
 666     sub.Fill(t, v);
 667   }
 668   //__________________________________________________________________
 669   void Sector::Fill(UShort_t t, Double_t v)
 670   {
 671     fFull->Fill(t, v);
 672     Strip& sub = GetOrAdd(t);
 673     //Info("Fill", "%s: Filling %d with %f", GetName(), t, v);
 674     sub.Fill(v);
 675   }
 676 }
 677 //======================================================================
 678 #include <AliFMDDigit.h>
 679 #include <AliFMDSDigit.h>
 680 #include <AliFMDRecPoint.h>
 681 #include <AliFMDHit.h>
 682 #include <AliESDFMD.h>
 683 #include <AliFMDInput.h>
 684 #include <AliFMDRawReader.h>
 685 #include <AliFMDParameters.h>
 686
 687 struct SpectraCollector : public AliFMDInput
 688 {
 689   //__________________________________________________________________
 690   SpectraCollector(UShort_t what, const char* src=0)
 691     : AliFMDInput("galice.root")
 692   {
 693     Setup(what, src);
 694   }
 695   //__________________________________________________________________
 696   SpectraCollector(const char* what, const char* src=0)
 697     : AliFMDInput("galice.root")
 698   {
 699     Setup(ParseLoad(what), src);
 700   }
 701   //__________________________________________________________________
 702   void Setup(UShort_t what, const char* src)
 703   {
 704     switch (what) {
 705     case kHits:         AddLoad(kHits);         break;
 706     case kDigits:       AddLoad(kDigits);       break;
 707     case kSDigits:      AddLoad(kSDigits);      break;
 708     case kRaw:          AddLoad(kRaw);          break;
 709     case kRawCalib:     AddLoad(kRawCalib);     break;
 710     case kRecPoints:    AddLoad(kRecPoints);    break;
 711     case kESD:          AddLoad(kESD);          break;
 712     default:
 713       Fatal("Spectra", "Unknown type %d reguested", what);
 714       return;
 715     }
 716     if ((what == kRaw || what == kRawCalib)) {
 717       if (!src || src[0] == '\0')
 718         Fatal("Spectra", "No raw input specified!");
 719       SetRawFile(src);
 720     }
 721
 722     switch (what) {
 723     case kHits:     fTop.SetAxis(500, 0, 1000); break;
 724     case kDigits:   // Fall-through
 725     case kSDigits:  // Fall-through
 726     case kRaw:      // Fall-through
 727     case kRawCalib: // This is the same for all kinds of digits
 728       fTop.SetAxis(1024, -.5, 1023.5);
 729       break;
 730     case kRecPoints:  // Fall-through
 731     case kESD:        // This is the same for ESD and rec-points
 732       fTop.SetAxis(200, -.05, 19.95);
 733       break;
 734     }
 735   }
 736   //____________________________________________________________________
 737   Bool_t ProcessHit(AliFMDHit* h, TParticle* /* p */)
 738   {
 739     if (h)
 740       fTop.Fill(h->Detector(), h->Ring(), h->Sector(), h->Strip(), h->Edep());
 741     return kTRUE;
 742   }
 743   //__________________________________________________________________
 744   Bool_t ProcessSDigit(AliFMDSDigit* d)
 745   {
 746     if (d)
 747       fTop.Fill(d->Detector(), d->Ring(), d->Sector(), d->Strip(), d->Counts());
 748     return kTRUE;
 749   }
 750   //__________________________________________________________________
 751   Bool_t ProcessRawDigit(AliFMDDigit* digit)
 752   {
 753     return ProcessDigit(digit);
 754   }
 755   //__________________________________________________________________
 756   Bool_t ProcessDigit(AliFMDDigit* d)
 757   {
 758     if (d)
 759       fTop.Fill(d->Detector(), d->Ring(), d->Sector(),d->Strip(), d->Counts());
 760     return kTRUE;
 761   }
 762   //__________________________________________________________________
 763   Bool_t ProcessRawCalibDigit(AliFMDDigit* digit)
 764   {
 765     AliFMDParameters* parm = AliFMDParameters::Instance();
 766     UShort_t d          =  digit->Detector();
 767     Char_t   r          =  digit->Ring();
 768     UShort_t s          =  digit->Sector();
 769     UShort_t t          =  digit->Strip();
 770     Double_t g          =  parm->GetPulseGain(d, r, s, t);
 771     Double_t p          =  parm->GetPedestal(d, r, s, t);
 772     Double_t w          =  parm->GetPedestalWidth(d, r, s, t);
 773     UShort_t c          =  digit->Counts();
 774     Double_t x          =  0;
 775     if (fFMDReader && fFMDReader->IsZeroSuppressed(d-1))
 776       x = c + fFMDReader->NoiseFactor(d-1) * w;
 777     else
 778       x = c - p;
 779
 780     Double_t m = x / (g * parm->GetDACPerMIP());
 781     if (g < 0.1 || g > 10) m = 0;
 782
 783     fTop.Fill(d, r, s, t, m);
 784     return kTRUE;
 785   }
 786   //__________________________________________________________________
 787   Bool_t ProcessRecPoint(AliFMDRecPoint* r)
 788   {
 789     if (r)
 790       fTop.Fill(r->Detector(),r->Ring(),r->Sector(),r->Strip(),r->Particles());
 791     return kTRUE;
 792   }
 793   //__________________________________________________________________
 794   Bool_t ProcessESD(UShort_t d, Char_t r, UShort_t s, UShort_t t,
 795                     Float_t, Float_t m)
 796   {
 797     if (m != AliESDFMD::kInvalidMult)
 798       fTop.Fill(d,r,s,t,m);
 799     return kTRUE;
 800   }
 801
 802   //__________________________________________________________________
 803   Bool_t Finish()
 804   {
 805     TBrowser* b = new TBrowser("b", "b");
 806     b->Add(&fTop);
 807
 808 #if 1
 809     TFile* file = TFile::Open("spectra.root", "RECREATE");
 810     // fTop.WriteOut(file);
 811     fTop.Write();
 812     file->Close();
 813 #endif
 814
 815     return kTRUE;
 816   }
 817   Spectra::Top fTop;
 818   ClassDef(SpectraCollector,0);
 819 };
 820
 821 #else
 822 //======================================================================
 823 void
 824 MakeSpectra(UShort_t what, Int_t n=1000, const char* src=0)
 825 {
 826   gROOT->LoadMacro("$ALICE_ROOT.trunk/FMD/scripts/Compile.C");
 827   gSystem->AddIncludePath("-DSPECTRA_BUILD");
 828   const char* script = "$ALICE_ROOT.trunk/FMD/scripts/MakeSpectra.C";
 829   const char* here   = gSystem->BaseName(script);
 830   Int_t ret = 0;
 831   if ((ret = gSystem->CopyFile(gSystem->ExpandPathName(script), here, true))) {
 832     Error("MakeSpectra", "Failed to copy %s to %s: %d", script, here, ret);
 833     return;
 834   }
 835   Compile(here, "+g");
 836   SpectraCollector* sd = new SpectraCollector(what, src);
 837   sd->Run(n);
 838 }
 839
 840 #endif // SPECTRA_BUILD
 841 //
 842 // EOF
 843 //