[u/mrichter/AliRoot.git] / MUON / MUONlist2digits.C

void MUONlist2digits (Int_t evNumber1=0,Int_t evNumber2=0,Int_t nCathode=1) 
{
/////////////////////////////////////////////////////////////////////////
//   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 
    }

// 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 ("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("AliRun object found on file\n");
      if (!gAlice) gAlice = new AliRun("gAlice","Alice test program");
   }
   printf ("I'm after gAlice \n");

       AliMUON *MUON  = gAlice->GetDetector("MUON");

   //   AliMUONchamber&  iChamber;
   AliMUONchamber*  iChamber;
   AliMUONsegmentation*  segmentation;

   const Float_t rmin[10] = {17.5,17.5,23.5,23.5,33.5,33.5,43,43,50,50};
   const Float_t rmax[10] = {91.5,91.5,122.5,122.5,173,173,221,221,256.5,256.5};


   Int_t Npx[10];
   Int_t Npy[10];


   for ( int i=0;i<10;i++) {
        Npx[i]=(Int_t)rmax[i]/0.75;
        Npy[i]=(Int_t)rmax[i]/0.5;
   }

       Int_t nxmax=1026;
       Int_t nymax=1026;
       AliMUONlist *elem[10][1026][1026];
       TObjArray *obj=new TObjArray;

       Int_t digits[3]; 
//
//   Loop over events 
//

   Int_t Nh=0;
   Int_t Nh1=0;
   for (int nev=0; nev<= evNumber2; nev++) {
       Int_t nparticles = gAlice->GetEvent(nev);
       cout << "nev         " << nev <<endl;
       cout << "nparticles  " << nparticles <<endl;
       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;

    for (int icat=0;icat<nCathode;icat++) { 

       //initialize the arrray of pointers !!!!!
       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;
           }
         }
       }

     
//
//   Loop over events
//
       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;
