Go to the documentation of this file.00001
00014 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
00015 #include "DataFormats/CLHEP/interface/AlgebraicObjects.h"
00016 #include "DataFormats/GeometrySurface/interface/TkRotation.h"
00017
00018 #include <iostream>
00019
00020
00021 namespace reco {
00022
00023 using namespace math;
00024
00025 BeamSpot::BeamSpot() {
00026
00027 position_ = Point(0.,0.,0.);
00028 sigmaZ_ = 0.;
00029 dxdz_ = 0.;
00030 dydz_ = 0.;
00031 BeamWidthX_ = 0.; BeamWidthY_ = 0;
00032 for (int j=0; j<7; j++) {
00033 for (int k=j; k<7; k++) {
00034 error_(j,k) = 0.;
00035 }
00036 }
00037 type_ = Unknown;
00038 emittanceX_ = 0;
00039 emittanceY_ = 0;
00040 betaStar_ = 0;
00041 }
00042
00043 void BeamSpot::print(std::stringstream& ss) const {
00044
00045 ss << "-----------------------------------------------------\n"
00046 << " Beam Spot Data\n\n"
00047 << " Beam type = " << type() << "\n"
00048 << " X0 = " << x0() << " +/- " << x0Error() << " [cm]\n"
00049 << " Y0 = " << y0() << " +/- " << y0Error() << " [cm]\n"
00050 << " Z0 = " << z0() << " +/- " << z0Error() << " [cm]\n"
00051 << " Sigma Z0 = " << sigmaZ() << " +/- " << sigmaZ0Error() << " [cm]\n"
00052 << " dxdz = " << dxdz() << " +/- " << dxdzError() << " [radians]\n"
00053 << " dydz = " << dydz() << " +/- " << dydzError() << " [radians]\n"
00054 << " Beam width X = " << BeamWidthX() << " +/- " << BeamWidthXError() << " [cm]\n"
00055 << " Beam width Y = " << BeamWidthY() << " +/- " << BeamWidthYError() << " [cm]\n"
00056 << " EmittanceX = " << emittanceX() << " [cm]\n"
00057 << " EmittanceY = " << emittanceY() << " [cm]\n"
00058 << " beta-star = " << betaStar() << " [cm]\n"
00059 << "-----------------------------------------------------\n\n";
00060
00061 }
00062
00063
00064 std::ostream& operator<< ( std::ostream& os, BeamSpot beam ) {
00065 std::stringstream ss;
00066 beam.print(ss);
00067 os << ss.str();
00068 return os;
00069 }
00070
00071 BeamSpot::Covariance3DMatrix BeamSpot::rotatedCovariance3D() const
00072 {
00073 AlgebraicVector3 newZ(dxdz(), dydz(), 1.);
00074 AlgebraicVector3 globalZ(0.,0.,1.);
00075 AlgebraicVector3 rotationAxis = ROOT::Math::Cross(globalZ.Unit(), newZ.Unit());
00076 float rotationAngle = -acos( ROOT::Math::Dot(globalZ.Unit(),newZ.Unit()));
00077 AlgebraicVector a = asHepVector(rotationAxis);
00078 Basic3DVector<float> aa(a[0], a[1], a[2]);
00079 TkRotation<float> rotation(aa ,rotationAngle);
00080 AlgebraicMatrix33 rotationMatrix;
00081 rotationMatrix(0,0) = rotation.xx();
00082 rotationMatrix(0,1) = rotation.xy();
00083 rotationMatrix(0,2) = rotation.xz();
00084 rotationMatrix(1,0) = rotation.yx();
00085 rotationMatrix(1,1) = rotation.yy();
00086 rotationMatrix(1,2) = rotation.yz();
00087 rotationMatrix(2,0) = rotation.zx();
00088 rotationMatrix(2,1) = rotation.zy();
00089 rotationMatrix(2,2) = rotation.zz();
00090
00091 AlgebraicSymMatrix33 diagError ;
00092 diagError(0,0) = pow(BeamWidthX(),2);
00093 diagError(1,1) = pow(BeamWidthY(),2);
00094 diagError(2,2) = pow(sigmaZ(),2);
00095
00096 Covariance3DMatrix matrix;
00097 matrix = ROOT::Math::Similarity(rotationMatrix, diagError) + covariance3D();
00098 return matrix;
00099 }
00100
00101 }