]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - MUON/MUONdigit.C
Minor syntax for the Alpha OSF
[u/mrichter/AliRoot.git] / MUON / MUONdigit.C
index ccef496a3e3a07e77af5795260fb8a57c4fb5c72..a1e3e8a14e190640a63db46051a799f87d1bc1c1 100644 (file)
@@ -1,69 +1,46 @@
-void MUONdigit (Int_t evNumber1=0,Int_t evNumber2=0,Int_t nCathode=1) 
+#include "iostream.h"
+
+void MUONdigit (Int_t evNumber1=0,Int_t evNumber2=0,Int_t nsignal  =25) 
 {
 /////////////////////////////////////////////////////////////////////////
-//                                                                     //
-//   This macros converts the pad-hit information into digits          //
-//                                                                     //
+//   This macro is a small example of a ROOT macro
+//   illustrating how to read the output of GALICE
+//   and do some analysis.
+//   
 /////////////////////////////////////////////////////////////////////////
 
 // Dynamically link some shared libs
 
    if (gClassTable->GetID("AliRun") < 0) {
-       gSystem->Load("$ALITOP/cern.so/lib/libpdfDUMMY.so");
-       gSystem->Load("$ALITOP/cern.so/lib/libPythia.so");
-       gSystem->Load("$ROOTSYS/lib/libEG.so");       
-       gSystem->Load("$ROOTSYS/lib/libEGPythia.so");    
-       gSystem->Load("libGeant3Dummy.so");        //a dummy version of Geant3
-       gSystem->Load("PHOS/libPHOSdummy.so");     //the standard Alice classes 
-       gSystem->Load("libgalice.so");             //the standard Alice classes 
+      gROOT->LoadMacro("loadlibs.C");
+      loadlibs();
    }
 
+
 // Connect the Root Galice file containing Geometry, Kine and Hits
 
    TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject("galice.root");
    if (file) file->Close(); 
    file = new TFile("galice.root","UPDATE");
    file->ls();
-//   file->Map();
 
-   printf ("\n File loaded !");
+   printf ("I'm after Map \n");
 
 // Get AliRun object from file or create it if not on file
 
    if (!gAlice) {
        gAlice = (AliRun*)file->Get("gAlice");
-       if (gAlice) printf("\n AliRun object found on file\n");
+       if (gAlice) printf("AliRun object found on file\n");
        if (!gAlice) gAlice = new AliRun("gAlice","Alice test program");
    }
-   printf ("\n gAlice created !\n");
-
-   AliMUON *MUON  = gAlice->GetDetector("MUON");
-
-   AliMUONchamber*  iChamber;
-   AliMUONsegmentation*  segmentation;
-
-   Int_t Npx[10];
-   Int_t Npy[10];
-
-   Int_t nxmax=1026;
-   Int_t nymax=1026;
-   AliMUONlist *elem[10][1026][1026];
-   TObjArray *obj=new TObjArray;
-
-   Int_t trk[500];
-   Int_t chtrk[500];
-
-
-   Int_t digits[3]; 
-   TVector *trinfo;
-   trinfo=new TVector(2);
-
+   printf ("I'm after gAlice \n");
+   
+   AliMUON *MUON  = (AliMUON*) gAlice->GetModule("MUON");
 //
-//   Loop over events 
+// Event Loop
 //
-
-   Int_t Nh=0;
-   Int_t Nh1=0;
+   Int_t nbgr_ev=0;
+   
    for (int nev=0; nev<= evNumber2; nev++) {
        Int_t nparticles = gAlice->GetEvent(nev);
        cout << "nev         " <<nev<<endl;
@@ -71,217 +48,27 @@ void MUONdigit (Int_t evNumber1=0,Int_t evNumber2=0,Int_t nCathode=1)
        if (nev < evNumber1) continue;
        if (nparticles <= 0) return;
 
-       TTree *TH = gAlice->TreeH();
-       Int_t ntracks = TH->GetEntries();
-       Int_t Nc=0;
-
-       Int_t counter=0;
-       //
-       // loop over cathodes
-       for (int icat=0; icat<nCathode; icat++) { 
-          //
-          // initialize the array of pointers 
-          printf("\n Initialize array of pointers \n");
-          for (int i=0; i<10; i++) {
-              for (int j=0; j<nymax; j++) {
-                  for (int k=0; k<nxmax; k++) {        
-                      elem[i][j][k]=0;
-                  }
-              }
-          }
-          printf("\n Finished !\n");
-     
-//
-//   Loop over tracks
-//
-          printf("\n ntracks: %d",ntracks);
-          
-          for (Int_t track=0; track<ntracks;track++) {
-              gAlice->ResetHits();
-              Int_t nbytes += TH->GetEvent(track);
-              if (MUON)  {
-                  for(AliMUONhit* mHit=(AliMUONhit*)MUON->FirstHit(-1); 
-                      mHit;
-                      mHit=(AliMUONhit*)MUON->NextHit()) 
-                  {
-                      Int_t   nch   = mHit->fChamber;  // chamber number
-                      Float_t x     = mHit->fX;        // x-pos of hit
-                      Float_t y     = mHit->fY;        // y-pos
-                      
-
-                      iChamber = &(MUON->Chamber(nch-1));
-                      response=iChamber->GetResponseModel();
-                      
-                      Int_t nsec=iChamber->Nsec();
-                      Int_t rmin = (Int_t)iChamber->frMin;
-                      Int_t rmax = (Int_t)iChamber->frMax;
-                      if (nch > 10) continue;
-                      
-//
-//
-                      for (AliMUONcluster* mPad=
-                               (AliMUONcluster*)MUON->FirstPad(mHit);
-                           mPad;
-                           mPad=(AliMUONcluster*)MUON->NextPad())
-                      {
-                          Int_t cathode  = mPad->fCathode;   // cathode number
-                          Int_t nhit     = mPad->fHitNumber; // hit number
-                          Int_t qtot     = mPad->fQ;         // charge
-                          Int_t ipx      = mPad->fPadX;      // pad  X
-                          Int_t ipy      = mPad->fPadY;      // pad  Y
-                          Int_t iqpad    = mPad->fQpad;      // charge per pad
-                          Int_t rsec     = mPad->fRSec;      // sector#
-//
-//
-                          if (cathode != (icat+1)) continue;
-                          segmentation=iChamber-> 
-                              GetSegmentationModel(cathode);
-                          Int_t Npx[nch-1]  = segmentation->Npx();
-                          Int_t Npy[nch-1]  = segmentation->Npy();
-                          
-                          
-                          Int_t npx=Npx[nch-1];
-                          Int_t npy=Npy[nch-1];
-                          Float_t thex, they;
-                          segmentation->GetPadCxy(ipx,ipy,thex,they);
-                          
-                          Float_t rpad=TMath::Sqrt(thex*thex+they*they);
-                          
-                          
-                       // check boundaries
-                          if (rpad < rmin || iqpad ==0 || rpad > rmax) continue;
-
-                          // fill the info array
+       Int_t nbgr_ev=Int_t(nev/nsignal);
+//       printf("nbgr_ev %d\n",nbgr_ev);
+       //if (MUON) MUON->Digitise(nev,nbgr_ev,"Add"," ","galice_bgr.root");
+       if (MUON) MUON->Digitise(nev,nbgr_ev,"rien"," ","galice_bgr.root");  
+//       char hname[30];
+//       sprintf(hname,"TreeD%d",nev);
+//       file->ls();
+   } // event loop 
+   file->Close();
+}
 
-                          trinfo(0)=(Float_t)track;
-                          trinfo(1)=(Float_t)iqpad;
 
-                          digits[0]=ipx;
-                          digits[1]=ipy;
-                          digits[2]=iqpad;
 
-                          // build the list of fired pads and update the info
-                          
-                          AliMUONlist *pdigit=elem[nch-1][ipy+npy][ipx+npx];
 
-                          if (pdigit==0) {
-                              obj->AddAtAndExpand(
-                                  new AliMUONlist(rpad,digits),counter);
-                              counter++;
-                              Int_t last=obj->GetLast();
-                              pdigit=(AliMUONlist*)obj->At(last);
-                              elem[nch-1][ipy+npy][ipx+npx]=pdigit;
-                              // list of tracks
-                              TObjArray *trlist=
-                                  (TObjArray*)pdigit->TrackList();
-                              trlist->Add(trinfo);
-                          } else {
-                              // update charge
-                              (*pdigit).fSignal+=iqpad;
-                              // update list of tracks
-                              TObjArray* trlist=(TObjArray*)pdigit
-                                  ->TrackList();
-                              Int_t nptracks=trlist->GetEntriesFast();
-                              
-                              for (Int_t tr=0;tr<nptracks;tr++) {
-                                  TVector *ptrk=(TVector*)trlist->At(tr);
-//                                TVector &vtrk = *ptrk;
-                                  Int_t entries=ptrk->GetNrows();
-                                  trk[tr]=Int_t(ptrk(0));
-                                  chtrk[tr]=Int_t(ptrk(1));
-                                  if (trk[tr]==track) {
-                                      chtrk[tr]+=iqpad;
-                                      trlist->RemoveAt(tr);
-                                      trinfo(0)=trk[tr];
-                                      trinfo(1)=chtrk[tr];
-                                      trlist->AddAt(trinfo,tr);
-                                  } else {
-                                      trlist->Add(trinfo);
-                                  }
-                              } // end loop over list of tracks for one pad
-                          } // end if pdigit
-                      } //end loop over clust
-                  } // hit loop
-              } // if MUON
-          } // track loop
 
-          Int_t tracks[10];
-          Int_t charges[10];
-          cout<<"Start filling digits now !"<<endl;
-          const Int_t zero_supm = 6;
-          const Int_t adc_satm = 1024;
-          Int_t nd=0;
-          
 
-          // start filling the digits
-          for (Int_t id=0; id<10; id++) {
-              Int_t nd=0;
-              for (Int_t iy=0;iy<nymax;iy++) {
-                  for (Int_t ix=0;ix<nxmax;ix++) {
-                      AliMUONlist *address=elem[id][iy][ix];
-                      if (address==0) continue; 
-                      // do zero-suppression and signal truncation
-                      Int_t q=address->fSignal;        
-                      // 
-                      if ( q <= zero_supm ) continue;
-                      if ( q > adc_satm)  q=adc_satm;
-                      digits[0]=address->fPadX;
-                      digits[1]=address->fPadY;
-                      digits[2]=address->fSignal;
 
-                      TObjArray* trlist=(TObjArray*)address->TrackList();
-                      Int_t nptracks=trlist->GetEntriesFast();
-                      // this was changed to accomodate the real number of tracks
-                      if (nptracks > 10) {
-                          cout<<"Attention - nptracks > 10 "<<nptracks<<endl;
-                          nptracks=10;
-                      }
-                      for (Int_t tr=0;tr<nptracks;tr++) {
-                          TVector *pp=(TVector*)trlist->At(tr);
-                          Int_t entries=pp->GetNrows();
-                          TVector &v2 = *pp;
-                          tracks[tr]=(Int_t)v2(0);
-                          charges[tr]=(Int_t)v2(1);
-                      } //end loop over list of tracks for one pad
 
-                      // fill digits
-                      MUON->AddDigits(id,tracks,charges,digits);
-                      nd++;
-                  } // end loop over ix
-              } // end loop over iy
-              cout<<"I'm out of the loops over pads"<<endl;
-              cout<<" id  nd  "<<id<<" "<<nd<<endl;
-          } //loop over chambers
-          cout<<"I'm out of the loops for digitisation"<<endl;
 
-          gAlice->TreeD()->Fill();
-          TTree *TD=gAlice->TreeD();
-          int ndig=TD->GetEntries();
-          cout<<"number of digits  "<<ndig<<endl;
 
-          TClonesArray *fDch;
-          for (int i=0;i<10;i++) {
-              fDch= MUON->DigitsAddress(i);
-              int ndig=fDch->GetEntriesFast();
-              printf (" i, ndig %d %d \n",i,ndig);
-          }
 
 
-          MUON->ResetDigits();
-          obj->Clear();
-       } //end loop over cathodes
 
-       char hname[30];
-       sprintf(hname,"TreeD%d",nev);
-       gAlice->TreeD()->Write(hname);
-       
-       file->ls();
-//       file->Map();
-   } // event loop 
-   obj->Delete();
-   delete obj;
 
-   file->Close();
-   cout<<"END  digitisation "<<endl;
-   
-}