1 /**************************************************************************
2 * Copyright(c) 2004, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 //__________________________________________________________________
20 // Utility class to help implement collection of FMD modules into
21 // rings. This is used by AliFMDSubDetector and AliFMD.
23 // The AliFMD object owns the AliFMDRing objects, and the
24 // AliFMDSubDetector objects reference these. That is, the AliFMDRing
25 // objects are share amoung the AliFMDSubDetector objects.
27 // Latest changes by Christian Holm Christensen
29 #include "AliFMDRing.h" // ALIFMDRING_H
30 #include "AliLog.h" // ALILOG_H
31 #include "TMath.h" // ROOT_TMath
32 #include "TH2.h" // ROOT_TH2
33 #include "TVirtualMC.h" // ROOT_TVirtualMC
34 #include "TVector2.h" // ROOT_TVector2
35 #include "TBrowser.h" // ROOT_TBrowser
36 #include "TString.h" // ROOT_TString
37 #include "TArc.h" // ROOT_TArc
38 #include "TObjArray.h" // ROOT_TObjArray
39 #include "TXTRU.h" // ROOT_TXTRU
40 #include "TNode.h" // ROOT_TNode
41 #include "TRotMatrix.h" // ROOT_TRotMatrix
42 #include "TList.h" // ROOT_TList
44 const Char_t* AliFMDRing::fgkRingFormat = "FRG%c";
45 const Char_t* AliFMDRing::fgkVirtualFormat = "FV%c%c";
46 const Char_t* AliFMDRing::fgkActiveFormat = "FAC%c";
47 const Char_t* AliFMDRing::fgkSectorFormat = "FSE%c";
48 const Char_t* AliFMDRing::fgkStripFormat = "FST%c";
49 const Char_t* AliFMDRing::fgkPrintboardFormat = "FP%c%c";
52 //____________________________________________________________________
55 //____________________________________________________________________
56 AliFMDRing::AliFMDRing(Char_t id, Bool_t detailed)
60 fPrintboardBottomId(0),
79 fPrintboardThickness(0),
83 // Construct a alifmdring.
85 // id Id of the ring (either 'i' or 'o').
86 // detailed Whether the strips are made or not.
90 //____________________________________________________________________
91 AliFMDRing::AliFMDRing(const AliFMDRing& other)
94 fDetailed(other.fDetailed),
95 fActiveId(other.fActiveId),
96 fPrintboardBottomId(other.fPrintboardBottomId),
97 fPrintboardTopId(other.fPrintboardTopId),
98 fRingId(other.fRingId),
99 fSectionId(other.fSectionId),
100 fStripId(other.fStripId),
101 fVirtualBackId(other.fVirtualBackId),
102 fVirtualFrontId(other.fVirtualFrontId),
103 fBondingWidth(other.fBondingWidth),
104 fWaferRadius(other.fWaferRadius),
105 fSiThickness(other.fSiThickness),
107 fHighR(other.fHighR),
108 fTheta(other.fTheta),
109 fNStrips(other.fNStrips),
110 fRingDepth(other.fRingDepth),
111 fLegRadius(other.fLegRadius),
112 fLegLength(other.fLegLength),
113 fLegOffset(other.fLegOffset),
114 fModuleSpacing(other.fModuleSpacing),
115 fPrintboardThickness(other.fPrintboardThickness),
116 fRotations(other.fRotations),
117 fShape(other.fShape),
118 fRotMatricies(other.fRotMatricies)
120 // Copy constructor of a AliFMDRing.
123 //____________________________________________________________________
125 AliFMDRing::operator=(const AliFMDRing& other)
127 // Assignment operator
130 fDetailed = other.fDetailed;
131 fActiveId = other.fActiveId;
132 fPrintboardBottomId = other.fPrintboardBottomId;
133 fPrintboardTopId = other.fPrintboardTopId;
134 fRingId = other.fRingId;
135 fSectionId = other.fSectionId;
136 fStripId = other.fStripId;
137 fVirtualBackId = other.fVirtualBackId;
138 fVirtualFrontId = other.fVirtualFrontId;
139 fBondingWidth = other.fBondingWidth;
140 fWaferRadius = other.fWaferRadius;
141 fSiThickness = other.fSiThickness;
143 fHighR = other.fHighR;
144 fTheta = other.fTheta;
145 fNStrips = other.fNStrips;
146 fRingDepth = other.fRingDepth;
147 fLegRadius = other.fLegRadius;
148 fLegLength = other.fLegLength;
149 fLegOffset = other.fLegOffset;
150 fModuleSpacing = other.fModuleSpacing;
151 fPrintboardThickness = other.fPrintboardThickness;
152 fRotations = other.fRotations;
154 if (other.fShape->IsA() == TXTRU::Class())
155 ((TXTRU*)other.fShape)->Copy(*fShape);
159 if (other.fRotMatricies) {
160 Int_t n = other.fRotMatricies->GetEntries();
161 if (!fRotMatricies) fRotMatricies = new TObjArray(n);
162 else fRotMatricies->Expand(n);
163 TIter next(other.fRotMatricies);
165 while ((o = next())) fRotMatricies->Add(o);
172 //____________________________________________________________________
176 // Initialize the ring object.
177 // DebugGuard guard("AliFMDRing::Init");
178 AliDebug(10, "AliFMDRing::Init");
183 //____________________________________________________________________
184 AliFMDRing::~AliFMDRing()
186 // Destructor - deletes shape and rotation matricies
187 if (fShape) delete fShape;
188 if (fRotMatricies) delete fRotMatricies;
192 //____________________________________________________________________
194 AliFMDRing::Browse(TBrowser* /* b */)
196 // DebugGuard guard("AliFMDRing::Browse");
197 AliDebug(10, "AliFMDRing::Browse");
201 //____________________________________________________________________
203 AliFMDRing::SetupCoordinates()
205 // Calculates the parameters of the polygon shape.
207 // DebugGuard guard("AliFMDRing::SetupCoordinates");
208 AliDebug(10, "AliFMDRing::SetupCoordinates");
209 // Get out immediately if we have already done all this
210 if (fPolygon.GetNVerticies() > 1) return;
212 double tanTheta = TMath::Tan(fTheta * TMath::Pi() / 180.);
213 double tanTheta2 = TMath::Power(tanTheta,2);
214 double r2 = TMath::Power(fWaferRadius,2);
215 double yA = tanTheta * fLowR;
216 double lr2 = TMath::Power(fLowR, 2);
217 double hr2 = TMath::Power(fHighR,2);
218 double xD = fLowR + TMath::Sqrt(r2 - tanTheta2 * lr2);
219 double xD2 = TMath::Power(xD,2);
220 //double xD_2 = fLowR - TMath::Sqrt(r2 - tanTheta2 * lr2);
221 double yB = sqrt(r2 - hr2 + 2 * fHighR * xD - xD2);
222 double xC = ((xD + TMath::Sqrt(-tanTheta2 * xD2 + r2
225 double yC = tanTheta * xC;
227 fPolygon.AddVertex(fLowR, -yA);
228 fPolygon.AddVertex(xC, -yC);
229 fPolygon.AddVertex(fHighR, -yB);
230 fPolygon.AddVertex(fHighR, yB);
231 fPolygon.AddVertex(xC, yC);
232 fPolygon.AddVertex(fLowR, yA);
235 //____________________________________________________________________
237 AliFMDRing::IsWithin(size_t moduleNo, double x, double y) const
239 // Checks if a point (x,y) is inside the module with number moduleNo
241 // DebugGuard guard("AliFMDRing::IsWithin");
242 AliDebug(10, "AliFMDRing::IsWithin");
244 double r2 = x * x + y * y;
245 if (r2 < fHighR * fHighR && r2 > fLowR * fLowR) {
246 // double point_angle = TMath::ATan2(y, x);
247 // int n_modules = 360 / Int_t(fTheta * 2);
248 double m_angle = (.5 + moduleNo) * 2 * fTheta;
249 double m_radians = TMath::Pi() * m_angle / 180.;
252 double xr = x * TMath::Cos(-m_radians) - y * TMath::Sin(-m_radians);
253 double yr = x * TMath::Sin(-m_radians) + y * TMath::Cos(-m_radians);
255 ret = fPolygon.Contains(xr,yr);
263 //____________________________________________________________________
265 AliFMDRing::Draw(Option_t* option) const
267 // Draw a the shape of the ring into a 2D histogram. Useful for
268 // superimposing the actual shape of the ring onto a scatter plot of
269 // hits in the detector.
271 // DebugGuard guard("AliFMDRing::Draw");
272 AliDebug(10, "AliFMDRing::Draw");
273 // The unrotated coordinates of the polygon verticies
274 if (fPolygon.GetNVerticies() < 1) return;
277 for (size_t i = 0; i < fPolygon.GetNVerticies(); i++)
278 v[i] = fPolygon.GetVertex(i);
280 Int_t nModules = 360 / Int_t(fTheta * 2);
281 Double_t dTheta = fTheta * 2;
284 if (opt.Contains("B", TString::kIgnoreCase)) {
285 opt.Remove(opt.Index("B", 1, TString::kIgnoreCase),1);
286 TH1* null = new TH2F("null", "Null",
287 100, -fHighR * 1.1, fHighR * 1.1,
288 100, -fHighR * 1.1, fHighR * 1.1);
290 null->Draw(opt.Data());
293 for (int i = 0; i < nModules; i++) {
294 Double_t theta = (i + .5) * dTheta;
296 for (int j = 0; j < 6; j++) {
297 TVector2 vr(v[j].Rotate(TMath::Pi() * theta / 180.));
298 if (!p.AddVertex(vr.X(),vr.Y())) {
299 // std::cerr << "Draw of polygon " << i << " failed" << std::endl;
303 p.Draw(opt.Data(), Form("MOD%c_%d", fId, i));
305 if (opt.Contains("0", TString::kIgnoreCase)) {
306 TArc* arcH = new TArc(0,0, fHighR);
307 arcH->SetLineStyle(2);
308 arcH->SetLineColor(4);
311 TArc* arcL = new TArc(0,0, fLowR);
312 arcL->SetLineStyle(2);
313 arcL->SetLineColor(4);
318 //____________________________________________________________________
320 AliFMDRing::SetupGeometry(Int_t vacuumId, Int_t siId, Int_t pcbId,
321 Int_t pbRotId, Int_t idRotId)
323 // Setup the geometry of the ring. It defines the volumes
324 // RNGI or RNGO which can later be positioned in a sub-detector
327 // The hieracy of the RNGx volume is
329 // FRGx // Ring volume
330 // FVFx // Container of hybrid + legs
331 // FACx // Active volume (si sensor approx)
332 // FSEx // Section division
333 // FSTx // Strip division
334 // FPTx // Print board (bottom)
335 // FPBx // Print board (top)
336 // FLL // Support leg (long version)
337 // FVBx // Container of hybrid + legs
338 // FACx // Active volume (si sensor approx)
339 // FSEx // Section division
340 // FSTx // Strip division
341 // FPTx // Print board (bottom)
342 // FPBx // Print board (top)
343 // FSL // Support leg (long version)
347 // vacuumId Medium of inactive virtual volumes
348 // siId Medium of Silicon sensor (active)
349 // pcbId Medium of print boards
350 // pbRotId Print board rotation matrix
351 // idRotId Identity rotation matrix
353 // DebugGuard guard("AliFMDRing::SetupGeometry");
354 AliDebug(10, "AliFMDRing::SetupGeometry");
356 const TVector2& bCorner = fPolygon.GetVertex(3); // Third corner
357 const TVector2& aCorner = fPolygon.GetVertex(5); // First corner
358 const TVector2& cCorner = fPolygon.GetVertex(4); // Second corner
361 Double_t dStrip = (bCorner.Mod() - aCorner.Mod()) / fNStrips;
362 Double_t stripOff = aCorner.Mod();
363 Double_t rmin = fLowR;
364 Double_t rmax = bCorner.Mod();
366 fRingDepth = (fSiThickness
367 + fPrintboardThickness
371 // Ring virtual volume
374 pars[2] = fRingDepth / 2;
375 name = Form(fgkRingFormat, fId);
376 fRingId = gMC->Gsvolu(name.Data(), "TUBE", vacuumId, pars, 3);
378 // Virtual volume for modules with long legs
382 name = Form(fgkVirtualFormat, 'F', fId);
383 fVirtualFrontId = gMC->Gsvolu(name.Data(), "TUBS", vacuumId, pars, 5);
385 // Virtual volume for modules with long legs
386 pars[2] = (fRingDepth - fModuleSpacing) / 2;
387 name = Form(fgkVirtualFormat, 'B', fId);
388 fVirtualBackId = gMC->Gsvolu(name.Data(), "TUBS", vacuumId, pars, 5);
390 // Virtual mother volume for silicon
391 pars[2] = fSiThickness/2;
393 name = Form(fgkActiveFormat, fId);
394 fActiveId = gMC->Gsvolu(name.Data(), "TUBS", vacuumId , pars, 5);
397 // Virtual sector volumes
399 name = Form(fgkSectorFormat, fId);
400 gMC->Gsdvn2(name.Data(), name2.Data(), 2, 2, -fTheta, vacuumId);
401 fSectionId = gMC->VolId(name.Data());
403 // Active strip volumes
405 name = Form(fgkStripFormat, fId);
406 gMC->Gsdvt2(name.Data(), name2.Data(), dStrip, 1,stripOff, siId, fNStrips);
407 fStripId = gMC->VolId(name.Data());
410 // Print-board on back of module
411 pars[4] = TMath::Tan(TMath::Pi() * fTheta / 180) * fBondingWidth;
412 // Top of the print board
413 pars[0] = cCorner.Y() - pars[4];
414 pars[1] = bCorner.Y() - pars[4];
415 pars[2] = fPrintboardThickness / 2; // PCB half thickness
416 pars[3] = (bCorner.X() - cCorner.X()) / 2;
417 name = Form(fgkPrintboardFormat, 'T', fId);
418 fPrintboardTopId = gMC->Gsvolu(name.Data(), "TRD1", pcbId, pars, 4);
420 // Bottom of the print board
421 pars[0] = aCorner.Y() - pars[4];
422 pars[1] = cCorner.Y() - pars[4];
423 pars[3] = (cCorner.X() - aCorner.X()) / 2;
424 name = Form(fgkPrintboardFormat, 'B', fId);
425 fPrintboardBottomId = gMC->Gsvolu(name.Data(), "TRD1", pcbId, pars, 4);
427 // Define rotation matricies
428 Int_t nModules = 360 / Int_t(fTheta * 2);
429 Double_t dTheta = fTheta * 2;
430 fRotations.Set(nModules);
431 for (int i = 0; i < nModules; i++) {
432 Double_t theta = (i + .5) * dTheta;
434 // Rotation matrix for virtual module volumes
435 gMC->Matrix(idrot, 90, theta, 90, fmod(90 + theta, 360), 0, 0);
436 fRotations[i] = idrot;
440 // Int_t nModules = 360 / Int_t(fTheta * 2);
441 // Double_t dTheta = fTheta * 2;
442 Double_t pbTopL = (bCorner.X() - cCorner.X());
443 Double_t pbBotL = (cCorner.X() - aCorner.X());
444 Double_t yoffset = ((TMath::Tan(TMath::Pi() * fTheta / 180)
447 for (int i = 0; i < nModules; i++) {
448 TString name2 = Form(fgkRingFormat, fId);
451 // Double_t theta = (i + .5) * dTheta;
452 Bool_t isFront = (i % 2 == 1);
454 Double_t w = fRingDepth - (isFront ? 0 : fModuleSpacing);
456 // Place virtual module volume
457 name = Form(fgkVirtualFormat, (isFront ? 'F' : 'B'), fId);
458 dz = (w - fRingDepth) / 2;
459 gMC->Gspos(name.Data(), id, name2.Data(), 0., 0., dz,fRotations[i]);
461 // We only need to place the children once, they are copied when
462 // we place the other virtual volumes.
466 // Place active silicon wafer - this is put so that the front of
467 // the silicon is on the edge of the virtual volume.
468 name = Form(fgkActiveFormat, fId);
469 dz = (w - fSiThickness) / 2;
470 gMC->Gspos(name.Data(), id, name2.Data(),0.,0.,dz,idRotId);
472 // Place print board. This is put immediately behind the silicon
473 name = Form(fgkPrintboardFormat, 'T', fId);
474 dz = w / 2 - fSiThickness - fPrintboardThickness / 2;
475 gMC->Gspos(name.Data(), id, name2.Data(),
476 fLowR + pbBotL + pbTopL / 2, 0, dz, pbRotId, "ONLY");
477 name = Form(fgkPrintboardFormat, 'B', fId);
478 gMC->Gspos(name.Data(), id, name2.Data(),
479 fLowR + pbBotL / 2, 0, dz, pbRotId, "ONLY");
482 // This is put immediately behind the pringboard.
483 dz = (w / 2 - fSiThickness - fPrintboardThickness
484 - (fLegLength + (isFront ? fModuleSpacing : 0)) /2);
485 name = (isFront ? "FLL" : "FSL");
486 gMC->Gspos(name.Data(), id*10 + 1, name2.Data(),
487 aCorner.X() + fLegOffset + fLegRadius, 0., dz, idRotId, "");
488 Double_t y = cCorner.Y() - yoffset - fLegOffset - fLegRadius;
489 gMC->Gspos(name.Data(),id*10+2,name2.Data(),cCorner.X(), y,dz,idRotId,"");
490 gMC->Gspos(name.Data(),id*10+3,name2.Data(),cCorner.X(), -y,dz,idRotId,"");
493 //____________________________________________________________________
495 AliFMDRing::Geometry(const char* mother, Int_t baseId, Double_t z,
496 Int_t /* pbRotId */, Int_t idRotId)
498 // Positions a RNGx volume inside a mother.
502 // mother Mother volume to position the RNGx volume in
503 // baseId Base copy number
504 // z Z coordinate where the front of the active silicon
505 // should be in the mother volume, so we need to
506 // subtract half the ring width.
507 // idRotId Identity rotation matrix
509 // DebugGuard guard("AliFMDRing::Geometry");
510 AliDebug(10, "AliFMDRing::Geometry");
512 Double_t offsetZ = (fSiThickness
513 + fPrintboardThickness
514 + fLegLength + fModuleSpacing) / 2;
515 name = Form(fgkRingFormat, fId);
516 gMC->Gspos(name.Data(), baseId, mother, 0., 0., z - offsetZ, idRotId, "");
519 //____________________________________________________________________
521 AliFMDRing::SimpleGeometry(TList* nodes,
527 // Make a simple geometry of the ring for event display.
529 // The simple geometry is made from ROOT TNode and TShape objects.
530 // Note, that we cache the TShape and TRotMatrix objects used for
535 // nodes List of nodes to register all create nodes in
536 // mother Mother node to put the ring in.
537 // colour Colour of the nodes
538 // z Z position of the node in the mother volume
541 // DebugGuard guard("AliFMDRing::SimpleGeometry");
542 AliDebug(10, "AliFMDRing::SimpleGeometry");
545 // If the shape hasn't been defined yet, we define it here.
547 TString name(Form(fgkActiveFormat, fId));
548 TString title(Form("Shape of modules in %c Rings", fId));
549 Int_t n = fPolygon.GetNVerticies();
550 TXTRU* shape = new TXTRU(name.Data(), title.Data(), "void", n, 2);
551 for (Int_t i = 0; i < n; i++) {
552 const TVector2& v = fPolygon.GetVertex(i);
553 shape->DefineVertex(i, v.X(), v.Y());
555 shape->DefineSection(0, - fSiThickness / 2, 1, 0, 0);
556 shape->DefineSection(1, + fSiThickness / 2, 1, 0, 0);
558 fShape->SetLineColor(colour);
561 Int_t nModules = 360 / Int_t(fTheta * 2);
562 Double_t dTheta = fTheta * 2;
564 // If the roation matricies hasn't been defined yet, we do so here
565 if (!fRotMatricies) {
566 fRotMatricies = new TObjArray(nModules);
567 for (int i = 0; i < nModules; i++) {
568 Double_t theta = (i + .5) * dTheta;
569 TString name(Form("FMD_ring_%c_rot", fId));
570 TString title(Form("FMD Ring %c Rotation", fId));
572 new TRotMatrix(name.Data(), title.Data(),
573 90, theta, 90, fmod(90 + theta, 360), 0, 0);
574 fRotMatricies->AddAt(rot, i);
578 Double_t offsetZ = (fSiThickness
579 + fPrintboardThickness
580 + fLegLength + fModuleSpacing) / 2;
582 // Make all the nodes
583 for (int i = 0; i < nModules; i++) {
584 Bool_t isFront = (i % 2 == 1);
586 TRotMatrix* rot = static_cast<TRotMatrix*>(fRotMatricies->At(i));
587 TString name(Form("FAC%c_%d_%d", fId, n, i));
588 TString title(Form("Active FMD%d volume in %c Ring", n, fId));
589 TNode* node = new TNode(name.Data(), title.Data(), fShape,
591 z - offsetZ + (isFront ? fModuleSpacing : 0),
593 node->SetLineColor(colour);
600 //____________________________________________________________________
604 // Set drawing attributes for the RING
606 // DebugGuard guard("AliFMDRing::Gsatt");
607 AliDebug(10, "AliFMDRing::Gsatt");
609 name = Form(fgkRingFormat,fId);
610 gMC->Gsatt(name.Data(), "SEEN", 0);
612 name = Form(fgkVirtualFormat, 'T', fId);
613 gMC->Gsatt(name.Data(), "SEEN", 0);
615 name = Form(fgkVirtualFormat, 'B', fId);
616 gMC->Gsatt(name.Data(), "SEEN", 0);
618 name = Form(fgkActiveFormat,fId);
619 gMC->Gsatt(name.Data(), "SEEN", 1);
621 name = Form(fgkPrintboardFormat, 'T', fId);
622 gMC->Gsatt(name.Data(), "SEEN", 1);
624 name = Form(fgkPrintboardFormat, 'B',fId);
625 gMC->Gsatt(name.Data(), "SEEN", 1);