5530fae2a54bccfc85978a0b741f71a1daf45ad9
[u/mrichter/AliRoot.git] / FMD / AliFMDRing.cxx
1 /**************************************************************************
2  * Copyright(c) 2004, 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 /* $Id$ */
17
18 //__________________________________________________________________
19 //
20 // Utility class to help implement collection of FMD modules into
21 // rings.  This is used by AliFMDSubDetector and AliFMD.  
22 //
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.
26 //
27 // Latest changes by Christian Holm Christensen
28 //
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
43
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";
50
51
52 //____________________________________________________________________
53 ClassImp(AliFMDRing);
54
55 //____________________________________________________________________
56 AliFMDRing::AliFMDRing(Char_t id, Bool_t detailed) 
57   : fId(id), 
58     fDetailed(detailed),
59     fActiveId(0),
60     fPrintboardBottomId(0),
61     fPrintboardTopId(0),
62     fRingId(0),
63     fSectionId(0),
64     fStripId(0),
65     fVirtualBackId(0),
66     fVirtualFrontId(0),
67     fBondingWidth(0),
68     fWaferRadius(0), 
69     fSiThickness(0),
70     fLowR(0), 
71     fHighR(0), 
72     fTheta(0), 
73     fNStrips(0), 
74     fRingDepth(0),
75     fLegRadius(0),
76     fLegLength(0),
77     fLegOffset(0),
78     fModuleSpacing(0),
79     fPrintboardThickness(0),
80     fShape(0),
81     fRotMatricies(0)
82 {
83   // Construct a alifmdring. 
84   // 
85   //     id             Id of the ring (either 'i' or 'o').
86   //     detailed       Whether the strips are made or not.
87   // 
88 }
89
90 //____________________________________________________________________
91 AliFMDRing::AliFMDRing(const AliFMDRing& other) 
92   : TObject(other),
93     fId(other.fId), 
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),
106     fLowR(other.fLowR), 
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)
119 {
120   // Copy constructor of a AliFMDRing. 
121 }
122
123 //____________________________________________________________________
124 AliFMDRing&
125 AliFMDRing::operator=(const AliFMDRing& other) 
126 {
127   // Assignment operator 
128   // 
129   fId                   = other.fId; 
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;
142   fLowR                 = other.fLowR; 
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;
153   if (other.fShape)  {
154     if (other.fShape->IsA() == TXTRU::Class()) 
155       ((TXTRU*)other.fShape)->Copy(*fShape);
156     else 
157       fShape = 0;
158   }
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);
164     TObject* o = 0;
165     while ((o = next())) fRotMatricies->Add(o);
166   }
167   return *this;
168 }
169
170   
171
172 //____________________________________________________________________
173 void 
174 AliFMDRing::Init() 
175 {
176   // Initialize the ring object.
177   // DebugGuard guard("AliFMDRing::Init");
178   AliDebug(10, "AliFMDRing::Init");
179   fPolygon.Clear();
180   SetupCoordinates();  
181 }
182
183 //____________________________________________________________________
184 AliFMDRing::~AliFMDRing() 
185 {
186   // Destructor - deletes shape and rotation matricies 
187   if (fShape) delete fShape;
188   if (fRotMatricies) delete fRotMatricies;
189 }
190
191
192 //____________________________________________________________________
193 void 
194 AliFMDRing::Browse(TBrowser* /* b */)
195 {
196   // DebugGuard guard("AliFMDRing::Browse");
197   AliDebug(10, "AliFMDRing::Browse");
198 }
199
200   
201 //____________________________________________________________________
202 void 
203 AliFMDRing::SetupCoordinates() 
204 {
205   // Calculates the parameters of the polygon shape. 
206   // 
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;
211
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 
223                                         + r2 * tanTheta2)) 
224                        / (1 + tanTheta2));
225   double yC        = tanTheta * xC;
226
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);
233 }
234
235 //____________________________________________________________________
236 bool
237 AliFMDRing::IsWithin(size_t moduleNo, double x, double y) const
238 {
239   // Checks if a point (x,y) is inside the module with number moduleNo 
240   //
241   // DebugGuard guard("AliFMDRing::IsWithin");
242   AliDebug(10, "AliFMDRing::IsWithin");
243   bool   ret            = false;
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.;
250     
251     // Rotate the point.
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);
254   
255     ret = fPolygon.Contains(xr,yr);
256   }
257   return ret;
258 }
259
260
261
262
263 //____________________________________________________________________
264 void
265 AliFMDRing::Draw(Option_t* option) const
266 {
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. 
270   // 
271   // DebugGuard guard("AliFMDRing::Draw");
272   AliDebug(10, "AliFMDRing::Draw");
273   // The unrotated coordinates of the polygon verticies
274   if (fPolygon.GetNVerticies() < 1) return;
275   
276   TVector2 v[6];
277   for (size_t i = 0; i < fPolygon.GetNVerticies(); i++) 
278     v[i] = fPolygon.GetVertex(i);
279   
280   Int_t    nModules  = 360 / Int_t(fTheta * 2);
281   Double_t dTheta    = fTheta * 2;
282   
283   TString opt(option);
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);
289     null->SetStats(0);
290     null->Draw(opt.Data());
291   }
292    
293   for (int i = 0; i < nModules; i++) {
294     Double_t theta = (i + .5) * dTheta;
295     AliFMDPolygon p;
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;
300         break;
301       }
302     }
303     p.Draw(opt.Data(), Form("MOD%c_%d", fId, i));
304   }
305   if (opt.Contains("0", TString::kIgnoreCase)) {
306     TArc* arcH = new TArc(0,0, fHighR);
307     arcH->SetLineStyle(2);
308     arcH->SetLineColor(4);
309     arcH->Draw();
310
311     TArc* arcL = new TArc(0,0, fLowR);
312     arcL->SetLineStyle(2);
313     arcL->SetLineColor(4);
314     arcL->Draw();
315   }
316 }
317
318 //____________________________________________________________________
319 void 
320 AliFMDRing::SetupGeometry(Int_t vacuumId, Int_t siId, Int_t pcbId, 
321                           Int_t pbRotId, Int_t idRotId)
322 {
323   // Setup the geometry of the ring.  It defines the volumes 
324   // RNGI or RNGO which can later be positioned in a sub-detector
325   // volume. 
326   // 
327   // The hieracy of the RNGx volume is 
328   // 
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)
344   //        
345   // Parameters: 
346   //
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 
352   //
353   // DebugGuard guard("AliFMDRing::SetupGeometry");
354   AliDebug(10, "AliFMDRing::SetupGeometry");
355
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
359   TString name;
360   TString name2;
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();
365   Double_t pars[10];
366   fRingDepth          = (fSiThickness 
367                          + fPrintboardThickness 
368                          + fLegLength 
369                          + fModuleSpacing);
370
371   // Ring virtual volume 
372   pars[0]             = rmin;
373   pars[1]             = rmax;
374   pars[2]             = fRingDepth / 2;
375   name                = Form(fgkRingFormat, fId);
376   fRingId             = gMC->Gsvolu(name.Data(), "TUBE", vacuumId, pars, 3);
377   
378   // Virtual volume for modules with long legs 
379   pars[1]             = rmax;
380   pars[3]             = -fTheta;
381   pars[4]             =  fTheta;
382   name                = Form(fgkVirtualFormat, 'F', fId);
383   fVirtualFrontId     = gMC->Gsvolu(name.Data(), "TUBS", vacuumId, pars, 5);
384
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);
389   
390   // Virtual mother volume for silicon
391   pars[2]             =  fSiThickness/2;
392   name2               =  name;
393   name                =  Form(fgkActiveFormat, fId);
394   fActiveId           =  gMC->Gsvolu(name.Data(), "TUBS", vacuumId , pars, 5);
395
396   if (fDetailed) {
397     // Virtual sector volumes 
398     name2               = name;
399     name                = Form(fgkSectorFormat, fId);
400     gMC->Gsdvn2(name.Data(), name2.Data(), 2, 2, -fTheta, vacuumId);
401     fSectionId          = gMC->VolId(name.Data());
402     
403     // Active strip volumes 
404     name2               = name;
405     name                = Form(fgkStripFormat, fId);
406     gMC->Gsdvt2(name.Data(), name2.Data(), dStrip, 1,stripOff, siId, fNStrips);
407     fStripId            = gMC->VolId(name.Data());
408   }
409   
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);
419
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);
426
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;
433     Int_t    idrot  = 0;
434     // Rotation matrix for virtual module volumes
435     gMC->Matrix(idrot, 90, theta, 90, fmod(90 + theta, 360), 0, 0);
436     fRotations[i] = idrot;
437   }
438
439
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) 
445                          * fBondingWidth)); 
446   
447   for (int i = 0; i < nModules; i++) {
448     TString  name2    = Form(fgkRingFormat, fId);
449
450     Int_t     id      = i;
451     // Double_t  theta   = (i + .5) * dTheta;
452     Bool_t    isFront = (i % 2 == 1);
453     Double_t  dz      = 0;
454     Double_t  w       = fRingDepth - (isFront ? 0 : fModuleSpacing);
455
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]);
460
461     // We only need to place the children once, they are copied when
462     // we place the other virtual volumes. 
463     if (i > 1) continue;
464     name2 = name;
465
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);
471
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");
480
481     // Support legs 
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,"");
491   }
492 }
493 //____________________________________________________________________
494 void 
495 AliFMDRing::Geometry(const char* mother, Int_t baseId, Double_t z, 
496                      Int_t /* pbRotId */, Int_t idRotId)
497 {
498   // Positions a RNGx volume inside a mother. 
499   // 
500   // Parameters
501   //  
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 
508   // 
509   // DebugGuard guard("AliFMDRing::Geometry");
510   AliDebug(10, "AliFMDRing::Geometry");
511   TString  name;
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, "");
517 }
518
519 //____________________________________________________________________
520 void 
521 AliFMDRing::SimpleGeometry(TList* nodes, 
522                            TNode* mother, 
523                            Int_t colour, 
524                            Double_t z, 
525                            Int_t n) 
526 {
527   // Make a simple geometry of the ring for event display. 
528   // 
529   // The simple geometry is made from ROOT TNode and TShape objects. 
530   // Note, that we cache the TShape and TRotMatrix objects used for
531   // this. 
532   // 
533   // Parameters
534   // 
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 
539   //    n         Detector number
540   //
541   // DebugGuard guard("AliFMDRing::SimpleGeometry");
542   AliDebug(10, "AliFMDRing::SimpleGeometry");
543   SetupCoordinates();
544
545   // If the shape hasn't been defined yet, we define it here. 
546   if (!fShape) {
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());
554     }
555     shape->DefineSection(0, - fSiThickness / 2, 1, 0, 0);
556     shape->DefineSection(1, + fSiThickness / 2, 1, 0, 0);
557     fShape = shape;
558     fShape->SetLineColor(colour);
559   }
560   
561   Int_t    nModules  = 360 / Int_t(fTheta * 2);
562   Double_t dTheta    = fTheta * 2;
563
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));
571       TRotMatrix* rot = 
572         new TRotMatrix(name.Data(), title.Data(), 
573                        90, theta, 90, fmod(90 + theta, 360), 0, 0);
574       fRotMatricies->AddAt(rot, i);
575     }
576   }
577
578   Double_t offsetZ   = (fSiThickness 
579                         + fPrintboardThickness
580                         + fLegLength + fModuleSpacing) / 2;
581
582   // Make all the nodes
583   for (int i = 0; i < nModules; i++) {
584     Bool_t    isFront = (i % 2 == 1);
585     mother->cd();
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, 
590                             0, 0, 
591                             z - offsetZ + (isFront ? fModuleSpacing : 0), 
592                             rot);
593     node->SetLineColor(colour);
594     nodes->Add(node);
595   }
596 }
597
598   
599
600 //____________________________________________________________________
601 void 
602 AliFMDRing::Gsatt() 
603 {
604   // Set drawing attributes for the RING 
605   // 
606   // DebugGuard guard("AliFMDRing::Gsatt");
607   AliDebug(10, "AliFMDRing::Gsatt");
608   TString name;
609   name = Form(fgkRingFormat,fId);
610   gMC->Gsatt(name.Data(), "SEEN", 0);
611
612   name = Form(fgkVirtualFormat, 'T', fId);
613   gMC->Gsatt(name.Data(), "SEEN", 0);
614
615   name = Form(fgkVirtualFormat, 'B', fId);
616   gMC->Gsatt(name.Data(), "SEEN", 0);
617
618   name = Form(fgkActiveFormat,fId);
619   gMC->Gsatt(name.Data(), "SEEN", 1);
620
621   name = Form(fgkPrintboardFormat, 'T', fId);
622   gMC->Gsatt(name.Data(), "SEEN", 1);
623
624   name = Form(fgkPrintboardFormat, 'B',fId);
625   gMC->Gsatt(name.Data(), "SEEN", 1);
626 }
627
628 //
629 // EOF
630 //