# Example of ROOT I/O with a TTree

    // @(#)root/test:$Name:  $:$Id: Event.cxx,v 1.27 2004/01/25 20:33:32 brun Exp $
    // Author: Rene Brun   19/08/96

    ////////////////////////////////////////////////////////////////////////
    //
    //                       Event and Track classes
    //                       =======================
    //
    //  The Event class is a naive/simple example of an event structure.
    //     public:
    //        char           fType[20];
    //        char          *fEventName;         //run+event number in character format
    //        Int_t          fNtrack;
    //        Int_t          fNseg;
    //        Int_t          fNvertex;
    //        UInt_t         fFlag;
    //        Double32_t     fTemperature;
    //        Int_t          fMeasures[10];
    //        Double32_t     fMatrix[4][4];
    //        Double32_t    *fClosestDistance; //[fNvertex] indexed array!
    //        EventHeader    fEvtHdr;
    //        TClonesArray  *fTracks;
    //        TRefArray     *fHighPt;            //array of High Pt tracks only
    //        TRefArray     *fMuons;             //array of Muon tracks only
    //        TRef           fLastTrack;         //pointer to last track
    //        TRef           fHistoWeb;          //EXEC:GetHistoWeb reference to an histogram in a TWebFile
    //        TH1F          *fH;
    //        TBits          fTriggerBits;       //Bits triggered by this event.
    //
    //   The EventHeader class has 3 data members (integers):
    //     public:
    //        Int_t          fEvtNum;
    //        Int_t          fRun;
    //        Int_t          fDate;
    //
    //
    //   The Event data member fTracks is a pointer to a TClonesArray.
    //   It is an array of a variable number of tracks per event.
    //   Each element of the array is an object of class Track with the members:
    //     private:
    //        Float_t      fPx;           //X component of the momentum
    //        Float_t      fPy;           //Y component of the momentum
    //        Float_t      fPz;           //Z component of the momentum
    //        Float_t      fRandom;       //A random track quantity
    //        Float_t      fMass2;        //The mass square of this particle
    //        Float_t      fBx;           //X intercept at the vertex
    //        Float_t      fBy;           //Y intercept at the vertex
    //        Float_t      fMeanCharge;   //Mean charge deposition of all hits of this track
    //        Float_t      fXfirst;       //X coordinate of the first point
    //        Float_t      fXlast;        //X coordinate of the last point
    //        Float_t      fYfirst;       //Y coordinate of the first point
    //        Float_t      fYlast;        //Y coordinate of the last point
    //        Float_t      fZfirst;       //Z coordinate of the first point
    //        Float_t      fZlast;        //Z coordinate of the last point
    //        Double32_t   fCharge;       //Charge of this track
    //        Double32_t   fVertex[3];    //Track vertex position
    //        Int_t        fNpoint;       //Number of points for this track
    //        Short_t      fValid;        //Validity criterion
    //        Int_t        fNsp;          //Number of points for this track with a special value
    //        Double32_t  *fPointValue;   //[fNsp] a special quantity for some point.
    //        TBits        fTriggerBits;  //Bits triggered by this track.
    //
    //   An example of a batch program to use the Event/Track classes is given
    //   in this directory: MainEvent.
    //   Look also in the same directory at the following macros:
    //     - eventa.C  an example how to read the tree
    //     - eventb.C  how to read events conditionally
    //
    //   During the processing of the event (optionally) also a large number
    //   of histograms can be filled. The creation and handling of the
    //   histograms is taken care of by the HistogramManager class.
    //
    ////////////////////////////////////////////////////////////////////////

    #include "TRandom.h"
    #include "TDirectory.h"
    #include "TProcessID.h"

    #include "Event.h"

    ClassImp(EventHeader)
    ClassImp(Event)
    ClassImp(Track)
    ClassImp(HistogramManager)

    TClonesArray *Event::fgTracks = 0;
    TH1F *Event::fgHist = 0;

    //______________________________________________________________________________
    Event::Event()
    {
       // Create an Event object.
       // When the constructor is invoked for the first time, the class static
       // variable fgTracks is 0 and the TClonesArray fgTracks is created.

       if (!fgTracks) fgTracks = new TClonesArray("Track", 1000);
       fTracks = fgTracks;
       fHighPt = new TRefArray;
       fMuons  = new TRefArray;
       fNtrack = 0;
       fH      = 0;
       Int_t i0,i1;
       for (i0 = 0; i0 < 4; i0++) {
          for (i1 = 0; i1 < 4; i1++) {
             fMatrix[i0][i1] = 0.0;
          }
       }
       for (i0 = 0; i0 <10; i0++) fMeasures[i0] = 0;
       for (i0 = 0; i0 <20; i0++) fType[i0] = 0;
       fClosestDistance = 0;
       fEventName = 0;
       fWebHistogram.SetAction(this);
    }

    //______________________________________________________________________________
    Event::~Event()
    {
       Clear();
       if (fH == fgHist) fgHist = 0;
       delete fH; fH = 0;
       delete fHighPt; fHighPt = 0;
       delete fMuons;  fMuons = 0;
       delete [] fClosestDistance;
       if (fEventName) delete [] fEventName;
    }

    //______________________________________________________________________________
    void Event::Build(Int_t ev, Int_t arg5, Float_t ptmin) {
      char etype[20];
      Float_t sigmat, sigmas;
      gRandom->Rannor(sigmat,sigmas);
      Int_t ntrack   = Int_t(arg5 +arg5*sigmat/120.);
      Float_t random = gRandom->Rndm(1);

      //Save current Object count
      Int_t ObjectNumber = TProcessID::GetObjectCount();
      Clear();
      fHighPt->Delete();
      fMuons->Delete();

      Int_t nch = 15;
      if (ev >= 100)   nch += 3;
      if (ev >= 10000) nch += 3;
      if (fEventName) delete [] fEventName;
      fEventName = new char[nch];
      sprintf(fEventName,"Event%d_Run%d",ev,200);
      sprintf(etype,"type%d",ev%5);
      SetType(etype);
      SetHeader(ev, 200, 960312, random);
      SetNseg(Int_t(10*ntrack+20*sigmas));
      SetNvertex(Int_t(1+20*gRandom->Rndm()));
      SetFlag(UInt_t(random+0.5));
      SetTemperature(random+20.);

      for(UChar_t m = 0; m < 10; m++) {
         SetMeasure(m, Int_t(gRandom->Gaus(m,m+1)));
      }
      for(UChar_t i0 = 0; i0 < 4; i0++) {
        for(UChar_t i1 = 0; i1 < 4; i1++) {
           SetMatrix(i0,i1,gRandom->Gaus(i0*i1,1));
        }
      }

      fTriggerBits.SetBitNumber((UInt_t)(64*gRandom->Rndm(1)));
      fTriggerBits.SetBitNumber((UInt_t)(64*gRandom->Rndm(1)));
      fTriggerBits.SetBitNumber((UInt_t)(64*gRandom->Rndm(1)));

      //  Create and Fill the Track objects
      for (Int_t t = 0; t < ntrack; t++) AddTrack(random,ptmin);

      //Restore Object count
      //To save space in the table keeping track of all referenced objects
      //we assume that our events do not address each other. We reset the
      //object count to what it was at the beginning of the event.
      TProcessID::SetObjectCount(ObjectNumber);
    }

    //______________________________________________________________________________
    Track *Event::AddTrack(Float_t random, Float_t ptmin)
    {
       // Add a new track to the list of tracks for this event.
       // To avoid calling the very time consuming operator new for each track,
       // the standard but not well know C++ operator "new with placement"
       // is called. If tracks[i] is 0, a new Track object will be created
       // otherwise the previous Track[i] will be overwritten.

       TClonesArray &tracks = *fTracks;
       Track *track = new(tracks[fNtrack++]) Track(random);
       //Save reference to last Track in the collection of Tracks
       fLastTrack = track;
       //Save reference in fHighPt if track is a high Pt track
       if (track->GetPt() > ptmin)   fHighPt->Add(track);
       //Save reference in fMuons if track is a muon candidate
       if (track->GetMass2() < 0.11) fMuons->Add(track);
       return track;
    }

    //______________________________________________________________________________
    void Event::Clear(Option_t * /*option*/)
    {
       fTracks->Clear("C"); //will also call Track::Clear
       fHighPt->Delete();
       fMuons->Delete();
    }

    //______________________________________________________________________________
    void Event::Reset(Option_t * /*option*/)
    {
    // Static function to reset all static objects for this event
    //   fgTracks->Delete(option);

       delete fgTracks; fgTracks = 0;
       fgHist   = 0;
    }

    //______________________________________________________________________________
    void Event::SetHeader(Int_t i, Int_t run, Int_t date, Float_t random)
    {
       fNtrack = 0;
       fEvtHdr.Set(i, run, date);
       if (!fgHist) fgHist = new TH1F("hstat","Event Histogram",100,0,1);
       fH = fgHist;
       fH->Fill(random);
    }

    //______________________________________________________________________________
    void Event::SetMeasure(UChar_t which, Int_t what) {
       if (which<10) fMeasures[which] = what;
    }

    //______________________________________________________________________________
    void Event::SetRandomVertex() {
       // This delete is to test the relocation of variable length array
       if (fClosestDistance) delete [] fClosestDistance;
       if (!fNvertex) {
          fClosestDistance = 0;
          return;
       }
       fClosestDistance = new Double32_t[fNvertex];
       for (Int_t k = 0; k < fNvertex; k++ ) {
          fClosestDistance[k] = gRandom->Gaus(1,1);
       }
    }

    //______________________________________________________________________________
    Track::Track(const Track &orig) : TObject(orig)
    {
       // Copy a track object

       fPx = orig.fPx;
       fPy = orig.fPy;
       fPz = orig.fPx;
       fRandom = orig.fRandom;
       fMass2 = orig.fMass2;
       fBx = orig.fBx;
       fBy = orig.fBy;
       fMeanCharge = orig.fMeanCharge;
       fXfirst = orig.fXfirst;
       fXlast  = orig.fXlast;
       fYfirst = orig.fYfirst;
       fYlast  = orig.fYlast;
       fZfirst = orig.fZfirst;
       fZlast  = orig.fZlast;
       fCharge = orig.fCharge;

       fVertex[0] = orig.fVertex[0];
       fVertex[1] = orig.fVertex[1];
       fVertex[2] = orig.fVertex[2];
       fNpoint = orig.fNpoint;
       fNsp = orig.fNsp;
       if (fNsp) {
          fPointValue = new Double32_t[fNsp];
          for(int i=0; i<fNsp; i++) {
             fPointValue[i] = orig.fPointValue[i];
          }
       } else {
          fPointValue = 0;
       }
       fValid  = orig.fValid;

       fTriggerBits = orig.fTriggerBits;

    }

    //______________________________________________________________________________
    Track::Track(Float_t random) : TObject(),fTriggerBits(64)
    {
       // Create a track object.
       // Note that in this example, data members do not have any physical meaning.

       Float_t a,b,px,py;
       gRandom->Rannor(px,py);
       fPx = px;
       fPy = py;
       fPz = TMath::Sqrt(px*px+py*py);
       fRandom = 1000*random;
       if (fRandom < 10) fMass2 = 0.106;
       else if (fRandom < 100) fMass2 = 0.8;
       else if (fRandom < 500) fMass2 = 4.5;
       else if (fRandom < 900) fMass2 = 8.9;
       else  fMass2 = 9.8;
       gRandom->Rannor(a,b);
       fBx = 0.1*a;
       fBy = 0.1*b;
       fMeanCharge = 0.01*gRandom->Rndm(1);
       gRandom->Rannor(a,b);
       fXfirst = a*10;
       fXlast  = b*10;
       gRandom->Rannor(a,b);
       fYfirst = a*12;
       fYlast  = b*16;
       gRandom->Rannor(a,b);
       fZfirst = 50 + 5*a;
       fZlast  = 200 + 10*b;
       fCharge = Double32_t(Int_t(3*gRandom->Rndm(1)) - 1);

       fTriggerBits.SetBitNumber((UInt_t)(64*gRandom->Rndm(1)));
       fTriggerBits.SetBitNumber((UInt_t)(64*gRandom->Rndm(1)));
       fTriggerBits.SetBitNumber((UInt_t)(64*gRandom->Rndm(1)));

       fVertex[0] = gRandom->Gaus(0,0.1);
       fVertex[1] = gRandom->Gaus(0,0.2);
       fVertex[2] = gRandom->Gaus(0,10);
       fNpoint = Int_t(60+10*gRandom->Rndm(1));
       fNsp = Int_t(3*gRandom->Rndm(1));
       if (fNsp) {
          fPointValue = new Double32_t[fNsp];
          for(int i=0; i<fNsp; i++) {
             fPointValue[i] = i+1;
          }
       } else {
          fPointValue = 0;
       }
       fValid  = Int_t(0.6+gRandom->Rndm(1));
    }

    //______________________________________________________________________________
    void Track::Clear(Option_t * /*option*/)
    {
       fTriggerBits.Clear();
       delete [] fPointValue;
       fPointValue=0;
    }

    //______________________________________________________________________________
    HistogramManager::HistogramManager(TDirectory *dir)
    {
       // Create histogram manager object. Histograms will be created
       // in the "dir" directory.

       // Save current directory and cd to "dir".
       TDirectory *saved = gDirectory;
       dir->cd();

       fNtrack      = new TH1F("hNtrack",    "Ntrack",100,575,625);
       fNseg        = new TH1F("hNseg",      "Nseg",100,5800,6200);
       fTemperature = new TH1F("hTemperature","Temperature",100,19.5,20.5);
       fPx          = new TH1F("hPx",        "Px",100,-4,4);
       fPy          = new TH1F("hPy",        "Py",100,-4,4);
       fPz          = new TH1F("hPz",        "Pz",100,0,5);
       fRandom      = new TH1F("hRandom",    "Random",100,0,1000);
       fMass2       = new TH1F("hMass2",     "Mass2",100,0,12);
       fBx          = new TH1F("hBx",        "Bx",100,-0.5,0.5);
       fBy          = new TH1F("hBy",        "By",100,-0.5,0.5);
       fMeanCharge  = new TH1F("hMeanCharge","MeanCharge",100,0,0.01);
       fXfirst      = new TH1F("hXfirst",    "Xfirst",100,-40,40);
       fXlast       = new TH1F("hXlast",     "Xlast",100,-40,40);
       fYfirst      = new TH1F("hYfirst",    "Yfirst",100,-40,40);
       fYlast       = new TH1F("hYlast",     "Ylast",100,-40,40);
       fZfirst      = new TH1F("hZfirst",    "Zfirst",100,0,80);
       fZlast       = new TH1F("hZlast",     "Zlast",100,0,250);
       fCharge      = new TH1F("hCharge",    "Charge",100,-1.5,1.5);
       fNpoint      = new TH1F("hNpoint",    "Npoint",100,50,80);
       fValid       = new TH1F("hValid",     "Valid",100,0,1.2);

       // cd back to original directory
       saved->cd();
    }

    //______________________________________________________________________________
    HistogramManager::~HistogramManager()
    {
       // Clean up all histograms.

       // Nothing to do. Histograms will be deleted when the directory
       // in which tey are stored is closed.
    }

    //______________________________________________________________________________
    void HistogramManager::Hfill(Event *event)
    {
       // Fill histograms.

       fNtrack->Fill(event->GetNtrack());
       fNseg->Fill(event->GetNseg());
       fTemperature->Fill(event->GetTemperature());

       for (Int_t itrack = 0; itrack < event->GetNtrack(); itrack++) {
          Track *track = (Track*)event->GetTracks()->UncheckedAt(itrack);
          fPx->Fill(track->GetPx());
          fPy->Fill(track->GetPy());
          fPz->Fill(track->GetPz());
          fRandom->Fill(track->GetRandom());
          fMass2->Fill(track->GetMass2());
          fBx->Fill(track->GetBx());
          fBy->Fill(track->GetBy());
          fMeanCharge->Fill(track->GetMeanCharge());
          fXfirst->Fill(track->GetXfirst());
          fXlast->Fill(track->GetXlast());
          fYfirst->Fill(track->GetYfirst());
          fYlast->Fill(track->GetYlast());
          fZfirst->Fill(track->GetZfirst());
          fZlast->Fill(track->GetZlast());
          fCharge->Fill(track->GetCharge());
          fNpoint->Fill(track->GetNpoint());
          fValid->Fill(track->GetValid());
       }
    }