CMSSW Reference Manual

TODO: calculate the parameter errors include the beam parameters calculate the global parameters

LASEndcapAlignmentParameterSet LASEndcapAlgorithm::CalculateParameters	(	LASGlobalData< LASCoordinateSet > &	measuredCoordinates,
		LASGlobalData< LASCoordinateSet > &	nominalCoordinates
	)

implementation of the analytical solution for the endcap; described in bruno's thesis (Appendix B): https://darwin.bth.rwth-aachen.de/opus3/volltexte/2002/348 but extended with the beams' phi positions

Definition at line 22 of file LASEndcapAlgorithm.cc.

References funct::cos(), GenMuonPlsPt100GeV_cfg::cout, muonGeometry::disk, lat::endl(), LASEndcapAlignmentParameterSet::GetParameter(), LASCoordinateSet::GetPhi(), LASCoordinateSet::GetR(), LASGlobalData< T >::GetTECEntry(), funct::pow(), radius(), funct::sin(), and LASGlobalLoop::TECLoop().

Referenced by LaserAlignment::endJob().

00023                                                                                                                               {
00024   
00025   std::cout << " [LASEndcapAlgorithm::CalculateParameters] -- Starting." << std::endl;
00026 
00027   // loop object
00028   LASGlobalLoop globalLoop;
00029   int det, ring, beam, disk;
00030 
00031   // vector containing the z positions of the disks in mm;
00032   // outer vector: TEC+/-, inner vector: 9 disks
00033   double zPositions[9] = { 1250., 1390., 1530., 1670., 1810., 1985., 2175., 2380., 2595. };
00034   std::vector<std::vector<double> > diskZPositions( 2, std::vector<double>( 9, 0. ) );
00035   for( det = 0; det < 2; ++ det ) {
00036     for( disk = 0; disk < 9; ++disk ) {
00037       diskZPositions.at( det ).at( disk ) = ( det==0 ? zPositions[disk] : -1. * zPositions[disk] );
00038     }
00039   }
00040   
00041   // constants
00042   const double endcapLength = 1345.; // mm
00043 
00044   // now come some sums which are later used in the formulas for the parameters.
00045   // the rings are implicitly summed over, however, this brings a little complication:
00046   // to make the calculation of the parameters independent of the ring (=radius),
00047   // we define some of the sums twice, once for phi and once for y=r*phi
00048 
00049   // sum over r*phi or phi for each endcap and for each disk (both rings)
00050   // outer vector: TEC+/-, inner vector: 9 disks
00051   std::vector<std::vector<double> > sumOverY( 2, std::vector<double>( 9, 0. ) );
00052   std::vector<std::vector<double> > sumOverPhi( 2, std::vector<double>( 9, 0. ) );
00053 
00054   // double sum over r*phi or phi, for each endcap (both rings)
00055   // outer vector: TEC+/-
00056   std::vector<double> doubleSumOverY( 2, 0. );
00057   std::vector<double> doubleSumOverPhi( 2, 0. );
00058 
00059   // sum over r*phi*z or phi*z, for each endcap (both rings)
00060   // outer vector: TEC+/-
00061   std::vector<double> sumOverYZ( 2, 0. );
00062   std::vector<double> sumOverPhiZ( 2, 0. );
00063 
00064   // sum over sin(phi_nominal)*R*phi for each endcap and for each disk (both rings)
00065   std::vector<std::vector<double> > sumOverSinThetaY( 2, std::vector<double>( 9, 0. ) );
00066 
00067   // sum over cos(phi_nominal)*R*phi for each endcap and for each disk (both rings)
00068   std::vector<std::vector<double> > sumOverCosThetaY( 2, std::vector<double>( 9, 0. ) );
00069 
00070   // double sum over sin or cos(phi_nominal)*phi, for each endcap
00071   std::vector<double> doubleSumOverSinThetaY( 2, 0. );
00072   std::vector<double> doubleSumOverCosThetaY( 2, 0. );
00073 
00074   // double sum over sin or cos(phi_nominal)*phi*z, for each endcap
00075   std::vector<double> doubleSumOverSinThetaYZ( 2, 0. );
00076   std::vector<double> doubleSumOverCosThetaYZ( 2, 0. );
00077 
00078   // sum over z values / sum over z^2 values
00079   std::vector<double> sumOverZ( 2, 0. );
00080   std::vector<double> sumOverZSquared( 2, 0. );
00081 
00082 
00083   // now calculate them
00084   det = 0; ring = 0; beam = 0; disk = 0;
00085   do {
00086 
00087     // current radius, depends on the ring
00088     const double radius = nominalCoordinates.GetTECEntry( det, ring, beam, disk ).GetR();
00089 
00090     // residual in r*phi (in the formulas this corresponds to y_ik)
00091     const double residual = ( measuredCoordinates.GetTECEntry( det, ring, beam, disk ).GetPhi() - nominalCoordinates.GetTECEntry( det, ring, beam, disk ).GetPhi() ) * radius;
00092     
00093     sumOverY.at( det ).at( disk ) += residual;
00094     sumOverPhi.at( det ).at( disk ) += residual / radius;
00095     
00096     doubleSumOverY.at( det ) += residual;
00097     doubleSumOverPhi.at( det ) += residual / radius;
00098     
00099     sumOverYZ.at( det ) += diskZPositions.at( det ).at( disk ) * residual;
00100     sumOverPhiZ.at( det ) +=  diskZPositions.at( det ).at( disk ) * residual / radius;
00101 
00102     sumOverSinThetaY.at( det ).at( disk ) += sin( nominalCoordinates.GetTECEntry( det, ring, beam, disk ).GetPhi() ) * residual;
00103     sumOverCosThetaY.at( det ).at( disk ) += cos( nominalCoordinates.GetTECEntry( det, ring, beam, disk ).GetPhi() ) * residual;
00104 
00105     doubleSumOverSinThetaY.at( det ) += sin( nominalCoordinates.GetTECEntry( det, ring, beam, disk ).GetPhi() ) * residual;
00106     doubleSumOverCosThetaY.at( det ) += cos( nominalCoordinates.GetTECEntry( det, ring, beam, disk ).GetPhi() ) * residual;
00107 
00108     doubleSumOverSinThetaYZ.at( det ) += sin( nominalCoordinates.GetTECEntry( det, ring, beam, disk ).GetPhi() ) * diskZPositions.at( det ).at( disk ) * residual;
00109     doubleSumOverCosThetaYZ.at( det ) += cos( nominalCoordinates.GetTECEntry( det, ring, beam, disk ).GetPhi() ) * diskZPositions.at( det ).at( disk ) * residual;
00110 
00111   } while( globalLoop.TECLoop( det, ring, beam, disk ) );
00112 
00113 
00114   // disk-wise sums
00115   for( disk = 0; disk < 9; ++disk ) {
00116     sumOverZ.at( 0 ) += diskZPositions.at( 0 ).at( disk ); sumOverZ.at( 1 ) += diskZPositions.at( 0 ).at( disk );
00117     sumOverZSquared.at( 0 ) += pow( diskZPositions.at( 0 ).at( disk ), 2 ); sumOverZSquared.at( 1 ) += pow( diskZPositions.at( 0 ).at( disk ), 2 );
00118   }
00119   
00120   
00121 
00122   // now we can calculate the parameters for both TECs simultaneously,
00123   // so they're all vectors( 2 ) for TEC+/- (global parameters), or dim 2*9 (disk parameters)
00124 
00125   // define them..
00126 
00127   // deltaPhi_0
00128   std::vector<double> deltaPhi0( 2, 0. );
00129   
00130   // deltaPhi_t
00131   std::vector<double> deltaPhiT( 2, 0. );
00132   
00133   // deltaPhi_k (k=0..8)
00134   std::vector<std::vector<double> > deltaPhiK( 2, std::vector<double>( 9, 0. ) );
00135 
00136   // deltaX_0
00137   std::vector<double> deltaX0( 2, 0. );
00138   
00139   // deltaX_t
00140   std::vector<double> deltaXT( 2, 0. );
00141   
00142   // deltaX_k (k=0..8)
00143   std::vector<std::vector<double> > deltaXK( 2, std::vector<double>( 9, 0. ) );
00144 
00145   // deltaY_0
00146   std::vector<double> deltaY0( 2, 0. );
00147   
00148   // deltaY_t
00149   std::vector<double> deltaYT( 2, 0. );
00150   
00151   // deltaY_k (k=0..8)
00152   std::vector<std::vector<double> > deltaYK( 2, std::vector<double>( 9, 0. ) );
00153 
00154 
00155   // ..and fill them
00156   for( det = 0; det < 2; ++det ) { // TEC+/- loop
00157     
00158     // deltaPhi_0
00159     // here we use the phi-sums to eliminate the radius
00160     deltaPhi0.at( det ) = ( sumOverZSquared.at( det ) * doubleSumOverPhi.at( det ) - sumOverZ.at( det ) * sumOverPhiZ.at( det ) )
00161       / ( 8. * ( pow( sumOverZ.at( det ), 2 ) - 9. * sumOverZSquared.at( det ) ) );
00162 
00163     // deltaPhi_t
00164     // again use the phi-sums
00165     deltaPhiT.at( det ) = endcapLength * ( 9. * sumOverPhiZ.at( det ) - sumOverZ.at( det ) * doubleSumOverPhi.at( det ) )
00166       / ( 8. * ( pow( sumOverZ.at( det ), 2 ) - 9. * sumOverZSquared.at( det ) ) );
00167 
00168     // deltaPhi_k (k=0..8)
00169     // again use the phi-sums
00170     for( disk = 0; disk < 9; ++disk ) {
00171       deltaPhiK.at( det ).at( disk ) = ( -1. * diskZPositions.at( det ).at( disk ) * deltaPhiT.at( det ) / endcapLength )
00172         -  ( deltaPhi0.at( det ) )  -  sumOverPhi.at( det ).at( disk ) / 8.;
00173     }
00174     
00175     // deltaX_0
00176     deltaX0.at( det ) = 2. * ( sumOverZ.at( det ) * doubleSumOverSinThetaYZ.at( det ) - sumOverZSquared.at( det ) * doubleSumOverSinThetaY.at( det ) )
00177       / ( 8. * ( pow( sumOverZ.at( det ), 2 ) - 9. * sumOverZSquared.at( det ) ) );
00178 
00179     // deltaX_t
00180     deltaXT.at( det ) = 2. * endcapLength * ( sumOverZ.at( det ) * doubleSumOverSinThetaY.at( det ) - 9. * doubleSumOverSinThetaYZ.at( det ) )
00181       / ( 8. * ( pow( sumOverZ.at( det ), 2 ) - 9. * sumOverZSquared.at( det ) ) );
00182 
00183     // deltaX_k (k=0..8)
00184     for( disk = 0; disk < 9; ++disk ) {
00185       deltaXK.at( det ).at( disk ) =  ( -1. * diskZPositions.at( det ).at( disk ) * deltaXT.at( det ) / endcapLength )
00186         -  ( deltaX0.at( det ) )  +  2. * sumOverSinThetaY.at( det ).at( disk ) / 8.;
00187     }
00188 
00189     // deltaY_0
00190     deltaY0.at( det ) = 2. * ( sumOverZSquared.at( det ) * doubleSumOverCosThetaY.at( det ) - sumOverZ.at( det ) * doubleSumOverCosThetaYZ.at( det ) )
00191       / ( 8. * ( pow( sumOverZ.at( det ), 2 ) - 9. * sumOverZSquared.at( det ) ) );
00192 
00193     // deltaY_t
00194     deltaYT.at( det ) = 2. * endcapLength * ( 9. * doubleSumOverCosThetaYZ.at( det ) - sumOverZ.at( det ) * doubleSumOverCosThetaY.at( det ) )
00195       / ( 8. * ( pow( sumOverZ.at( det ), 2 ) - 9. * sumOverZSquared.at( det ) ) );
00196 
00197     // deltaY_k (k=0..8)
00198     for( disk = 0; disk < 9; ++disk ) {
00199       deltaYK.at( det ).at( disk ) =  ( -1. * diskZPositions.at( det ).at( disk ) * deltaYT.at( det ) / endcapLength )
00200         -  ( deltaY0.at( det ) )  -  2. * sumOverCosThetaY.at( det ).at( disk ) / 8.;
00201     }
00202   
00203   }
00204 
00205 
00206   // fill the result
00207   LASEndcapAlignmentParameterSet theResult;
00208 
00209   // for the moment we fill only the values, not the errors
00210 
00211   for( det = 0; det < 2; ++det ) {
00212     for( disk = 0; disk < 9; ++disk ) {
00213       
00214       // the rotation parameters: deltaPhi_k
00215       theResult.GetParameter( det, disk, 0 ).first = deltaPhiK.at( det ).at( disk );
00216 
00217       // the x offsets: deltaX_k
00218       theResult.GetParameter( det, disk, 1 ).first = deltaXK.at( det ).at( disk );
00219       
00220       // the y offsets: deltaY_k
00221       theResult.GetParameter( det, disk, 2 ).first = deltaYK.at( det ).at( disk );
00222 
00223     }
00224   }
00225 
00226   std::cout << " [LASEndcapAlgorithm::CalculateParameters] -- Done." << std::endl;
00227 
00228   return( theResult );
00229 
00230 }


Public Member Functions
LASEndcapAlignmentParameterSet	CalculateParameters (LASGlobalData< LASCoordinateSet > &, LASGlobalData< LASCoordinateSet > &)
	implementation of the analytical solution for the endcap; described in bruno's thesis (Appendix B): https://darwin.bth.rwth-aachen.de/opus3/volltexte/2002/348 but extended with the beams' phi positions
	LASEndcapAlgorithm ()

LASEndcapAlgorithm Class Reference

Public Member Functions

Detailed Description

Constructor & Destructor Documentation

Member Function Documentation