G4ParameterisationPara.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: G4ParameterisationPara.cc 69784 2013-05-15 09:16:06Z gcosmo $
00028 //
00029 // class G4ParameterisationPara Implementation file
00030 //
00031 // 26.05.03 - P.Arce, Initial version
00032 // 08.04.04 - I.Hrivnacova, Implemented reflection
00033 // 21.04.10 - M.Asai, Added gaps
00034 // --------------------------------------------------------------------
00035 
00036 #include "G4ParameterisationPara.hh"
00037 
00038 #include <iomanip>
00039 
00040 #include "G4PhysicalConstants.hh"
00041 #include "G4ThreeVector.hh"
00042 #include "G4Transform3D.hh"
00043 #include "G4RotationMatrix.hh"
00044 #include "G4VPhysicalVolume.hh"
00045 #include "G4ReflectedSolid.hh"
00046 #include "G4Para.hh"
00047 
00048 //--------------------------------------------------------------------------
00049 G4VParameterisationPara::
00050 G4VParameterisationPara( EAxis axis, G4int nDiv, G4double width,
00051                          G4double offset, G4VSolid* msolid,
00052                          DivisionType divType )
00053   :  G4VDivisionParameterisation( axis, nDiv, width, offset, divType, msolid )
00054 {
00055   G4Para* msol = (G4Para*)(msolid);
00056   if (msolid->GetEntityType() == "G4ReflectedSolid")
00057   {
00058     // Get constituent solid  
00059     G4VSolid* mConstituentSolid 
00060        = ((G4ReflectedSolid*)msolid)->GetConstituentMovedSolid();
00061     msol = (G4Para*)(mConstituentSolid);
00062     fmotherSolid = msol;
00063 
00064     // Create a new solid with inversed parameters
00065     G4Para* newSolid
00066       = new G4Para(msol->GetName(),
00067                    msol->GetXHalfLength(), 
00068                    msol->GetYHalfLength(),
00069                    msol->GetZHalfLength(),
00070                    std::atan(msol->GetTanAlpha()),
00071                    pi - msol->GetSymAxis().theta(),
00072                    msol->GetSymAxis().phi());
00073    
00074     msol = newSolid;
00075     fmotherSolid = newSolid;
00076     fReflectedSolid = true;
00077     fDeleteSolid = true;
00078   }    
00079 }
00080 
00081 //------------------------------------------------------------------------
00082 G4VParameterisationPara::~G4VParameterisationPara()
00083 {
00084 }
00085 
00086 //------------------------------------------------------------------------
00087 G4ParameterisationParaX::
00088 G4ParameterisationParaX( EAxis axis, G4int nDiv,
00089                          G4double width, G4double offset,
00090                          G4VSolid* msolid, DivisionType divType )
00091   :  G4VParameterisationPara( axis, nDiv, width, offset, msolid, divType )
00092 {
00093   CheckParametersValidity();
00094   SetType( "DivisionParaX" );
00095 
00096   G4Para* mpara = (G4Para*)(fmotherSolid);
00097   if( divType == DivWIDTH )
00098   {
00099     fnDiv = CalculateNDiv( 2*mpara->GetXHalfLength(), width, offset );
00100   }
00101   else if( divType == DivNDIV )
00102   {
00103     fwidth = CalculateWidth( 2*mpara->GetXHalfLength(), nDiv, offset );
00104   }
00105 
00106 #ifdef G4DIVDEBUG
00107   if( verbose >= 1 )
00108   {
00109     G4cout << " G4ParameterisationParaX - # divisions " << fnDiv
00110            << " = " << nDiv << G4endl
00111            << " Offset " << foffset << " = " << offset << G4endl
00112            << " Width " << fwidth << " = " << width << G4endl;
00113   }
00114 #endif
00115 }
00116 
00117 //------------------------------------------------------------------------
00118 G4double G4ParameterisationParaX::GetMaxParameter() const
00119 {
00120   G4Para* msol = (G4Para*)(fmotherSolid);
00121   return 2*msol->GetXHalfLength();
00122 }
00123 
00124 //------------------------------------------------------------------------
00125 G4ParameterisationParaX::~G4ParameterisationParaX()
00126 {
00127 }
00128 
00129 //------------------------------------------------------------------------
00130 void
00131 G4ParameterisationParaX::
00132 ComputeTransformation( const G4int copyNo, G4VPhysicalVolume *physVol ) const
00133 {
00134   G4Para* msol = (G4Para*)(fmotherSolid );
00135   G4double mdx = msol->GetXHalfLength( );
00136 
00137   //----- translation 
00138   G4ThreeVector origin(0.,0.,0.); 
00139   G4double posi = -mdx + foffset+(copyNo+0.5)*fwidth;
00140   origin.setX( posi ); 
00141   
00142 #ifdef G4DIVDEBUG
00143   if( verbose >= 2 )
00144   {
00145     G4cout << std::setprecision(8) << " G4ParameterisationParaX "
00146            << copyNo << G4endl
00147            << " Position: " << origin << " - Axis: " << faxis << G4endl;
00148   }
00149 #endif
00150 
00151   //----- set translation 
00152   physVol->SetTranslation( origin );
00153 }
00154 
00155 //--------------------------------------------------------------------------
00156 void
00157 G4ParameterisationParaX::
00158 ComputeDimensions(G4Para& para, const G4int,
00159                   const G4VPhysicalVolume*) const
00160 {
00161   //---- The division along X of a Para will result a Para
00162   G4Para* msol = (G4Para*)(fmotherSolid);
00163 
00164   //---- Get
00165   G4double pDx = fwidth/2. - fhgap;
00166   G4double pDy = msol->GetYHalfLength();
00167   G4double pDz = msol->GetZHalfLength();
00168   G4double pAlpha = std::atan(msol->GetTanAlpha());
00169   G4double pTheta = msol->GetSymAxis().theta();
00170   G4double pPhi = msol->GetSymAxis().phi();
00171  
00172   para.SetAllParameters ( pDx, pDy, pDz, pAlpha, pTheta, pPhi );
00173 
00174 #ifdef G4DIVDEBUG
00175   if( verbose >= 1 )
00176   {
00177     G4cout << " G4ParameterisationParaX::ComputeDimensions()"
00178            << " - Mother PARA " << G4endl;
00179     msol->DumpInfo();
00180     G4cout << " - Parameterised PARA: " << G4endl;
00181     para.DumpInfo();
00182   }
00183 #endif
00184 }
00185 
00186 //------------------------------------------------------------------------
00187 G4ParameterisationParaY::
00188 G4ParameterisationParaY( EAxis axis, G4int nDiv,
00189                          G4double width, G4double offset,
00190                          G4VSolid* msolid, DivisionType divType )
00191   :  G4VParameterisationPara( axis, nDiv, width, offset, msolid, divType )
00192 {
00193   CheckParametersValidity();
00194   SetType( "DivisionParaY" );
00195 
00196   G4Para* mpara = (G4Para*)(fmotherSolid);
00197   if( divType == DivWIDTH )
00198   {
00199     fnDiv = CalculateNDiv( 2*mpara->GetYHalfLength(), width, offset );
00200   }
00201   else if( divType == DivNDIV )
00202   {
00203     fwidth = CalculateWidth( 2*mpara->GetYHalfLength(), nDiv, offset );
00204   }
00205 
00206 #ifdef G4DIVDEBUG
00207   if( verbose >= 1 )
00208   {
00209     G4cout << " G4ParameterisationParaY - # divisions " << fnDiv
00210            << " = " << nDiv << G4endl
00211            << " Offset " << foffset << " = " << offset << G4endl
00212            << " Width " << fwidth << " = " << width << G4endl;
00213   }
00214 #endif
00215 }
00216 
00217 //------------------------------------------------------------------------
00218 G4ParameterisationParaY::~G4ParameterisationParaY()
00219 {
00220 }
00221 
00222 //------------------------------------------------------------------------
00223 G4double G4ParameterisationParaY::GetMaxParameter() const
00224 {
00225   G4Para* msol = (G4Para*)(fmotherSolid);
00226   return 2*msol->GetYHalfLength();
00227 }
00228 
00229 //------------------------------------------------------------------------
00230 void
00231 G4ParameterisationParaY::
00232 ComputeTransformation( const G4int copyNo, G4VPhysicalVolume *physVol ) const
00233 {
00234   G4Para* msol = (G4Para*)(fmotherSolid );
00235   G4double mdy = msol->GetYHalfLength( );
00236 
00237   //----- translation 
00238   G4ThreeVector origin(0.,0.,0.); 
00239   G4double posiY = -mdy + foffset+(copyNo+0.5)*fwidth;
00240   origin.setY( posiY );
00241   G4double posiX = posiY * msol->GetTanAlpha();
00242   origin.setX( posiX );
00243 
00244 #ifdef G4DIVDEBUG
00245   if( verbose >= 2 )
00246   {
00247     G4cout << std::setprecision(8) << " G4ParameterisationParaY "
00248            << copyNo << G4endl
00249            << " Position: " << origin << " - Axis: " << faxis << G4endl;
00250   }
00251 #endif
00252 
00253   //----- set translation 
00254   physVol->SetTranslation( origin );
00255 }
00256 
00257 //--------------------------------------------------------------------------
00258 void
00259 G4ParameterisationParaY::
00260 ComputeDimensions(G4Para& para, const G4int,
00261                   const G4VPhysicalVolume*) const
00262 {
00263   //---- The division along Y of a Para will result a Para
00264   G4Para* msol = (G4Para*)(fmotherSolid);
00265 
00266   //---- Get
00267   G4double pDx = msol->GetXHalfLength();
00268   G4double pDy = fwidth/2. - fhgap;
00269   G4double pDz = msol->GetZHalfLength();
00270   G4double pAlpha = std::atan(msol->GetTanAlpha());
00271   G4double pTheta = msol->GetSymAxis().theta();
00272   G4double pPhi = msol->GetSymAxis().phi();
00273  
00274   para.SetAllParameters ( pDx, pDy, pDz, pAlpha, pTheta, pPhi );
00275 
00276 #ifdef G4DIVDEBUG
00277   if( verbose >= -1 )
00278   {
00279     G4cout << " G4ParameterisationParaY::ComputeDimensions()"
00280            << " - Mother PARA " << G4endl;
00281     msol->DumpInfo();
00282     G4cout << " - Parameterised PARA: " << G4endl;
00283     para.DumpInfo();
00284   }
00285 #endif
00286 }
00287 
00288 //------------------------------------------------------------------------
00289 G4ParameterisationParaZ::
00290 G4ParameterisationParaZ( EAxis axis, G4int nDiv,
00291                          G4double width, G4double offset,
00292                          G4VSolid* msolid, DivisionType divType )
00293   :  G4VParameterisationPara( axis, nDiv, width, offset, msolid, divType )
00294 {
00295   CheckParametersValidity();
00296   SetType( "DivisionParaZ" );
00297 
00298   G4Para* mpara = (G4Para*)(fmotherSolid);
00299   if( divType == DivWIDTH )
00300   {
00301     fnDiv = CalculateNDiv( 2*mpara->GetZHalfLength(), width, offset );
00302   }
00303   else if( divType == DivNDIV )
00304   {
00305     fwidth = CalculateWidth( 2*mpara->GetZHalfLength(), nDiv, offset );
00306   }
00307 
00308 #ifdef G4DIVDEBUG
00309   if( verbose >= -1 )
00310   {
00311     G4cout << " G4ParameterisationParaZ - # divisions " << fnDiv
00312            << " = " << nDiv << G4endl
00313            << " Offset " << foffset << " = " << offset << G4endl
00314            << " Width " << fwidth << " = " << width << G4endl;
00315   }
00316 #endif
00317 }
00318 
00319 //------------------------------------------------------------------------
00320 G4ParameterisationParaZ::~G4ParameterisationParaZ()
00321 {
00322 }
00323 
00324 //------------------------------------------------------------------------
00325 G4double G4ParameterisationParaZ::GetMaxParameter() const
00326 {
00327   G4Para* msol = (G4Para*)(fmotherSolid);
00328   return 2*msol->GetZHalfLength();
00329 }
00330 
00331 //------------------------------------------------------------------------
00332 void
00333 G4ParameterisationParaZ::
00334 ComputeTransformation( const G4int copyNo, G4VPhysicalVolume *physVol ) const
00335 {
00336   G4Para* msol = (G4Para*)(fmotherSolid );
00337   G4double mdz = msol->GetZHalfLength( );
00338 
00339   //----- translation 
00340   G4double posi = -mdz + OffsetZ() + (copyNo+0.5)*fwidth;
00341   G4ThreeVector symAxis = msol->GetSymAxis();
00342   G4ThreeVector origin( symAxis * posi / symAxis.z() ); 
00343   
00344 #ifdef G4DIVDEBUG
00345   if( verbose >= 2 )
00346   {
00347     G4cout << std::setprecision(8) << " G4ParameterisationParaZ "
00348            << copyNo << G4endl
00349            << " Position: " << origin << " - Axis: " << faxis << G4endl;
00350   }
00351 #endif
00352 
00353   //----- set translation 
00354   physVol->SetTranslation( origin );
00355 }
00356 
00357 //--------------------------------------------------------------------------
00358 void
00359 G4ParameterisationParaZ::
00360 ComputeDimensions(G4Para& para, const G4int,
00361                   const G4VPhysicalVolume*) const
00362 {
00363   //---- The division along Z of a Para will result a Para
00364   G4Para* msol = (G4Para*)(fmotherSolid);
00365 
00366   //---- Get
00367   G4double pDx = msol->GetXHalfLength();
00368   G4double pDy = msol->GetYHalfLength();
00369   G4double pDz = fwidth/2. - fhgap;
00370   G4double pAlpha = std::atan(msol->GetTanAlpha());
00371   G4double pTheta = msol->GetSymAxis().theta();
00372   G4double pPhi = msol->GetSymAxis().phi();
00373  
00374   para.SetAllParameters ( pDx, pDy, pDz, pAlpha, pTheta, pPhi );
00375 
00376 #ifdef G4DIVDEBUG
00377   if( verbose >= -1 )
00378   {
00379     G4cout << " G4ParameterisationParaZ::ComputeDimensions()"
00380            << " - Mother PARA " << G4endl;
00381     msol->DumpInfo();
00382     G4cout << " - Parameterised PARA: " << G4endl;
00383     para.DumpInfo();
00384   }
00385 #endif
00386 }

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