]> git.uio.no Git - u/mrichter/AliRoot.git/blame - MUON/AliMUONproto.cxx
Double inclusion of AliResponse removed.
[u/mrichter/AliRoot.git] / MUON / AliMUONproto.cxx
CommitLineData
a9e2aefa 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15/*
16$Log$
17Revision 1.1.2.1 2000/04/18 09:11:15 morsch
18Implementation of MUON Chamber Prototype Class
19Both read digits from raw data or use the Monte-Carlo.
20Rachid GUERNANE, IPN Lyon guernane@ipnl.in2p3.fr
21
22*/
23
24/*
25Implementation of MUON Chamber Prototype Class
26Both read digits from raw data or use the Monte-Carlo.
271-Feb-2000 Rachid GUERNANE, IPN Lyon guernane@ipnl.in2p3.fr
28*/
29
30
31////////////////////////////////////////////////
32// Manager and hits classes for set:PROTO //
33////////////////////////////////////////////////
34
35#include <TTUBE.h>
36#include <TBRIK.h>
37#include <TRotMatrix.h>
38#include <TNode.h>
39#include <TTree.h>
40#include <TRandom.h>
41#include <TObject.h>
42#include <TVector.h>
43#include <TObjArray.h>
44#include <TMinuit.h>
45#include <TParticle.h>
46#include <TROOT.h>
47#include <TFile.h>
48#include <TNtuple.h>
49#include <TCanvas.h>
50#include <TPad.h>
51#include <TDirectory.h>
52#include <TObjectTable.h>
53#include <AliPDG.h>
54#include <TArrayI.h>
55#include <AliDetector.h>
56
57#include "AliMUONproto.h"
58#include "TTUBE.h"
59#include "AliMUONClusterFinder.h"
60#include "AliRun.h"
61#include "AliMC.h"
62#include "iostream.h"
63#include "AliCallf77.h"
64#include "AliConst.h"
65#include "chainalice2.h"
66#include "AliMUONSegResV1.h"
67#include "AliMUONSegResV11.h"
68
69ClassImp(AliMUONproto)
70
71//___________________________________________
72AliMUONproto::AliMUONproto()
73 : AliMUON()
74{
75 cout << "\n Calling AliMUONproto constructor..." << endl;
76
77 //
78 //
79 //
80}
81
82//___________________________________________
83AliMUONproto::AliMUONproto(const char *name, const char *title)
84 : AliMUON(name,title)
85// AliMUON defines everything, but the chamber for NCH=1
86{
87
88// z-Positions of Chambers
89 const Float_t zch[1] = {975.};
90//
91// inner diameter
92 const Float_t dmi[1] = {0.};
93//
94// outer diameter
95 const Float_t dma[1] = {29.2};
96//
97//
98// Default Parameters for ALICE2 prototype
99
100//
101 (*fChambers)[0] = new AliMUONChamber();
102 AliMUONChamber* chamber = (AliMUONChamber*) (*fChambers)[0];
103 chamber->SetGid(0);
104 chamber->SetZPOS(zch[0]);
105//
106 chamber->InitGeo(zch[0]);
107 chamber->SetRInner(dmi[0]/2);
108 chamber->SetROuter(dma[0]/2);
109
110 for (int i = 0; i <= 99; i++) {
111 fThreshold[i] = 0.;
112 }
113
114}
115
116//___________________________________________
117void AliMUONproto::GetRawDigits(Int_t evnb, Int_t *lptr, Int_t ilen) {
118
119 Int_t ip = 0;
120 Int_t equip = 0;
121 Int_t nochnl;
122 Int_t loop;
123 Int_t val;
124 Int_t itype;
125 Int_t id;
126 Int_t serial;
127 Int_t equiplen;
128
129 Int_t digits[5];
130
131 AliMUON* MUON = (AliMUON*)gAlice->GetModule("MUON");
132 AliMUONSegmentationV1* seg = (AliMUONSegmentationV1*) Chamber(0).GetSegmentationModel(1);
133
134 Int_t tracks[10];
135 Int_t charges[10];
136
137 for (Int_t i = 0; i < 10; i++) {
138 tracks[i] = 0;
139 charges[i] = 0;
140 }
141
142 Int_t ich = 0;
143
144 nwtype:
145
146 itype = lptr[ip++];
147 id = lptr[ip++];
148 equiplen = lptr[ip++];
149
150 if (equiplen < 0 ) return;
151
152 if (itype != (int)(0XCACCA008)) {
153 ip += equiplen;
154 if (ip < ilen-2) goto nwtype;
155 }
156 else {
157 serial = id >> 16;
158 equip = id & 0x1;
159 if ((serial == 190) && (equip == 1)) {
160 for (loop = 0; loop < equiplen; loop++) {
161 nochnl = (lptr[ip] & 0x7ff000 ) >> 12;
162 val = lptr[ip] & 0x3ff;
163 // fill digits from raw data according to cathode connexions
164 if (group[nochnl][2][1]!=0)
165 digits[0] = group[nochnl][2][1] - 12;
166 else if (group[nochnl][1][1]!=0)
167 digits[0] = group[nochnl][1][1] - 12;
168 else
169 digits[0] = group[nochnl][0][1] - 12;
170 digits[1] = group[nochnl][0][0];
171 if (digits[0] != seg->Ix(digits[0], digits[1]))
172 printf("Pb pour ix=%d,iy=%d\n", digits[0], digits[1]);
173 digits[2] = val;
174 digits[3] = 0;
175 digits[4] = 0;
176 if (digits[2] >= fThreshold[nochnl])
177 MUON->AddDigits(ich, tracks, charges, digits);
178
179 ip++;
180 }
181 }
182 else
183 ip += equiplen;
184
185 if (ip < ilen-2) goto nwtype;
186 }
187
188 gAlice->TreeD()->Fill();
189 MUON->ResetDigits();
190
191 gAlice->TreeD()->Fill();
192 MUON->ResetDigits();
193
194 char hname[30];
195 sprintf(hname, "TreeD%d", evnb);
196 gAlice->TreeD()->Write(hname);
197 // reset tree
198 gAlice->TreeD()->Reset();
199
200}
201
202//___________________________________________
203void AliMUONproto::SetPadSize(Int_t id, Int_t isec, Float_t p1, Float_t p2)
204{
205 Int_t i=2*(id-1);
206 ((AliMUONChamber*)(*fChambers)[i])->SetPadSize(isec,p1,p2);
207}
208
209//___________________________________________
210void AliMUONproto::SetChargeSlope(Int_t id, Float_t p1)
211{
212 Int_t i=2*(id-1);
213 ((AliMUONChamber*) (*fChambers)[i])->SetChargeSlope(p1);
214}
215
216//___________________________________________
217void AliMUONproto::SetChargeSpread(Int_t id, Float_t p1, Float_t p2)
218{
219 Int_t i=2*(id-1);
220 ((AliMUONChamber*) (*fChambers)[i])->SetChargeSpread(p1,p2);
221}
222
223//___________________________________________
224void AliMUONproto::SetMaxAdc(Int_t id, Float_t p1)
225{
226 Int_t i=2*(id-1);
227 ((AliMUONChamber*) (*fChambers)[i])->SetMaxAdc(p1);
228}
229
230//___________________________________________
231void AliMUONproto::CreateGeometry()
232{
233 Int_t *idtmed = fIdtmed->GetArray()-1099;
234//
235 Float_t zpos;
236 Float_t bpar[3];
237 Int_t idGas=idtmed[1100];
238
239 AliMUONChamber *iChamber;
240//********************************************************************
241// Prototype ALICE2 **
242//********************************************************************
243 iChamber=(AliMUONChamber*) (*fChambers)[0];
244
245 zpos=iChamber->ZPosition();
246//
247 const Float_t X_POS = 11*.975/2; //half x
248 const Float_t Y_POS = 18*.55/2; //half y
249 const Float_t Z_POS = 0.325;
250
251 bpar[0] = X_POS;
252 bpar[1] = Y_POS;
253 bpar[2] = Z_POS;
254
255 gMC->Gsvolu("C01G", "BOX", idGas, bpar, 3);
256 gMC->Gspos("C01G", 1, "ALIC", -bpar[0], bpar[1], zpos, 0, "ONLY");
257}
258
259//___________________________________________
260void AliMUONproto::CreateMaterials()
261{
262 // *** DEFINITION OF AVAILABLE MUON MATERIALS ***
263 //
264 // Ar-CO2 gas
265 Float_t ag1[3] = { 39.95,12.01,16. };
266 Float_t zg1[3] = { 18.,6.,8. };
267 Float_t wg1[3] = { .8,.0667,.13333 };
268 Float_t dg1 = .001821;
269
270 Float_t epsil = .001; // Tracking precision,
271 Float_t tmaxfd = -20.; // Maximum angle due to field deflection
272 Float_t stmin = -.8;
273 Int_t ISXFLD = gAlice->Field()->Integ();
274 Float_t SXMGMX = gAlice->Field()->Max();
275
276 //
277 // --- Define the various materials for GEANT ---
278 AliMixture(22, "ArCO2 80%$", ag1, zg1, dg1, 3, wg1);
279
280 //
281 // Ar/CO2
282 AliMedium(1, "ARG_CO2 ", 22, 1, ISXFLD, SXMGMX, tmaxfd, fMaxStepGas,
283 fMaxDestepAlu, epsil, stmin);
284 // Air
285 //AliMedium(1, "AIR_CH_US ", 15, 1, ISXFLD, SXMGMX, tmaxfd, -1., -.3, epsil, stmin);
286}
287
288//___________________________________________
289void AliMUONproto::Init()
290{
291 printf("\n\n\n Start Init for Prototype ALICE2 - CPC chamber type\n\n\n");
292
293 //
294 // Initialize Tracking Chambers
295 //
296
297 ( (AliMUONChamber*) (*fChambers)[0])->Init();
298
299 //
300 // Set the chamber (sensitive region) GEANT identifier
301 AliMC* gMC = AliMC::GetMC();
302 ((AliMUONChamber*)(*fChambers)[0])->SetGid(gMC->VolId("C01G"));
303
304 printf("\n\n\n Finished Init for Prototype ALICE2 - CPC chamber type\n\n\n");
305
306}
307
308//___________________________________________
309void AliMUONproto::StepManager()
310{
311 Int_t copy, id;
312 static Int_t idvol;
313 static Int_t vol[2];
314 Int_t ipart;
315 TLorentzVector pos;
316 TLorentzVector mom;
317 Float_t theta,phi;
318 Float_t destep, step;
319
320 static Float_t eloss, eloss2, xhit, yhit, tlength;
321 const Float_t big=1.e10;
322
323 // modifs perso
324 static Float_t hits[14];
325
326 TClonesArray &lhits = *fHits;
327
328 //
329 // Set maximum step size for gas
330 // numed=gMC->GetMedium();
331 //
332 // Only charged tracks
333 if( !(gMC->TrackCharge()) ) return;
334 //
335 // Only gas gap inside chamber
336 // Tag chambers and record hits when track enters
337 idvol=-1;
338 id=gMC->CurrentVolID(copy);
339
340 for (Int_t i=1; i<=NCH; i++) {
341 if(id==((AliMUONChamber*)(*fChambers)[i-1])->GetGid()){
342 vol[0]=i;
343 idvol=i-1;
344 }
345 }
346 if (idvol == -1) return;
347 //
348 // Get current particle id (ipart), track position (pos) and momentum (mom)
349 gMC->TrackPosition(pos);
350 gMC->TrackMomentum(mom);
351
352 ipart = gMC->TrackPid();
353 //Int_t ipart1 = gMC->IdFromPDG(ipart);
354 //printf("ich, ipart %d %d \n",vol[0],ipart1);
355
356 //
357 // momentum loss and steplength in last step
358 destep = gMC->Edep();
359 step = gMC->TrackStep();
360
361 //
362 // record hits when track enters ...
363 if( gMC->IsTrackEntering()) {
364 gMC->SetMaxStep(fMaxStepGas);
365 Double_t tc = mom[0]*mom[0]+mom[1]*mom[1];
366 Double_t rt = TMath::Sqrt(tc);
367 Double_t pmom = TMath::Sqrt(tc+mom[2]*mom[2]);
368 Double_t tx=mom[0]/pmom;
369 Double_t ty=mom[1]/pmom;
370 Double_t tz=mom[2]/pmom;
371 Double_t s=((AliMUONChamber*)(*fChambers)[idvol])
372 ->GetResponseModel()
373 ->Pitch()/tz;
374 theta = Float_t(TMath::ATan2(rt,Double_t(mom[2])))*kRaddeg;
375 phi = Float_t(TMath::ATan2(Double_t(mom[1]),Double_t(mom[0])))*kRaddeg;
376 hits[0] = Float_t(ipart); // Geant3 particle type
377 hits[1] = pos[0]+s*tx; // X-position for hit
378 hits[2] = pos[1]+s*ty; // Y-position for hit
379 hits[3] = pos[2]+s*tz; // Z-position for hit
380 hits[4] = theta; // theta angle of incidence
381 hits[5] = phi; // phi angle of incidence
382 hits[8] = (Float_t) fNclusters; // first padhit
383 hits[9] = -1; // last pad hit
384
385 // modifs perso
386 hits[10] = mom[3]; // hit momentum P
387 hits[11] = mom[0]; // Px/P
388 hits[12] = mom[1]; // Py/P
389 hits[13] = mom[2]; // Pz/P
390 // fin modifs perso
391
392 // phi angle of incidence
393 tlength = 0;
394 eloss = 0;
395 eloss2 = 0;
396 xhit = pos[0];
397 yhit = pos[1];
398 // Only if not trigger chamber
399 if(idvol<10) {
400 //
401 // Initialize hit position (cursor) in the segmentation model
402 ((AliMUONChamber*) (*fChambers)[idvol])
403 ->SigGenInit(pos[0], pos[1], pos[2]);
404 } else {
405 //geant3->Gpcxyz();
406 //printf("In the Trigger Chamber #%d\n",idvol-9);
407 }
408 }
409 eloss2+=destep;
410
411 //
412 // Calculate the charge induced on a pad (disintegration) in case
413 //
414 // Mip left chamber ...
415 if( gMC->IsTrackExiting() || gMC->IsTrackStop() || gMC->IsTrackDisappeared()){
416 gMC->SetMaxStep(big);
417 eloss += destep;
418 tlength += step;
419
420 // Only if not trigger chamber
421 if(idvol<10) {
422 if (eloss > 0)
423 MakePadHits(0.5*(xhit+pos[0]),0.5*(yhit+pos[1]),eloss,idvol);
424 }
425
426 hits[6]=tlength;
427 hits[7]=eloss2;
428 if (fNclusters > (Int_t)hits[8]) {
429 hits[8]= hits[8]+1;
430 hits[9]= (Float_t) fNclusters;
431 }
432
433 new(lhits[fNhits++])
434 AliMUONHit(fIshunt,gAlice->CurrentTrack(),vol,hits);
435 eloss = 0;
436 //
437 // Check additional signal generation conditions
438 // defined by the segmentation
439 // model (boundary crossing conditions)
440 } else if
441 (((AliMUONChamber*) (*fChambers)[idvol])
442 ->SigGenCond(pos[0], pos[1], pos[2]))
443 {
444 ((AliMUONChamber*) (*fChambers)[idvol])
445 ->SigGenInit(pos[0], pos[1], pos[2]);
446// printf("\n-> MakePadHits, reason special %d",ipart);
447 if (eloss > 0)
448 MakePadHits(0.5*(xhit+pos[0]),0.5*(yhit+pos[1]),eloss,idvol);
449 xhit = pos[0];
450 yhit = pos[1];
451 eloss = destep;
452 tlength += step ;
453 //
454 // nothing special happened, add up energy loss
455 } else {
456 eloss += destep;
457 tlength += step ;
458 }
459}
460
461//___________________________________________
462void AliMUONproto::BuildGeometry()
463{
464 TNode *Node;
465 TNode *Top;
466
467 const int kColorMUON = kBlue;
468 //
469 Top=gAlice->GetGeometry()->GetNode("alice");
470 //
471 //
472 //
473 float dx, dy, dz, zpos;
474
475 const Float_t cz[1]={975.};
476
477 zpos=cz[0];
478
479 dx=11*.975/2;
480 dy=18*.55/2;
481 dz=0.325;
482
483 new TBRIK("C_MUON101","Mother volume for Proto.","void",dx*2,dy*2,dz*2);
484 TBRIK* PROTO = new TBRIK("PROT", "Proto. Sens. Region","void",dx,dy,dz);
485 Top->cd();
486 Node = new TNode("MUON101","ChamberNode","C_MUON101",0,0,zpos,"");
487 Node->SetLineColor(kColorMUON);
488 Node->SetVisibility(0);
489 fNodes->Add(Node);
490 Node->cd();
491 Node = new TNode("MUON201", "Proto. Sens. Region Node", PROTO, -dx, dy, dz);
492 Node->SetLineColor(kColorMUON);
493
494}
495
496//_____________________________________________________________________________
497void AliMUONproto::FindClusters(Int_t nev, Int_t last_entry)
498{
499
500//
501// Loop on chambers and on cathode planes
502//
503 for (Int_t icat = 0; icat < 2; icat++) {
504 gAlice->ResetDigits();
505 gAlice->TreeD()->GetEvent(last_entry+icat);
506
507 for (Int_t ich = 0; ich < NTrackingCh; ich++) {
508 AliMUONChamber* iChamber=(AliMUONChamber*) (*fChambers)[ich];
509 TClonesArray *MUONdigits = this->DigitsAddress(ich);
510 if (MUONdigits == 0) continue;
511 //
512 // Get ready the current chamber stuff
513 //
514
515 AliMUONResponse* response = iChamber->GetResponseModel();
516 AliMUONSegmentation* seg = iChamber->GetSegmentationModel(icat+1);
517 AliMUONClusterFinder* rec = iChamber->GetReconstructionModel();
518
519 if (seg) {
520 rec->SetSegmentation(seg);
521 rec->SetResponse(response);
522 rec->SetDigits(MUONdigits);
523 rec->SetChamber(ich);
524 rec->FindRawClusters();
525 }
526
527 TClonesArray *fRch;
528 fRch=RawClustAddress(ich);
529 fRch->Sort();
530 // it seems to work
531
532 } // for ich
533 // fill the tree
534 gAlice->TreeR()->Fill();
535
536 ResetRawClusters();
537 } // for icat
538
539 char hname[30];
540 sprintf(hname,"TreeR%d",nev);
541 gAlice->TreeR()->Write(hname);
542 gAlice->TreeR()->Reset();
543}
544
545//_____________________________________________________________________________
546void AliMUONproto::SetThreshold()
547{
548
549 ifstream inputFile("/home/alice/guernane/aliroot/pro/MUON/crped190.dat", ios::in);
550
551 if (inputFile.fail()) {
552 cout << "Error opening file" << endl;
553 exit(2);
554 }
555
556 char buff[32];
557 Int_t Serial;
558 Int_t Ntrigger;
559 Int_t Nchannel;
560 Int_t i1;
561 Int_t i2;
562
563 inputFile >> buff;
564
565 inputFile >> Serial;
566 inputFile >> Ntrigger;
567 inputFile >> Nchannel;
568 inputFile >> i1;
569 inputFile >> i2;
570
571 Float_t ped0[Nchannel];
572 Float_t sig0[Nchannel];
573 Float_t ped1[Nchannel];
574 Float_t sig1[Nchannel];
575 Int_t ichannel;
576
577 for (Int_t i = 0; i < Nchannel-1; i++) {
578 ped0[i] = 0;
579 sig0[i] = 0;
580 ped1[i] = 0;
581 sig1[i] = 0;
582 }
583
584 for (Int_t i = 0; i < Nchannel-1; i++) {
585 inputFile >> ichannel;
586 inputFile >> ped0[i];
587 inputFile >> sig0[i];
588 inputFile >> ped1[i];
589 inputFile >> sig1[i];
590 fThreshold[i] = fNsigma*sig1[i];
591 }
592
593 inputFile.close();
594}