Back to home page

EIC code displayed by LXR
 
 

    

File indexing completed on 2025-01-31 09:21:53 

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025 //
0026 /*
0027  * ============================================================================
0028  *
0029  *       Filename:  CexmcReconstructor.cc
0030  *
0031  *    Description:  reconstructor base class
0032  *
0033  *        Version:  1.0
0034  *        Created:  02.12.2009 16:21:15
0035  *       Revision:  none
0036  *       Compiler:  gcc
0037  *
0038  *         Author:  Alexey Radkov (), 
0039  *        Company:  PNPI
0040  *
0041  * ============================================================================
0042  */
0043 
0044 #include "CexmcReconstructor.hh"
0045 #include "CexmcReconstructorMessenger.hh"
0046 #include "CexmcEnergyDepositStore.hh"
0047 #include "CexmcRunManager.hh"
0048 
0049 
0050 CexmcReconstructor::CexmcReconstructor() : hasBasicTrigger( false ),
0051     epDefinitionAlgorithm( CexmcEntryPointBySqrtEDWeights ),
0052     epDepthDefinitionAlgorithm( CexmcEntryPointDepthPlain ),
0053     csAlgorithm( CexmcSelectAllCrystals ), useInnerRefCrystal( false ),
0054     epDepth( 0 ), theAngle( 0 ), calorimeterEDLeftAdjacent( 0 ),
0055     calorimeterEDRightAdjacent( 0 ), collectEDInAdjacentCrystals( false ),
0056     targetEPInitialized( false ), messenger( NULL )
0057 {
0058     G4RunManager *      runManager( G4RunManager::GetRunManager() );
0059     const CexmcSetup *  setup( static_cast< const CexmcSetup * >(
0060                                 runManager->GetUserDetectorConstruction() ) );
0061     calorimeterGeometry = setup->GetCalorimeterGeometry();
0062     targetTransform = setup->GetTargetTransform();
0063     calorimeterLeftTransform = setup->GetCalorimeterLeftTransform();
0064     calorimeterRightTransform = setup->GetCalorimeterRightTransform();
0065 
0066     messenger = new CexmcReconstructorMessenger( this );
0067 }
0068 
0069 
0070 CexmcReconstructor::~CexmcReconstructor()
0071 {
0072     delete messenger;
0073 }
0074 
0075 
0076 void  CexmcReconstructor::Reconstruct(
0077                                     const CexmcEnergyDepositStore *  edStore )
0078 {
0079     ReconstructEntryPoints( edStore );
0080     if ( hasBasicTrigger )
0081         ReconstructTargetPoint();
0082     if ( hasBasicTrigger )
0083         ReconstructAngle();
0084 }
0085 
0086 
0087 G4bool  CexmcReconstructor::HasFullTrigger( void ) const
0088 {
0089     return hasBasicTrigger;
0090 }
0091 
0092 
0093 void  CexmcReconstructor::ReconstructEntryPoints(
0094                                     const CexmcEnergyDepositStore *  edStore )
0095 {
0096     G4int     columnLeft( edStore->calorimeterEDLeftMaxX );
0097     G4int     rowLeft( edStore->calorimeterEDLeftMaxY );
0098     G4int     columnRight( edStore->calorimeterEDRightMaxX );
0099     G4int     rowRight( edStore->calorimeterEDRightMaxY );
0100     G4double  crystalLength( calorimeterGeometry.crystalLength );
0101 
0102     if ( useInnerRefCrystal )
0103     {
0104         TransformToAdjacentInnerCrystal( columnLeft, rowLeft );
0105         TransformToAdjacentInnerCrystal( columnRight, rowRight );
0106     }
0107 
0108     calorimeterEPLeftPosition.setX( 0 );
0109     calorimeterEPLeftPosition.setY( 0 );
0110     calorimeterEPLeftPosition.setZ( -crystalLength / 2 + epDepth );
0111     calorimeterEPLeftDirection.setX( 0 );
0112     calorimeterEPLeftDirection.setY( 0 );
0113     calorimeterEPLeftDirection.setZ( 0 );
0114     calorimeterEPRightPosition.setX( 0 );
0115     calorimeterEPRightPosition.setY( 0 );
0116     calorimeterEPRightPosition.setZ( -crystalLength / 2 + epDepth );
0117     calorimeterEPRightDirection.setX( 0 );
0118     calorimeterEPRightDirection.setY( 0 );
0119     calorimeterEPRightDirection.setZ( 0 );
0120 
0121     G4bool  edInAdjacentCrystalsCollected( false );
0122 
0123     switch ( epDefinitionAlgorithm )
0124     {
0125     case CexmcEntryPointInTheCenter :
0126         break;
0127     case CexmcEntryPointInTheCenterOfCrystalWithMaxED :
0128         {
0129             G4int     nCrystalsInColumn(
0130                                     calorimeterGeometry.nCrystalsInColumn );
0131             G4int     nCrystalsInRow( calorimeterGeometry.nCrystalsInRow );
0132             G4double  crystalWidth( calorimeterGeometry.crystalWidth );
0133             G4double  crystalHeight( calorimeterGeometry.crystalHeight );
0134 
0135             calorimeterEPLeftPosition.setX(
0136                 ( G4double( columnLeft ) - G4double( nCrystalsInRow ) / 2 ) *
0137                                         crystalWidth  + crystalWidth / 2 );
0138             calorimeterEPLeftPosition.setY(
0139                 ( G4double( rowLeft ) - G4double( nCrystalsInColumn ) / 2 ) *
0140                                         crystalHeight + crystalHeight / 2 );
0141             calorimeterEPRightPosition.setX(
0142                 ( G4double( columnRight ) - G4double( nCrystalsInRow ) / 2 ) *
0143                                         crystalWidth  + crystalWidth / 2 );
0144             calorimeterEPRightPosition.setY(
0145                 ( G4double( rowRight ) - G4double( nCrystalsInColumn ) / 2 ) *
0146                                         crystalHeight + crystalHeight / 2 );
0147         }
0148         break;
0149     case CexmcEntryPointByLinearEDWeights :
0150     case CexmcEntryPointBySqrtEDWeights :
0151         {
0152             G4double  x( 0 );
0153             G4double  y( 0 );
0154 
0155             CalculateWeightedEPPosition( edStore->calorimeterEDLeftCollection,
0156                                          rowLeft, columnLeft, x, y,
0157                                          calorimeterEDLeftAdjacent );
0158             calorimeterEPLeftPosition.setX( x );
0159             calorimeterEPLeftPosition.setY( y );
0160             CalculateWeightedEPPosition( edStore->calorimeterEDRightCollection,
0161                                          rowRight, columnRight, x, y,
0162                                          calorimeterEDRightAdjacent );
0163             calorimeterEPRightPosition.setX( x );
0164             calorimeterEPRightPosition.setY( y );
0165 
0166             if ( csAlgorithm == CexmcSelectAdjacentCrystals )
0167                 edInAdjacentCrystalsCollected = true;
0168         }
0169         break;
0170     default :
0171         break;
0172     }
0173 
0174     switch ( epDepthDefinitionAlgorithm )
0175     {
0176     case CexmcEntryPointDepthPlain :
0177         break;
0178     case CexmcEntryPointDepthSphere :
0179         {
0180             G4double  calorimeterEPLeftRadiusOfTheSphere(
0181                              calorimeterLeftTransform.NetTranslation().mag() +
0182                              calorimeterEPLeftPosition.z() );
0183             G4double  calorimeterEPLeftRadiusOfTheSphere2(
0184                                       calorimeterEPLeftRadiusOfTheSphere *
0185                                       calorimeterEPLeftRadiusOfTheSphere );
0186             G4double  calorimeterEPLeftPositionX2(
0187                                             calorimeterEPLeftPosition.x() *
0188                                             calorimeterEPLeftPosition.x() );
0189             G4double  calorimeterEPLeftPositionY2(
0190                                             calorimeterEPLeftPosition.y() *
0191                                             calorimeterEPLeftPosition.y() );
0192             G4double  calorimeterEPLeftPositionZOffset(
0193                            calorimeterEPLeftRadiusOfTheSphere - std::sqrt(
0194                                   calorimeterEPLeftRadiusOfTheSphere2 -
0195                                   calorimeterEPLeftPositionX2 -
0196                                   calorimeterEPLeftPositionY2 ) );
0197             G4double  calorimeterEPRightRadiusOfTheSphere(
0198                               calorimeterRightTransform.NetTranslation().mag() +
0199                               calorimeterEPRightPosition.z() );
0200             G4double  calorimeterEPRightRadiusOfTheSphere2(
0201                                        calorimeterEPRightRadiusOfTheSphere *
0202                                        calorimeterEPRightRadiusOfTheSphere );
0203             G4double  calorimeterEPRightPositionX2(
0204                                             calorimeterEPRightPosition.x() *
0205                                             calorimeterEPRightPosition.x() );
0206             G4double  calorimeterEPRightPositionY2(
0207                                             calorimeterEPRightPosition.y() *
0208                                             calorimeterEPRightPosition.y() );
0209             G4double  calorimeterEPRightPositionZOffset(
0210                             calorimeterEPRightRadiusOfTheSphere - std::sqrt(
0211                                     calorimeterEPRightRadiusOfTheSphere2 -
0212                                     calorimeterEPRightPositionX2 -
0213                                     calorimeterEPRightPositionY2 ) );
0214             calorimeterEPLeftPosition.setZ( calorimeterEPLeftPosition.z() -
0215                                             calorimeterEPLeftPositionZOffset );
0216             calorimeterEPRightPosition.setZ( calorimeterEPRightPosition.z() -
0217                                          calorimeterEPRightPositionZOffset );
0218         }
0219         break;
0220     default :
0221         break;
0222     }
0223 
0224     calorimeterEPLeftWorldPosition = calorimeterLeftTransform.TransformPoint(
0225                                                  calorimeterEPLeftPosition );
0226     calorimeterEPLeftWorldDirection = calorimeterLeftTransform.TransformAxis(
0227                                                  calorimeterEPLeftDirection );
0228     calorimeterEPRightWorldPosition = calorimeterRightTransform.TransformPoint(
0229                                                  calorimeterEPRightPosition );
0230     calorimeterEPRightWorldDirection = calorimeterRightTransform.TransformAxis(
0231                                                  calorimeterEPRightDirection );
0232 
0233     if ( collectEDInAdjacentCrystals && ! edInAdjacentCrystalsCollected )
0234     {
0235         CollectEDInAdjacentCrystals( edStore->calorimeterEDLeftCollection,
0236                                      rowLeft, columnLeft,
0237                                      calorimeterEDLeftAdjacent );
0238         CollectEDInAdjacentCrystals( edStore->calorimeterEDRightCollection,
0239                                      rowRight, columnRight,
0240                                      calorimeterEDRightAdjacent );
0241     }
0242 
0243     hasBasicTrigger = true;
0244 }
0245 
0246 
0247 void  CexmcReconstructor::ReconstructTargetPoint( void )
0248 {
0249     if ( ! targetEPInitialized )
0250     {
0251         targetEPWorldPosition = targetTransform.TransformPoint(
0252                                                     G4ThreeVector( 0, 0, 0 ) );
0253         targetEPWorldDirection.setX( 0 );
0254         targetEPWorldDirection.setY( 0 );
0255         targetEPWorldDirection.setZ( 1 );
0256 
0257         targetEPPosition = targetTransform.Inverse().TransformPoint(
0258                                                     targetEPWorldPosition );
0259         targetEPDirection = targetTransform.Inverse().TransformAxis(
0260                                                     targetEPWorldDirection );
0261         targetEPInitialized = true;
0262     }
0263 
0264     hasBasicTrigger = true;
0265 }
0266 
0267 
0268 void  CexmcReconstructor::ReconstructAngle( void )
0269 {
0270     G4ThreeVector  epLeft( calorimeterEPLeftWorldPosition -
0271                            targetEPWorldPosition );
0272     G4ThreeVector  epRight( calorimeterEPRightWorldPosition -
0273                             targetEPWorldPosition );
0274     theAngle = epLeft.angle( epRight );
0275 
0276     hasBasicTrigger = true;
0277 }
0278 
0279 
0280 void  CexmcReconstructor::CollectEDInAdjacentCrystals(
0281                         const CexmcEnergyDepositCalorimeterCollection &  edHits,
0282                         G4int  row, G4int  column, G4double &  ed )
0283 {
0284     G4int  i( 0 );
0285 
0286     for ( CexmcEnergyDepositCalorimeterCollection::const_iterator
0287                                 k( edHits.begin() ); k != edHits.end(); ++k )
0288     {
0289         if ( i - row > 1 || i - row < -1 )
0290         {
0291             ++i;
0292             continue;
0293         }
0294 
0295         G4int  j( 0 );
0296         for ( CexmcEnergyDepositCrystalRowCollection::const_iterator
0297                   l( k->begin() ); l != k->end(); ++l )
0298         {
0299             if ( j - column > 1 || j - column < -1 )
0300             {
0301                 ++j;
0302                 continue;
0303             }
0304             ed += *l;
0305             ++j;
0306         }
0307         ++i;
0308     }
0309 }
0310 
0311 
0312 void  CexmcReconstructor::CalculateWeightedEPPosition(
0313                 const CexmcEnergyDepositCalorimeterCollection &  edHits,
0314                 G4int  row, G4int  column, G4double &  x, G4double &  y,
0315                 G4double &  ed )
0316 {
0317     G4int     nCrystalsInColumn( calorimeterGeometry.nCrystalsInColumn );
0318     G4int     nCrystalsInRow( calorimeterGeometry.nCrystalsInRow );
0319     G4double  crystalWidth( calorimeterGeometry.crystalWidth );
0320     G4double  crystalHeight( calorimeterGeometry.crystalHeight );
0321 
0322     G4int     i( 0 );
0323     G4double  xWeightsSum( 0 );
0324     G4double  yWeightsSum( 0 );
0325     G4double  energyWeightsSum( 0 );
0326 
0327     if ( csAlgorithm == CexmcSelectAdjacentCrystals )
0328         ed = 0.;
0329 
0330     for ( CexmcEnergyDepositCalorimeterCollection::const_iterator
0331                                 k( edHits.begin() ); k != edHits.end(); ++k )
0332     {
0333         if ( csAlgorithm == CexmcSelectAdjacentCrystals &&
0334                                          ( i - row > 1 || i - row < -1 ) )
0335         {
0336             ++i;
0337             continue;
0338         }
0339 
0340         G4int  j( 0 );
0341         for ( CexmcEnergyDepositCrystalRowCollection::const_iterator
0342                   l( k->begin() ); l != k->end(); ++l )
0343         {
0344             if ( csAlgorithm == CexmcSelectAdjacentCrystals &&
0345                                          ( j - column > 1 || j - column < -1 ) )
0346             {
0347                 ++j;
0348                 continue;
0349             }
0350 
0351             if ( csAlgorithm == CexmcSelectAdjacentCrystals )
0352                 ed += *l;
0353             
0354             G4double  xInCalorimeterOffset(
0355                         ( G4double( j ) - G4double( nCrystalsInRow ) / 2 ) *
0356                         crystalWidth  + crystalWidth / 2 );
0357             G4double  energyWeight(
0358                         epDefinitionAlgorithm ==
0359                                         CexmcEntryPointBySqrtEDWeights ?
0360                                                         std::sqrt( *l ) : *l );
0361             xWeightsSum += energyWeight * xInCalorimeterOffset;
0362             G4double  yInCalorimeterOffset(
0363                         ( G4double( i ) - G4double( nCrystalsInColumn ) / 2 ) *
0364                         crystalHeight  + crystalHeight / 2 );
0365             yWeightsSum += energyWeight * yInCalorimeterOffset;
0366             energyWeightsSum += energyWeight;
0367             ++j;
0368         }
0369         ++i;
0370     }
0371 
0372     x = xWeightsSum / energyWeightsSum;
0373     y = yWeightsSum / energyWeightsSum;
0374 }
0375
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!