G4NURBStubesector.cc

Go to the documentation of this file.
00001 //
00002 // ********************************************************************
00003 // * License and Disclaimer                                           *
00004 // *                                                                  *
00005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
00006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
00007 // * conditions of the Geant4 Software License,  included in the file *
00008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
00009 // * include a list of copyright holders.                             *
00010 // *                                                                  *
00011 // * Neither the authors of this software system, nor their employing *
00012 // * institutes,nor the agencies providing financial support for this *
00013 // * work  make  any representation or  warranty, express or implied, *
00014 // * regarding  this  software system or assume any liability for its *
00015 // * use.  Please see the license in the file  LICENSE  and URL above *
00016 // * for the full disclaimer and the limitation of liability.         *
00017 // *                                                                  *
00018 // * This  code  implementation is the result of  the  scientific and *
00019 // * technical work of the GEANT4 collaboration.                      *
00020 // * By using,  copying,  modifying or  distributing the software (or *
00021 // * any work based  on the software)  you  agree  to acknowledge its *
00022 // * use  in  resulting  scientific  publications,  and indicate your *
00023 // * acceptance of all terms of the Geant4 Software license.          *
00024 // ********************************************************************
00025 //
00026 //
00027 // $Id$
00028 //
00029 // 
00030 // Olivier Crumeyrolle  12 September 1996
00031 
00032 // Tubesector builder implementation
00033 // OC 290896
00034 
00035 #include <sstream>
00036 
00037 #include "G4NURBStubesector.hh"
00038 #include "G4PhysicalConstants.hh"
00039 
00040 G4NURBStubesector::G4NURBStubesector(G4double r, G4double R,
00041                                      G4double DZ, G4double PHI1, G4double PHI2)
00042   : G4NURBS( 2, 3,  // linear along U, quadratic along V
00043              5, DecideNbrCtrlPts(PHI1, PHI2),  
00044                     // rectangle along U,  required stuff along V
00045                     // we must use a static function which
00046                     // take the two angles because the
00047                     // mother constructor is initialised
00048                     // before everything
00049              Regular,     // the knot vector along U will be generated
00050              RegularRep ) // circular like knot vector also
00051 {  
00052   // check angles
00053   G4double deltaPHI = PHI2-PHI1;
00054   while (deltaPHI <= 0) { PHI2 += twopi; deltaPHI += twopi; };
00055 
00056   G4int f = (int)floor(deltaPHI / (halfpi));  //number of pi/2 arcs
00057 
00058   const G4double mr = (r+R)/2;
00059 
00060   const G4double cp1 = std::cos(PHI1);
00061   const G4double sp1 = std::sin(PHI1);
00062   const G4double cp2 = std::cos(PHI2);
00063   const G4double sp2 = std::sin(PHI2);
00064 
00065 
00066   // define control points
00067   CP(mpCtrlPts[ 0] ,  cp1*mr, sp1*mr,  0, 1 );
00068   CP(mpCtrlPts[ 1] ,  cp1*mr, sp1*mr,  0, 1 );
00069   CP(mpCtrlPts[ 2] ,  cp1*mr, sp1*mr,  0, 1 );
00070   CP(mpCtrlPts[ 3] ,  cp1*mr, sp1*mr,  0, 1 );
00071   CP(mpCtrlPts[ 4] ,  cp1*mr, sp1*mr,  0, 1 );
00072 
00073   CP(mpCtrlPts[ 5] ,  cp1*mr, sp1*mr,  0, 1 );
00074   CP(mpCtrlPts[ 6] ,  cp1*mr, sp1*mr,  0, 1 );
00075   CP(mpCtrlPts[ 7] ,  cp1*mr, sp1*mr,  0, 1 );
00076   CP(mpCtrlPts[ 8] ,  cp1*mr, sp1*mr,  0, 1 );
00077   CP(mpCtrlPts[ 9] ,  cp1*mr, sp1*mr,  0, 1 );
00078 
00079   CP(mpCtrlPts[10] ,  cp1*r, sp1*r,  DZ, 1 );
00080   CP(mpCtrlPts[11] ,  cp1*R, sp1*R,  DZ, 1 );
00081   CP(mpCtrlPts[12] ,  cp1*R, sp1*R, -DZ, 1 );
00082   CP(mpCtrlPts[13] ,  cp1*r, sp1*r, -DZ, 1 );
00083   CP(mpCtrlPts[14] ,  cp1*r, sp1*r,  DZ, 1 );
00084 
00085   t_indCtrlPt  i = 15;
00086   G4double  srcAngle = PHI1;
00087   G4double  deltaAngleo2;
00088 
00089   G4double destAngle = halfpi + PHI1;
00090 
00091   for(; f > 0; f--)
00092   {
00093     // the first arc CP is already Done
00094 
00095     deltaAngleo2 = (destAngle - srcAngle) / 2;
00096     const G4double csa = std::cos(srcAngle);
00097     const G4double ssa = std::sin(srcAngle);
00098     const G4double tdao2 = std::tan(deltaAngleo2); 
00099 
00100     // to calculate the intermediate CP :
00101     // rotate by srcAngle the (1, tdao2) point
00102     const t_Coord x = csa - ssa*tdao2;
00103     const t_Coord y = ssa + csa*tdao2;
00104 
00105     // weight of the CP
00106     const G4Float weight = (std::cos(deltaAngleo2));
00107 
00108     // initialization. postfix ++ because i initialized to 15
00109     CP(mpCtrlPts[i++], x*r, y*r,  DZ, 1, weight);
00110     CP(mpCtrlPts[i++], x*R, y*R,  DZ, 1, weight);
00111     CP(mpCtrlPts[i++], x*R, y*R, -DZ, 1, weight);
00112     CP(mpCtrlPts[i++], x*r, y*r, -DZ, 1, weight);
00113     CP(mpCtrlPts[i++], x*r, y*r,  DZ, 1, weight);
00114 
00115     // end CP (which is the first CP of the next arc)
00116     const G4double cda = std::cos(destAngle);
00117     const G4double sda = std::sin(destAngle);
00118     CP(mpCtrlPts[i++], cda*r, sda*r,  DZ, 1);
00119     CP(mpCtrlPts[i++], cda*R, sda*R,  DZ, 1);
00120     CP(mpCtrlPts[i++], cda*R, sda*R, -DZ, 1);
00121     CP(mpCtrlPts[i++], cda*r, sda*r, -DZ, 1);
00122     CP(mpCtrlPts[i++], cda*r, sda*r,  DZ, 1);
00123 
00124     // prepare next arc
00125     srcAngle = destAngle;
00126     destAngle += halfpi;
00127   }
00128 
00129   // f == 0, final Arc
00130   // could be handled in the loops
00131 
00132   destAngle = PHI2;
00133   deltaAngleo2 = (destAngle - srcAngle) / 2;
00134   const G4double csa = std::cos(srcAngle);
00135   const G4double ssa = std::sin(srcAngle);
00136   const G4double tdao2 = std::tan(deltaAngleo2); 
00137 
00138   // to calculate the intermediate CP :
00139   // rotate by srcAngle the (1, tdao2) point
00140   const t_Coord x = csa - ssa*tdao2;
00141   const t_Coord y = ssa + csa*tdao2;
00142 
00143   // weight of the CP
00144   const G4Float weight = (std::cos(deltaAngleo2));
00145 
00146   // initialization.
00147   CP(mpCtrlPts[i++], x*r, y*r,  DZ, 1, weight);
00148   CP(mpCtrlPts[i++], x*R, y*R,  DZ, 1, weight);
00149   CP(mpCtrlPts[i++], x*R, y*R, -DZ, 1, weight);
00150   CP(mpCtrlPts[i++], x*r, y*r, -DZ, 1, weight);
00151   CP(mpCtrlPts[i++], x*r, y*r,  DZ, 1, weight);
00152 
00153   // end CP
00154   const G4double cda = std::cos(destAngle);
00155   const G4double sda = std::sin(destAngle);
00156   CP(mpCtrlPts[i++], cda*r, sda*r,  DZ, 1);
00157   CP(mpCtrlPts[i++], cda*R, sda*R,  DZ, 1);
00158   CP(mpCtrlPts[i++], cda*R, sda*R, -DZ, 1);
00159   CP(mpCtrlPts[i++], cda*r, sda*r, -DZ, 1);
00160   CP(mpCtrlPts[i++], cda*r, sda*r,  DZ, 1);
00161 
00162   if (i != (mtotnbrCtrlPts - 10) ) 
00163   { 
00164     G4cerr << "\nERROR: G4NURBStubesector::G4NURBStubesector: wrong index,"
00165            << i << " instead of " << (mtotnbrCtrlPts - 10)
00166            << "\n\tThe tubesector won't be correct."
00167            << G4endl;
00168   }
00169 
00170   CP(mpCtrlPts[i++] ,  cp2*mr, sp2*mr, 0, 1);
00171   CP(mpCtrlPts[i++] ,  cp2*mr, sp2*mr, 0, 1);
00172   CP(mpCtrlPts[i++] ,  cp2*mr, sp2*mr, 0, 1);
00173   CP(mpCtrlPts[i++] ,  cp2*mr, sp2*mr, 0, 1);
00174   CP(mpCtrlPts[i++] ,  cp2*mr, sp2*mr, 0, 1);
00175 
00176   CP(mpCtrlPts[i++] ,  cp2*mr, sp2*mr, 0, 1);
00177   CP(mpCtrlPts[i++] ,  cp2*mr, sp2*mr, 0, 1);
00178   CP(mpCtrlPts[i++] ,  cp2*mr, sp2*mr, 0, 1);
00179   CP(mpCtrlPts[i++] ,  cp2*mr, sp2*mr, 0, 1);
00180   CP(mpCtrlPts[i++] ,  cp2*mr, sp2*mr, 0, 1);
00181 
00182   // possible to put a DZ DZ -DZ -DZ DZ column to scratch
00183   // to a line instead of a point
00184 
00185   // creating the nurbs identity
00186   std::ostringstream  tmpstr;
00187   tmpstr << "Tubs" << " \tPHI1=" << PHI1 << " ; PHI2=" << PHI2;
00188   mpwhoami = new char [tmpstr.str().length() + 1];
00189   mpwhoami = std::strcpy(mpwhoami, tmpstr.str().c_str());
00190 }
00191 
00192 const char* G4NURBStubesector::Whoami() const
00193 {
00194   return mpwhoami;
00195 }
00196 
00197 G4NURBStubesector::~G4NURBStubesector()
00198 {
00199   if (mpwhoami) { delete [] mpwhoami; mpwhoami = 0; }
00200 }
00201 
00202 G4NURBStubesector::t_inddCtrlPt
00203 G4NURBStubesector::DecideNbrCtrlPts(G4double PHI1, G4double PHI2)
00204 {
00205   // check angles
00206   G4double deltaPHI = PHI2-PHI1;
00207   while (deltaPHI <= 0) { PHI2 += twopi; deltaPHI += twopi; }
00208   G4double k = deltaPHI / (halfpi);
00209 
00210   //    G4cerr << " k " << k << G4endl;
00211   //    G4cerr << " fk " << std::floor(k) << G4endl;
00212   //    G4cerr <<  " ifk " << ((int)(std::floor(k))) << G4endl;
00213   //    G4cerr << " n " << (2*((int)(std::floor(k))) + 7) << G4endl;
00214 
00215   return ( 2*((int)(std::floor(k))) + 7 );     
00216 }

Generated on Mon May 27 17:49:08 2013 for Geant4 by  doxygen 1.4.7