//
//
		   for (AliMUONcluster* mPad=(AliMUONcluster*)MUON->FirstPad(mHit);
		        mPad;
		        mPad=(AliMUONcluster*)MUON->NextPad())
		   {
	        	Int_t cathode     = mPad->fCathode;   // chamber number
		        Int_t nhit     = mPad->fHitNumber; // hit number
		        Int_t qtot     = mPad->fQ;         // charge
		        Int_t ipx      = mPad->fPadX;      // pad number on X
		        Int_t ipy      = mPad->fPadY;      // pad number on Y
		        Int_t iqpad    = mPad->fQpad;      // charge per pad
		        Int_t rpad     = mPad->fRpad;      // r-pos of pad
//
//
        if (cathode != (icat+1)) continue;
	segmentation=iChamber->GetSegmentationModel(cathode);
	//        Int_t Npx[nch-1]  = segmentation->Npx();
	//        Int_t Npy[nch-1]  = segmentation->Npy();

	if (cathode==1) {
             Int_t npx=Npx[nch-1];
             Int_t npy=Npy[nch-1];
        } else {
             Int_t npx=Npy[nch-1];
             Int_t npy=Npx[nch-1];
        }

	//         printf("nch, cathode, rmin, rmax, npx, npy %d %d %d %d %d %d\n",
	//	                  nch,cathode,rmin,rmax,npx,npy);

	//	       printf("icat, iqpad, ipx, ipy %d %d %d %d \n",
	//			         icat, iqpad, ipx, ipy);


	// check boundaries
        if (rpad < rmin[nch-1] || iqpad ==0 || rpad > rmax[nch-1]) continue;
	//         if (rpad < rmin[nch-1] || iqpad ==0 || rpad > 81.65) continue;
	 //         if (rpad < rmin || iqpad ==0 || rpad > rmax) continue;

	  // fill the info array
          TVector *trinfo;
           trinfo=new TVector(2);

          trinfo(0)=(Float_t)track;
          trinfo(1)=(Float_t)iqpad;

	  //          Int_t digits[3]; 
          digits[0]=ipx;
          digits[1]=ipy;
          digits[2]=iqpad;


          Int_t trk[50];
          Int_t chtrk[50];

    // 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 nentrobj=obj->GetEntriesFast();
       Int_t last=obj->GetLast();
       pdigit=(AliMUONlist*)obj->At(last);
       elem[nch-1][ipy+npy][ipx+npx]=pdigit;
       //       Int_t q=pdigit->fSignal;
       // 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);
           Int_t entries=ptrk->GetNrows();
           TVector &vtrk = *ptrk;
           trk[tr]=(Int_t)vtrk(0);
           chtrk[tr]=(Int_t)vtrk(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 "<<endl;
    //    Int_t digits[3];
    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;
	     //             Int_t npx=Npx[id];
	     //             Int_t npy=Npy[id];

	     //	     Int_t nx=2*npx;
	     //	     Int_t ny=2*npx;
   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 rpad=address->fRpad;
        Int_t q=address->fSignal;        
	if ( q <= zero_supm  || rpad > 55) continue;
	//        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 should be changed to accomodate the real number of tracks
        if (nptracks > 10) {
	    cout<<"Attention - nptracks > 3 "<<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+1,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;

    //    CathodeIndex[id+icat*10] = nd;     // update cathodes index
   

 } //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;

        int ndig1=(MUON->Dch1())->GetEntriesFast();
        cout<<"number of digits  1   "<<ndig1<<endl;
        int ndig2=(MUON->Dch2())->GetEntriesFast();
        cout<<"number of digits  2   "<<ndig2<<endl;
        int ndig3=(MUON->Dch3())->GetEntriesFast();
        cout<<"number of digits  3   "<<ndig3<<endl;
        int ndig4=(MUON->Dch4())->GetEntriesFast();
        cout<<"number of digits  4   "<<ndig4<<endl;
        int ndig5=(MUON->Dch5())->GetEntriesFast();
        cout<<"number of digits  5   "<<ndig5<<endl;
        int ndig6=(MUON->Dch6())->GetEntriesFast();
        cout<<"number of digits  6   "<<ndig6<<endl;
        int ndig7=(MUON->Dch7())->GetEntriesFast();
        cout<<"number of digits  7   "<<ndig7<<endl;
        int ndig8=(MUON->Dch8())->GetEntriesFast();
        cout<<"number of digits  8   "<<ndig8<<endl;
        int ndig9=(MUON->Dch9())->GetEntriesFast();
        cout<<"number of digits  9   "<<ndig9<<endl;
        int ndig10=(MUON->Dch10())->GetEntriesFast();
        cout<<"number of digits  10   "<<ndig10<<endl;

        MUON->ResetDigits();


	//        char hname[30];
	//        sprintf(hname,"TreeD%d",nev);
	//        gAlice->TreeD()->Write(hname);

    } //end loop over cathodes

        char hname[30];
        sprintf(hname,"TreeD%d",nev);
        gAlice->TreeD()->Write(hname);
    //    MUON->ResetDigits();
    //    if (CathodeIndex) CathodeIndex->Clear();

    //    file->Write();

     file->ls();
     file->Map();
   } // event loop 

   //     delete [] CathodeIndex;

     file->Close();
     cout<<"END  digitisation "<<endl;
   
}
Commit	Line	Data
fe4da5cc	1	void MUONlist2digits (Int_t evNumber1=0,Int_t evNumber2=0,Int_t nCathode=1)
	2	{
	3	/////////////////////////////////////////////////////////////////////////
	4	// This macro is a small example of a ROOT macro
	5	// illustrating how to read the output of GALICE
	6	// and do some analysis.
	7	//
	8	/////////////////////////////////////////////////////////////////////////
	9
	10	// Dynamically link some shared libs
	11
	12
	13	if (gClassTable->GetID("AliRun") < 0) {
	14	gSystem->Load("$ALITOP/cern.so/lib/libpdfDUMMY.so");
	15	gSystem->Load("$ALITOP/cern.so/lib/libPythia.so");
	16	gSystem->Load("$ROOTSYS/lib/libEG.so");
	17	gSystem->Load("$ROOTSYS/lib/libEGPythia.so");
	18	gSystem->Load("libGeant3Dummy.so"); //a dummy version of Geant3
	19	gSystem->Load("PHOS/libPHOSdummy.so"); //the standard Alice classes
	20	gSystem->Load("libgalice.so"); // the standard Alice classes
	21	}
	22
	23	// Connect the Root Galice file containing Geometry, Kine and Hits
	24
	25	TFile file = (TFile)gROOT->GetListOfFiles()->FindObject("galice.root");
	26	if (file) file->Close();
	27	file = new TFile("galice.root","UPDATE");
	28	file->ls();
	29	file->Map();
	30
	31	printf ("I'm after Map \n");
	32
	33	// Get AliRun object from file or create it if not on file
	34
	35	if (!gAlice) {
	36	gAlice = (AliRun*)file->Get("gAlice");
	37	if (gAlice) printf("AliRun object found on file\n");
	38	if (!gAlice) gAlice = new AliRun("gAlice","Alice test program");
	39	}
	40	printf ("I'm after gAlice \n");
	41
	42	AliMUON *MUON = gAlice->GetDetector("MUON");
	43
	44	// AliMUONchamber& iChamber;
	45	AliMUONchamber* iChamber;
	46	AliMUONsegmentation* segmentation;
	47
	48	const Float_t rmin[10] = {17.5,17.5,23.5,23.5,33.5,33.5,43,43,50,50};
	49	const Float_t rmax[10] = {91.5,91.5,122.5,122.5,173,173,221,221,256.5,256.5};
	50
	51
	52	Int_t Npx[10];
	53	Int_t Npy[10];
	54
	55
	56	for ( int i=0;i<10;i++) {
	57	Npx[i]=(Int_t)rmax[i]/0.75;
	58	Npy[i]=(Int_t)rmax[i]/0.5;
	59	}
	60
	61	Int_t nxmax=1026;
	62	Int_t nymax=1026;
	63	AliMUONlist *elem[10][1026][1026];
	64	TObjArray *obj=new TObjArray;
65
66	Int_t digits[3];
67	//
68	// Loop over events
69	//
70
71	Int_t Nh=0;
72	Int_t Nh1=0;
73	for (int nev=0; nev<= evNumber2; nev++) {
74	Int_t nparticles = gAlice->GetEvent(nev);
75	cout << "nev " << nev <<endl;
76	cout << "nparticles " << nparticles <<endl;
77	if (nev < evNumber1) continue;
78	if (nparticles <= 0) return;
79
80	TTree *TH = gAlice->TreeH();
81	Int_t ntracks = TH->GetEntries();
82	Int_t Nc=0;
83	//
84	Int_t counter=0;
85
86	for (int icat=0;icat<nCathode;icat++) {
87
88	//initialize the arrray of pointers !!!!!
89	for (int i=0; i<10; i++) {
90	for (int j=0; j<nymax; j++) {
91	for (int k=0; k<nxmax; k++) {
92	elem[i][j][k]=0;
93	}
94	}
95	}
96
97
98	//
99	// Loop over events
100	//
101	for (Int_t track=0; track<ntracks;track++) {
102	gAlice->ResetHits();
103	Int_t nbytes += TH->GetEvent(track);
104	if (MUON) {
105	for(AliMUONhit* mHit=(AliMUONhit*)MUON->FirstHit(-1);
106	mHit;
107	mHit=(AliMUONhit*)MUON->NextHit())
108	{
109	Int_t nch = mHit->fChamber; // chamber number
110	Float_t x = mHit->fX; // x-pos of hit
111	Float_t y = mHit->fY; // y-pos
112
113	iChamber = &(MUON->Chamber(nch-1));
114	response=iChamber->GetResponseModel();
115
116	Int_t nsec=iChamber->Nsec();
117	// Int_t rmin = (Int_t)iChamber->frMin;
118	// Int_t rmax = (Int_t)iChamber->frMax;
119	//
120	//
121	for (AliMUONcluster* mPad=(AliMUONcluster*)MUON->FirstPad(mHit);
122	mPad;
123	mPad=(AliMUONcluster*)MUON->NextPad())
124	{
125	Int_t cathode = mPad->fCathode; // chamber number
126	Int_t nhit = mPad->fHitNumber; // hit number
127	Int_t qtot = mPad->fQ; // charge
128	Int_t ipx = mPad->fPadX; // pad number on X
129	Int_t ipy = mPad->fPadY; // pad number on Y
130	Int_t iqpad = mPad->fQpad; // charge per pad
131	Int_t rpad = mPad->fRpad; // r-pos of pad
132	//
133	//
134	if (cathode != (icat+1)) continue;
135	segmentation=iChamber->GetSegmentationModel(cathode);
136	// Int_t Npx[nch-1] = segmentation->Npx();
137	// Int_t Npy[nch-1] = segmentation->Npy();
138
139	if (cathode==1) {
140	Int_t npx=Npx[nch-1];
141	Int_t npy=Npy[nch-1];
142	} else {
143	Int_t npx=Npy[nch-1];
144	Int_t npy=Npx[nch-1];
145	}
146
147	// printf("nch, cathode, rmin, rmax, npx, npy %d %d %d %d %d %d\n",
148	// nch,cathode,rmin,rmax,npx,npy);
149
150	// printf("icat, iqpad, ipx, ipy %d %d %d %d \n",
151	// icat, iqpad, ipx, ipy);
152
153
154	// check boundaries
155	if (rpad < rmin[nch-1] \|\| iqpad ==0 \|\| rpad > rmax[nch-1]) continue;
156	// if (rpad < rmin[nch-1] \|\| iqpad ==0 \|\| rpad > 81.65) continue;
157	// if (rpad < rmin \|\| iqpad ==0 \|\| rpad > rmax) continue;
158
159	// fill the info array
160	TVector *trinfo;
161	trinfo=new TVector(2);
162
163	trinfo(0)=(Float_t)track;
164	trinfo(1)=(Float_t)iqpad;
165
166	// Int_t digits[3];
167	digits[0]=ipx;
168	digits[1]=ipy;
169	digits[2]=iqpad;
170
171
172	Int_t trk[50];
173	Int_t chtrk[50];
174
175	// build the list of fired pads and update the info
176	AliMUONlist *pdigit=elem[nch-1][ipy+npy][ipx+npx];
177	if (pdigit==0) {
178	obj->AddAtAndExpand(new AliMUONlist(rpad,digits),counter);
179	counter++;
180	// Int_t nentrobj=obj->GetEntriesFast();
181	Int_t last=obj->GetLast();
182	pdigit=(AliMUONlist*)obj->At(last);
183	elem[nch-1][ipy+npy][ipx+npx]=pdigit;
184	// Int_t q=pdigit->fSignal;
185	// list of tracks
186	TObjArray trlist=(TObjArray)pdigit->TrackList();
187	trlist->Add(trinfo);
188	} else {
189	// update charge
190	(*pdigit).fSignal+=iqpad;
191	// update list of tracks
192	TObjArray* trlist=(TObjArray*)pdigit->TrackList();
193	Int_t nptracks=trlist->GetEntriesFast();
194	for (Int_t tr=0;tr<nptracks;tr++) {
195	TVector ptrk=(TVector)trlist->At(tr);
196	Int_t entries=ptrk->GetNrows();
197	TVector &vtrk = *ptrk;
198	trk[tr]=(Int_t)vtrk(0);
199	chtrk[tr]=(Int_t)vtrk(1);
200	if (trk[tr]==track) {
201	chtrk[tr]+=iqpad;
202	trlist->RemoveAt(tr);
203	trinfo(0)=trk[tr];
204	trinfo(1)=chtrk[tr];
205	trlist->AddAt(trinfo,tr);
206	} else {
207	trlist->Add(trinfo);
208	}
209	} //end loop over list of tracks for one pad
210	} // end if pdigit
211
212
213	} //end loop over clust
214
215	} // hit loop
216	} // if MUON
217	} // track loop
218
219	Int_t tracks[10];
220	Int_t charges[10];
221	cout<<"start filling digits "<<endl;
222	// Int_t digits[3];
223	const Int_t zero_supm = 6;
224	const Int_t adc_satm = 1024;
225	Int_t nd=0;
226
227
228	// start filling the digits
229	for (Int_t id=0;id<10;id++) {
230	Int_t nd=0;
231	// Int_t npx=Npx[id];
232	// Int_t npy=Npy[id];
233
234	// Int_t nx=2*npx;
235	// Int_t ny=2*npx;
236	for (Int_t iy=0;iy<nymax;iy++) {
237	for (Int_t ix=0;ix<nxmax;ix++) {
238	AliMUONlist *address=elem[id][iy][ix];
239	if (address==0) continue;
240	// do zero-suppression and signal truncation
241	Int_t rpad=address->fRpad;
242	Int_t q=address->fSignal;
243	if ( q <= zero_supm \|\| rpad > 55) continue;
244	// if ( q <= zero_supm ) continue;
245	if ( q > adc_satm) q=adc_satm;
246	digits[0]=(*address).fPadX;
247	digits[1]=(*address).fPadY;
248	digits[2]=address->fSignal;
249
250	TObjArray* trlist=(TObjArray*)address->TrackList();
251	Int_t nptracks=trlist->GetEntriesFast();
252	// this should be changed to accomodate the real number of tracks
253	if (nptracks > 10) {
254	cout<<"Attention - nptracks > 3 "<<nptracks<<endl;
255	nptracks=10;
256	}
257	for (Int_t tr=0;tr<nptracks;tr++) {
258	TVector pp=(TVector)trlist->At(tr);
259	Int_t entries=pp->GetNrows();
260	TVector &v2 = *pp;
261	tracks[tr]=(Int_t)v2(0);
262	charges[tr]=(Int_t)v2(1);
263	} //end loop over list of tracks for one pad
264
265	// fill digits
266	MUON->AddDigits(id+1,tracks,charges,digits);
267	nd++;
268
269	} // end loop over ix
270	} // end loop over iy
271	cout<<"I'm out of the loops over pads"<<endl;
272	cout<<" id nd "<<id<<" "<<nd<<endl;
273
274	// CathodeIndex[id+icat*10] = nd; // update cathodes index
275
276
277	} //loop over chambers
278	cout<<"I'm out of the loops for digitisation"<<endl;
279
280
281	gAlice->TreeD()->Fill();
282	TTree *TD=gAlice->TreeD();
283	int ndig=TD->GetEntries();
284	cout<<"number of digits "<<ndig<<endl;
285
286	int ndig1=(MUON->Dch1())->GetEntriesFast();
287	cout<<"number of digits 1 "<<ndig1<<endl;
288	int ndig2=(MUON->Dch2())->GetEntriesFast();
289	cout<<"number of digits 2 "<<ndig2<<endl;
290	int ndig3=(MUON->Dch3())->GetEntriesFast();
291	cout<<"number of digits 3 "<<ndig3<<endl;
292	int ndig4=(MUON->Dch4())->GetEntriesFast();
293	cout<<"number of digits 4 "<<ndig4<<endl;
294	int ndig5=(MUON->Dch5())->GetEntriesFast();
295	cout<<"number of digits 5 "<<ndig5<<endl;
296	int ndig6=(MUON->Dch6())->GetEntriesFast();
297	cout<<"number of digits 6 "<<ndig6<<endl;
298	int ndig7=(MUON->Dch7())->GetEntriesFast();
299	cout<<"number of digits 7 "<<ndig7<<endl;
300	int ndig8=(MUON->Dch8())->GetEntriesFast();
301	cout<<"number of digits 8 "<<ndig8<<endl;
302	int ndig9=(MUON->Dch9())->GetEntriesFast();
303	cout<<"number of digits 9 "<<ndig9<<endl;
304	int ndig10=(MUON->Dch10())->GetEntriesFast();
305	cout<<"number of digits 10 "<<ndig10<<endl;
306
307	MUON->ResetDigits();
308
309
310	// char hname[30];
311	// sprintf(hname,"TreeD%d",nev);
312	// gAlice->TreeD()->Write(hname);
313
314	} //end loop over cathodes
315
316	char hname[30];
317	sprintf(hname,"TreeD%d",nev);
318	gAlice->TreeD()->Write(hname);
319	// MUON->ResetDigits();
320	// if (CathodeIndex) CathodeIndex->Clear();
321
322	// file->Write();
323
324	file->ls();
325	file->Map();
326	} // event loop
327
328	// delete [] CathodeIndex;
329
330	file->Close();
331	cout<<"END digitisation "<<endl;
332
333	}