FMD/scripts/NodeGeometry.C

   1 //____________________________________________________________________
   2 //
   3 //
   4 // $Id$
   5 //
   6 // Script I used for rapid prototyping of the FMD3 geometry - in
   7 // particular the support cone
   8 //
   9 /** @defgroup node_geom Simple geometry
  10     @ingroup FMD_script
  11 */
  12 #include <TGeometry.h>
  13 #include <TNode.h>
  14 #include <TXTRU.h>
  15 #include <TTUBE.h>
  16 #include <TTUBS.h>
  17 #include <TPCON.h>
  18 #include <TBRIK.h>
  19 #include <TCanvas.h>
  20 #include <vector>
  21 #include <algorithm>
  22 #include <cmath>
  23 #include <iostream>
  24
  25 //____________________________________________________________________
  26 /** @brief A 2D point
  27     @ingroup node_geom
  28  */
  29 struct point_t
  30 {
  31   point_t(double x=0, double y=0) : first(x), second(y) {}
  32   double first;
  33   double second;
  34 };
  35
  36 //____________________________________________________________________
  37 /** @brief Shape of a ring
  38     @ingroup node_geom
  39  */
  40 struct Ring
  41 {
  42   // typedef std::pair<double,double> point_t;
  43   typedef std::vector<point_t>     points_t;
  44   /** Constructor
  45       @param rL         Lower radius
  46       @param rH         Higer radius
  47       @param theta      Opening angle
  48       @param waferR     Wafer radius
  49       @param siThick    Silicon thickness
  50       @param staggering Staggering of modules */
  51   Ring(double rL, double rH, double theta,  double waferR,
  52        double siThick, double staggering)
  53     : fStaggering(staggering),
  54       fInnerRadius(rL),
  55       fOuterRadius(rH),
  56       fAngle(theta),
  57       fRadius(waferR),
  58       fThickness(siThick),
  59       fVerticies(6)
  60   {
  61     double tan_theta  = tan(fAngle * TMath::Pi() / 180.);
  62     double tan_theta2 = pow(tan_theta,2);
  63     double r2         = pow(fRadius,2);
  64     double ir2        = pow(fInnerRadius,2);
  65     double or2        = pow(fOuterRadius,2);
  66     double y_A        = tan_theta * fInnerRadius;
  67     double x_D        = fInnerRadius + sqrt(r2 - tan_theta2 * ir2);
  68     double x_D2       = pow(x_D,2);
  69     double y_B        = sqrt(r2 - or2 + 2 * fOuterRadius * x_D - x_D2);
  70     double x_C        = ((x_D + sqrt(-tan_theta2 * x_D2
  71                                      + r2 * (1 + tan_theta2)))
  72                          / (1 + tan_theta2));
  73     double y_C        = tan_theta * x_C;
  74
  75     fVerticies[0] = point_t(fInnerRadius,  y_A);
  76     fVerticies[1] = point_t(x_C,           y_C);
  77     fVerticies[2] = point_t(fOuterRadius,  y_B);
  78     fVerticies[3] = point_t(fOuterRadius, -y_B);
  79     fVerticies[4] = point_t(x_C,          -y_C);
  80     fVerticies[5] = point_t(fInnerRadius, -y_A);
  81   }
  82   /** Destructor */
  83   virtual ~Ring()
  84   {
  85     fVerticies.clear();
  86   }
  87   /** Create a shape
  88       @return pointer to new shape  */
  89   TShape* CreateShape()
  90   {
  91     std::cout << "Creating Module shape" << std::flush;
  92     TXTRU* moduleShape = new TXTRU("Module","Module", "", 6, 2);
  93     for (Int_t i = 0; i  < 6; i++) {
  94       std::cout << "." << std::flush;
  95       point_t& p = fVerticies[i];
  96       moduleShape->DefineVertex(i, p.first, p.second);
  97     }
  98     moduleShape->DefineSection(0, -fThickness/2, 1, 0, 0);
  99     moduleShape->DefineSection(1,  fThickness/2, 1, 0, 0);
 100     std::cout << std::endl;
 101     return (TShape*)moduleShape;
 102   }
 103   /** Create a node that represents a ring.
 104       @return Node */
 105   TNode* CreateRing(const char* name, double z)
 106   {
 107     std::cout << "Creating Ring node for " << name << std::flush;
 108     double bredth = fStaggering + fThickness;
 109     TShape* ringShape   = new TTUBE(Form("%sShape", name), "Ring Shape",
 110                                     "", fInnerRadius,
 111                                     fOuterRadius,bredth/2);
 112     TNode*  ringNode    = new TNode(Form("%sNode", name), "Ring Node",
 113                                     ringShape, 0, 0, z+bredth/2, 0);
 114     TShape* moduleShape = CreateShape();
 115     Int_t n = Int_t(360 / 2 / fAngle);
 116     for (Int_t i = 0; i < n; i++) {
 117       std::cout << "." << std::flush;
 118       ringNode->cd();
 119       Double_t theta  = 2  * fAngle * i;
 120       Double_t z      = -(bredth+fThickness)/2+(i%2?0:fStaggering);
 121       TRotMatrix* rot = new TRotMatrix(Form("%sRotation%02d", name, i),
 122                                        "Rotation", 90, theta, 90,
 123                                        fmod(90 + theta, 360), 0, 0);
 124       TNode* moduleNode = new TNode(Form("%sModule%02d", name, i),
 125                                     "Module", moduleShape, 0, 0, z,
 126                                     rot);
 127       moduleNode->SetFillColor(2);
 128       moduleNode->SetLineColor(2);
 129       moduleNode->SetLineWidth(2);
 130     }
 131     std::cout << std::endl;
 132     ringNode->SetVisibility(0);
 133     return ringNode;
 134   }
 135   double fStaggering;
 136   /** Inner radius */
 137   double fInnerRadius;
 138   /** Outer radius */
 139   double fOuterRadius;
 140   /** Opening angle (in degrees) */
 141   double fAngle;
 142   /** Radius (in centimeters) */
 143   double fRadius;
 144   /** Thickness */
 145   double fThickness;
 146   /** List of verticies */
 147   points_t fVerticies;
 148 };
 149
 150 //____________________________________________________________________
 151 /** @brief Shape of a detector
 152     @ingroup node_geom
 153  */
 154 struct Detector
 155 {
 156   /** Constructor
 157       @param id
 158       @param inner
 159       @param outer */
 160   Detector(Ring* inner, double iZ, Ring* outer=0, double oZ=0)
 161     : fInner(inner), fInnerZ(iZ), fOuter(outer), fOuterZ(oZ)
 162   {}
 163   /** Destructor */
 164   virtual ~Detector() {}
 165   /** Create rings */
 166   virtual void CreateRings()
 167   {
 168     if (fInner) fInner->CreateRing("inner", fInnerZ);
 169     if (fOuter) fOuter->CreateRing("outer", fOuterZ);
 170   }
 171   /** Create a node that represents the support */
 172   virtual void CreateSupport(double) { }
 173   /** Pointer to inner ring */
 174   Ring* fInner;
 175   /** Position in z of inner ring */
 176   double fInnerZ;
 177   /** Pointer to outer ring */
 178   Ring* fOuter;
 179   /** Position in z of inner ring */
 180   double fOuterZ;
 181 };
 182
 183 //____________________________________________________________________
 184 /** @brief FMD3 simple node geometry
 185     @ingroup node_geom
 186  */
 187 struct FMD3 : public Detector
 188 {
 189   /** Constructor
 190       @param inner Inner ring representation
 191       @param outer Outer ring representation */
 192   FMD3(Ring* inner, Ring* outer)
 193     : Detector(inner, -62.8,outer, -75.2)
 194   {
 195     fNoseRl     = 5.5;
 196     fNoseRh     = 6.7;
 197     fNoseDz     = 2.8 / 2;
 198     fNoseZ      = -46;
 199     fConeL      = 30.9;
 200     fBackRl     = 61 / 2;
 201     fBackRh     = 66.8 /2;
 202     fBackDz     = 1.4 / 2;
 203     fBeamDz     = .5 / 2;
 204     fBeamW      = 6;
 205     fFlangeR    = 49.25;
 206   }
 207   virtual ~FMD3() {}
 208   void CreateRings()
 209   {
 210     double zdist      = fConeL - 2 * fBackDz - 2 * fNoseDz;
 211     double tdist      = fBackRh - fNoseRh;
 212     double alpha      = tdist / zdist;
 213     double x, rl, rh, z;
 214     z  = fNoseZ - fConeL / 2;
 215     TPCON* fmd3Shape = new TPCON("fmd3Shape", "FMD 3 Shape", "", 0, 360, 7);
 216     x  = fNoseZ;
 217     rl = fNoseRl;
 218     rh = fNoseRh;
 219     fmd3Shape->DefineSection(0, x - z, rl, rh);
 220     x  = fNoseZ-2*fNoseDz;
 221     fmd3Shape->DefineSection(1, x - z, rl, rh);
 222     x  = fInnerZ - fInner->fStaggering - fInner->fThickness;
 223     rl = fInner->fInnerRadius;
 224     rh = fNoseRh + alpha * TMath::Abs(x-fNoseZ + 2 * fNoseDz);
 225     fmd3Shape->DefineSection(2, x - z, rl, rh);
 226     x  = fOuterZ;
 227     rl = fOuter->fInnerRadius;
 228     rh = fBackRh;
 229     fmd3Shape->DefineSection(3, x - z, rl, rh);
 230     x  = fNoseZ - zdist - 2 * fNoseDz;
 231     rl = fOuter->fInnerRadius;
 232     rh = fBackRh;
 233     fmd3Shape->DefineSection(4, x - z, rl, rh);
 234     x  = fNoseZ - zdist - 2 * fNoseDz;
 235     rl = fOuter->fInnerRadius;
 236     rh = fFlangeR;
 237     fmd3Shape->DefineSection(5, x - z, rl, rh);
 238     x  = fNoseZ - fConeL;
 239     rl = fOuter->fInnerRadius;
 240     rh = fFlangeR;
 241     fmd3Shape->DefineSection(6, x - z, rl, rh);
 242
 243     TNode* fmd3Node = new TNode("fmd3Node", "FMD3 Node", fmd3Shape,0,0,z,0);
 244     fmd3Node->SetLineColor(11);
 245     fmd3Node->SetFillColor(11);
 246     fmd3Node->SetVisibility(1);
 247     fmd3Node->cd();
 248     if (fInner) fInner->CreateRing("inner", fInnerZ-z);
 249     fmd3Node->cd();
 250     if (fOuter) fOuter->CreateRing("outer", fOuterZ-z);
 251     fmd3Node->cd();
 252     CreateSupport(fNoseZ - z);
 253   }
 254
 255   /** Create support volumes  */
 256   void CreateSupport(double noseZ)
 257   {
 258     TShape* noseShape = new TTUBE("noseShape", "Nose Shape", "",
 259                                   fNoseRl, fNoseRh, fNoseDz);
 260     TNode*  noseNode  = new TNode("noseNode", "noseNode", noseShape,
 261                                   0, 0, noseZ - fNoseDz, 0);
 262     noseNode->SetLineColor(0);
 263     double  zdist = fConeL - 2 * fBackDz - 2 * fNoseDz;
 264     double  tdist = fBackRh - fNoseRh;
 265     double  beamL = TMath::Sqrt(zdist * zdist + tdist * tdist);
 266     double  theta = -TMath::ATan2(tdist, zdist);
 267     TShape* backShape = new TTUBE("backShape", "Back Shape", "",
 268                                   fBackRl, fBackRh, fBackDz);
 269     TNode*  backNode  = new TNode("backNode", "backNode", backShape,
 270                                   0, 0, noseZ-2*fNoseDz-zdist-fBackDz, 0);
 271     backNode->SetLineColor(0);
 272     TShape* beamShape = new TBRIK("beamShape", "beamShape", "",
 273                                   fBeamDz, fBeamW / 2 , beamL / 2);
 274     Int_t    n = 8;
 275     Double_t r = fNoseRl + tdist / 2;
 276     for (Int_t i = 0; i < n; i++) {
 277       Double_t phi   = 360. / n * i;
 278       Double_t t     = 180. * theta / TMath::Pi();
 279       TRotMatrix* beamRotation = new TRotMatrix(Form("beamRotation%d", i),
 280                                                 Form("beamRotation%d", i),
 281                                                 180-t,phi,90,90+phi,t,phi);
 282       TNode* beamNode = new TNode(Form("beamNode%d", i),
 283                                   Form("beamNode%d", i), beamShape,
 284                                   r * TMath::Cos(phi / 180 * TMath::Pi()),
 285                                   r * TMath::Sin(phi / 180 * TMath::Pi()),
 286                                   noseZ-2*fNoseDz-zdist/2, beamRotation);
 287       beamNode->SetLineColor(0);
 288     }
 289     Double_t flangel    = (fFlangeR - fBackRh) / 2;
 290     TShape* flangeShape = new TBRIK("flangeShape", "FlangeShape", "",
 291                                     flangel, fBeamW / 2, fBackDz);
 292     n = 4;
 293     r = fBackRh + flangel;
 294     for (Int_t i = 0; i < n; i++) {
 295       Double_t phi = 360. / n * i + 180. / n;
 296       TRotMatrix* flangeRotation = new TRotMatrix(Form("flangeRotation%d", i),
 297                                                   Form("Flange Rotation %d",i),
 298                                                   90,phi,90,90+phi,0,0);
 299       TNode* flangeNode = new TNode(Form("flangeNode%d", i),
 300                                     Form("flangeNode%d", i),
 301                                     flangeShape,
 302                                     r * TMath::Cos(phi / 180 * TMath::Pi()),
 303                                     r * TMath::Sin(phi / 180 * TMath::Pi()),
 304                                     noseZ-2*fNoseDz-zdist-fBackDz,
 305                                     flangeRotation);
 306       flangeNode->SetLineColor(0);
 307
 308     }
 309   }
 310   /** Nose inner radius */
 311   double  fNoseRl;
 312   /** Nose outer radius */
 313   double  fNoseRh;
 314   /** Nose depth */
 315   double  fNoseDz;
 316   /** Nose start position */
 317   double  fNoseZ;
 318   /** Length of whole support structure */
 319   double  fConeL;
 320   /** Inner radius of back ring */
 321   double  fBackRl;
 322   /** Outer radius of back ring */
 323   double  fBackRh;
 324   /** Thickness of back ring */
 325   double  fBackDz;
 326   /** Thickness of beams */
 327   double  fBeamDz;
 328   /** Width of beams */
 329   double  fBeamW;
 330   /** Ending radius of flanges */
 331   double  fFlangeR;
 332 };
 333
 334 //____________________________________________________________________
 335 /** @brief Create a node geometry
 336     @ingroup node_geom
 337     @code
 338     .x NodeGeometry.C++
 339     @endcode
 340  */
 341 void
 342 NodeGeometry()
 343 {
 344   TGeometry* geometry = new TGeometry("geometry","geometry");
 345   TShape* topShape = new TBRIK("topShape", "topShape", "", 40, 40, 150);
 346   TNode* topNode = new TNode("topNode", "topNode", topShape, 0, 0, 0, 0);
 347   topNode->SetVisibility(0);
 348   topNode->cd();
 349
 350   Ring inner( 4.3, 17.2, 18, 13.4 / 2, .03, 1);
 351   Ring outer(15.6, 28.0,  9, 13.4 / 2, .03, 1);
 352   FMD3 fmd3(&inner, &outer);
 353   fmd3.CreateRings();
 354
 355   TCanvas* c = new TCanvas("c", "c", 800, 800);
 356   c->SetFillColor(1);
 357   geometry->Draw();
 358   // c->x3d("ogl");
 359 }
 360 //____________________________________________________________________
 361 //
 362 // EOF
 363 //