Skip to content

GitLab

Commit 135f53a2 authored by Claus Kleinwort's avatar Claus Kleinwort
Browse files

examples cleaned up 

git-svn-id: http://svnsrv.desy.de/public/GeneralBrokenLines/trunk@146 281f6f2b-e318-4fd1-8bce-1a4ba7aab212
parent 48c1f948
Hide whitespace changes
Inline Side-by-side
cpp/CMakeLists.txt
View file @ 135f53a2
......@@ -7,7 +7,7 @@ PROJECT(GBL)
# project version
SET( ${PROJECT_NAME}_VERSION_MAJOR 2 )
SET( ${PROJECT_NAME}_VERSION_MINOR 1 )
SET( ${PROJECT_NAME}_VERSION_PATCH 5 )
SET( ${PROJECT_NAME}_VERSION_PATCH 6 )
# make life easier and simply use the ilcsoft default settings
# load default ilcsoft settings (install prefix, build type, rpath, etc.)
......
cpp/examples/exampleDc.cpp
View file @ 135f53a2
......@@ -87,13 +87,13 @@ void exampleDc() {
double dist = 200.;
for (unsigned int iLayer = 0; iLayer < 6; ++iLayer) {
layers.push_back(
CreateLayerDc("CH1+", dist * sinTheta, 0.0, dist * cosTheta,
CreateLayerDc("CH1+", 1, dist * sinTheta, 0.0, dist * cosTheta,
thickness[iLayer], theta, 6., 0.030)); // +6 deg stereo layers
dist += 2.;
}
for (unsigned int iLayer = 6; iLayer < 12; ++iLayer) {
layers.push_back(
CreateLayerDc("CH1-", dist * sinTheta, 0.0, dist * cosTheta,
CreateLayerDc("CH1-", 1, dist * sinTheta, 0.0, dist * cosTheta,
thickness[iLayer], theta, -6., 0.030)); // -6 deg stereo layers
dist += 2.;
}
......@@ -101,13 +101,13 @@ void exampleDc() {
dist = 300.;
for (unsigned int iLayer = 0; iLayer < 6; ++iLayer) {
layers.push_back(
CreateLayerDc("CH2+", dist * sinTheta, 0.0, dist * cosTheta,
CreateLayerDc("CH2+", 2, dist * sinTheta, 0.0, dist * cosTheta,
thickness[iLayer], theta, 6., 0.030)); // +6 deg stereo layers
dist += 2.;
}
for (unsigned int iLayer = 6; iLayer < 12; ++iLayer) {
layers.push_back(
CreateLayerDc("CH2-", dist * sinTheta, 0.0, dist * cosTheta,
CreateLayerDc("CH2-", 2, dist * sinTheta, 0.0, dist * cosTheta,
thickness[iLayer], theta, -6., 0.030)); // -6 deg stereo layers
dist += 2.;
}
......@@ -115,13 +115,13 @@ void exampleDc() {
dist = 400.;
for (unsigned int iLayer = 0; iLayer < 6; ++iLayer) {
layers.push_back(
CreateLayerDc("CH3+", dist * sinTheta, 0.0, dist * cosTheta,
CreateLayerDc("CH3+", 3, dist * sinTheta, 0.0, dist * cosTheta,
thickness[iLayer], theta, 6., 0.030)); // +6 deg stereo layers
dist += 2.;
}
for (unsigned int iLayer = 6; iLayer < 12; ++iLayer) {
layers.push_back(
CreateLayerDc("CH3-", dist * sinTheta, 0.0, dist * cosTheta,
CreateLayerDc("CH3-", 3, dist * sinTheta, 0.0, dist * cosTheta,
thickness[iLayer], theta, -6., 0.030)); // -6 deg stereo layers
dist += 2.;
}
......@@ -240,8 +240,9 @@ void exampleDc() {
Matrix<double, 2, 6> derGlobal = layer.getRigidBodyDerLocal(pos,
dir); */
// Chamber alignment in global system (as common system for both stereo orientations)
unsigned int layerID = layer.getLayerID();
for (int p = 0; p < 6; p++)
labGlobal[p] = (iLayer / 12 + 1) * 1000 + p + 1;// chamber alignment
labGlobal[p] = layerID * 1000 + p + 1; // chamber alignment
Matrix<double, 2, 6> derGlobal = layer.getRigidBodyDerGlobal(pos,
dir).block<2, 6>(0, 0);
point.addGlobals(labGlobal, derGlobal);
......@@ -290,6 +291,7 @@ namespace gbl {
/**
* Create drift chamber layer with 1D measurement (u)
* \param [in] aName name
* \param [in] layer layer ID
* \param [in] xPos X-position (of center)
* \param [in] yPos Y-position (of center)
* \param [in] zPos Z-position (of center)
......@@ -298,8 +300,8 @@ namespace gbl {
* \param [in] stereoAngle stereo angle
* \param [in] uRes resolution in u-direction
*/
GblDetectorLayer CreateLayerDc(const std::string aName, double xPos,
double yPos, double zPos, double thickness, double xzAngle,
GblDetectorLayer CreateLayerDc(const std::string aName, unsigned int layer,
double xPos, double yPos, double zPos, double thickness, double xzAngle,
double stereoAngle, double uRes) {
Vector3d aCenter(xPos, yPos, zPos);
Vector2d aResolution(uRes, 0.);
......@@ -313,7 +315,7 @@ GblDetectorLayer CreateLayerDc(const std::string aName, double xPos,
* sinXz, sinXz, 0., cosXz; // U,V,N
Matrix3d alignTrafo;
alignTrafo << cosXz, 0., -sinXz, 0., 1., 0., sinXz, 0., cosXz; // I,J,K
return GblDetectorLayer(aName, 1, thickness, aCenter, aResolution,
return GblDetectorLayer(aName, layer, 1, thickness, aCenter, aResolution,
aPrecision, measTrafo, alignTrafo);
}
......
cpp/examples/exampleDc.h
View file @ 135f53a2
......@@ -35,8 +35,8 @@
//! Namespace for the general broken lines package
namespace gbl {
GblDetectorLayer CreateLayerDc(const std::string aName, double xPos,
double yPos, double zPos, double thickness, double xzAngle,
GblDetectorLayer CreateLayerDc(const std::string aName, unsigned int layer,
double xPos, double yPos, double zPos, double thickness, double xzAngle,
double stereoAngle, double uRes);
}
......
cpp/examples/exampleSit.cpp
View file @ 135f53a2
......@@ -82,27 +82,28 @@ void exampleSit() {
// name, position (x,y,z), thickness (X/X_0), (1 or 2) measurements (direction in YZ, resolution)
std::vector<GblDetectorLayer> layers;
layers.push_back(
CreateLayerSit("PIX1", 2.0, 0., 0., 0.0033, 0., 0.0010, 90.,
CreateLayerSit("PIX1", 0, 2.0, 0., 0., 0.0033, 0., 0.0010, 90.,
0.0020)); // pixel
layers.push_back(
CreateLayerSit("PIX2", 3.0, 0., 0., 0.0033, 0., 0.0010, 90.,
CreateLayerSit("PIX2", 1, 3.0, 0., 0., 0.0033, 0., 0.0010, 90.,
0.0020)); // pixel
layers.push_back(
CreateLayerSit("PIX3", 4.0, 0., 0., 0.0033, 0., 0.0010, 90.,
CreateLayerSit("PIX3", 2, 4.0, 0., 0., 0.0033, 0., 0.0010, 90.,
0.0020)); // pixel
layers.push_back(
CreateLayerSit("S2D4", 6.0, 0., 0., 0.0033, 0., 0.0025, 5.0,
CreateLayerSit("S2D4", 3, 6.0, 0., 0., 0.0033, 0., 0.0025, 5.0,
0.0025)); // strip 2D, +5 deg stereo angle
layers.push_back(
CreateLayerSit("S2D5", 8.0, 0., 0., 0.0033, 0., 0.0025, -5.,
CreateLayerSit("S2D5", 4, 8.0, 0., 0., 0.0033, 0., 0.0025, -5.,
0.0025)); // strip 2D, -5 deg stereo angle
layers.push_back(
CreateLayerSit("S2D6", 10., 0., 0., 0.0033, 0., 0.0025, 5.0,
CreateLayerSit("S2D6", 5, 10., 0., 0., 0.0033, 0., 0.0025, 5.0,
0.0025)); // strip 2D, +5 deg stereo angle
layers.push_back(
CreateLayerSit("S2D7", 12., 0., 0., 0.0033, 0., 0.0025, -5.,
CreateLayerSit("S2D7", 6, 12., 0., 0., 0.0033, 0., 0.0025, -5.,
0.0025)); // strip 2D, -5 deg stereo angle
layers.push_back(CreateLayerSit("S1D8", 15., 0., 0., 0.0033, 0., 0.0040)); // strip 1D, no sensitivity to Z
layers.push_back(
CreateLayerSit("S1D8", 7, 15., 0., 0., 0.0033, 0., 0.0040)); // strip 1D, no sensitivity to Z
/* print layers
for (unsigned int iLayer = 0; iLayer < layers.size(); ++iLayer) {
......@@ -210,8 +211,9 @@ void exampleSit() {
point.addMeasurement(proL2m, res, measPrecision);
// global labels and parameters for rigid body alignment
std::vector<int> labGlobal(6);
unsigned int layerID = layer.getLayerID();
for (int p = 0; p < 6; p++)
labGlobal[p] = iLayer * 10 + p + 1;
labGlobal[p] = layerID * 10 + p + 1;
Vector3d pos = pred.getPosition();
Vector3d dir = pred.getDirection();
Matrix<double, 2, 6> derGlobal = layer.getRigidBodyDerLocal(pos,
......@@ -263,6 +265,7 @@ namespace gbl {
* Create silicon layer with 1D measurement (u) at fixed X-position.
*
* \param [in] aName name
* \param [in] layer layer ID
* \param [in] xPos X-position (of center)
* \param [in] yPos Y-position (of center)
* \param [in] zPos Z-position (of center)
......@@ -270,8 +273,8 @@ namespace gbl {
* \param [in] uAngle angle of u-direction in YZ plane
* \param [in] uRes resolution in u-direction
*/
GblDetectorLayer CreateLayerSit(const std::string aName, double xPos,
double yPos, double zPos, double thickness, double uAngle,
GblDetectorLayer CreateLayerSit(const std::string aName, unsigned int layer,
double xPos, double yPos, double zPos, double thickness, double uAngle,
double uRes) {
Vector3d aCenter(xPos, yPos, zPos);
Vector2d aResolution(uRes, 0.);
......@@ -282,7 +285,7 @@ GblDetectorLayer CreateLayerSit(const std::string aName, double xPos,
measTrafo << 0., cosU, sinU, 0., -sinU, cosU, 1., 0., 0.; // U,V,N
Matrix3d alignTrafo;
alignTrafo << 0., 1., 0., 0., 0., 1., 1., 0., 0.; // Y,Z,X
return GblDetectorLayer(aName, 1, thickness, aCenter, aResolution,
return GblDetectorLayer(aName, layer, 1, thickness, aCenter, aResolution,
aPrecision, measTrafo, alignTrafo);
}
......@@ -293,6 +296,7 @@ GblDetectorLayer CreateLayerSit(const std::string aName, double xPos,
* (but must be different).
*
* \param [in] aName name
* \param [in] layer layer ID
* \param [in] xPos X-position (of center)
* \param [in] yPos Y-position (of center)
* \param [in] zPos Z-position (of center)
......@@ -302,9 +306,9 @@ GblDetectorLayer CreateLayerSit(const std::string aName, double xPos,
* \param [in] vAngle angle of v-direction in YZ plane
* \param [in] vRes resolution in v-direction
*/
GblDetectorLayer CreateLayerSit(const std::string aName, double xPos,
double yPos, double zPos, double thickness, double uAngle, double uRes,
double vAngle, double vRes) {
GblDetectorLayer CreateLayerSit(const std::string aName, unsigned int layer,
double xPos, double yPos, double zPos, double thickness, double uAngle,
double uRes, double vAngle, double vRes) {
Vector3d aCenter(xPos, yPos, zPos);
Vector2d aResolution(uRes, vRes);
Vector2d aPrecision(1. / (uRes * uRes), 1. / (vRes * vRes));
......@@ -316,7 +320,7 @@ GblDetectorLayer CreateLayerSit(const std::string aName, double xPos,
measTrafo << 0., cosU, sinU, 0., cosV, sinV, 1., 0., 0.; // U,V,N
Matrix3d alignTrafo;
alignTrafo << 0., 1., 0., 0., 0., 1., 1., 0., 0.; // Y,Z,X
return GblDetectorLayer(aName, 2, thickness, aCenter, aResolution,
return GblDetectorLayer(aName, layer, 2, thickness, aCenter, aResolution,
aPrecision, measTrafo, alignTrafo);
}
......
cpp/examples/exampleSit.h
View file @ 135f53a2
......@@ -35,11 +35,12 @@
//! Namespace for the general broken lines package
namespace gbl {
GblDetectorLayer CreateLayerSit(const std::string aName, double xPos,
double yPos, double zPos, double thickness, double uAngle, double uRes); // 1D measurement
GblDetectorLayer CreateLayerSit(const std::string aName, double xPos,
double yPos, double zPos, double thickness, double uAngle, double uRes,
double vAngle, double vRes); // 2D measurement
GblDetectorLayer CreateLayerSit(const std::string aName, unsigned int layer,
double xPos, double yPos, double zPos, double thickness, double uAngle,
double uRes); // 1D measurement
GblDetectorLayer CreateLayerSit(const std::string aName, unsigned int layer,
double xPos, double yPos, double zPos, double thickness, double uAngle,
double uRes, double vAngle, double vRes); // 2D measurement
}
......
cpp/examples/exampleUtil.cpp
View file @ 135f53a2
......@@ -174,6 +174,43 @@ GblSimpleHelix::GblSimpleHelix(double aRinv, double aPhi0, double aDca,
GblSimpleHelix::~GblSimpleHelix() {
}
/// Get phi (of point on circle) for given radius (to ref. point)
/**
* ( |dca| < radius < |rad-2*dca|, from H1/cjfphi, not restricted to -Pi .. +Pi )
*
* \param[in] aRadius radius
*/
double GblSimpleHelix::getPhi(double aRadius) const {
double arg = (0.5 * rinv * (aRadius * aRadius + dca * dca) - dca)
/ (aRadius * (1.0 - rinv * dca));
return asin(arg) + phi0;
}
/// Get (2D) arc length for given radius (to ref. point)
/**
* ( |dca| < radius < |rad-2*dca|, from H1/cjfsxy )
*
* \param[in] aRadius radius
*/
double GblSimpleHelix::getArcLengthR(double aRadius) const {
double arg = (0.5 * rinv * (aRadius * aRadius + dca * dca) - dca)
/ (aRadius * (1.0 - rinv * dca));
if (fabs(arg) >= 1.) {
std::cout << " bad arc " << aRadius << " " << rinv << " " << dca
<< std::endl;
return 0.;
}
// line
if (rinv == 0)
return sqrt(aRadius * aRadius - dca * dca);
// helix
double sxy = asin(aRadius * rinv * sqrt(1.0 - arg * arg)) / rinv;
if (0.5 * rinv * rinv * (aRadius * aRadius - dca * dca) - 1. + rinv * dca
> 0.)
sxy = M_PI / fabs(rinv) - sxy;
return sxy;
}
/// Get (2D) arc length for given point.
/**
* \param [in] xPos X Position
......@@ -299,6 +336,7 @@ GblHelixPrediction GblSimpleHelix::getPrediction(const Eigen::Vector3d& refPos,
/**
* Create detector layer with 1D or 2D measurement (u,v)
* \param [in] aName name
* \param [in] aLayer layer ID
* \param [in] aDim dimension (1,2)
* \param [in] thickness thickness / radiation_length
* \param [in] aCenter center of detector (origin of local systems)
......@@ -307,13 +345,14 @@ GblHelixPrediction GblSimpleHelix::getPrediction(const Eigen::Vector3d& refPos,
* \param [in] measTrafo matrix of row vectors defining local measurement system
* \param [in] alignTrafo matrix of row vectors defining local alignment system
*/
GblDetectorLayer::GblDetectorLayer(const std::string aName, const int aDim,
const double thickness, Eigen::Vector3d& aCenter,
Eigen::Vector2d& aResolution, Eigen::Vector2d& aPrecision,
Eigen::Matrix3d& measTrafo, Eigen::Matrix3d& alignTrafo) :
name(aName), measDim(aDim), xbyx0(thickness), center(aCenter), resolution(
aResolution), precision(aPrecision), global2meas(measTrafo), global2align(
alignTrafo) {
GblDetectorLayer::GblDetectorLayer(const std::string aName,
const unsigned int aLayer, const int aDim, const double thickness,
Eigen::Vector3d& aCenter, Eigen::Vector2d& aResolution,
Eigen::Vector2d& aPrecision, Eigen::Matrix3d& measTrafo,
Eigen::Matrix3d& alignTrafo) :
name(aName), layer(aLayer), measDim(aDim), xbyx0(thickness), center(
aCenter), resolution(aResolution), precision(aPrecision), global2meas(
measTrafo), global2align(alignTrafo) {
udir = global2meas.row(0);
vdir = global2meas.row(1);
ndir = global2meas.row(2);
......@@ -325,13 +364,18 @@ GblDetectorLayer::~GblDetectorLayer() {
/// Print GblSiliconlayer.
void GblDetectorLayer::print() const {
IOFormat CleanFmt(4, 0, ", ", "\n", "[", "]");
std::cout << " Layer " << name << " " << measDim << "D, " << xbyx0
<< " X0, @ " << center.transpose().format(CleanFmt) << ", res "
<< resolution.transpose().format(CleanFmt) << ", udir "
std::cout << " Layer " << name << " " << layer << " : " << measDim << "D, "
<< xbyx0 << " X0, @ " << center.transpose().format(CleanFmt)
<< ", res " << resolution.transpose().format(CleanFmt) << ", udir "
<< udir.transpose().format(CleanFmt) << ", vdir "
<< vdir.transpose().format(CleanFmt) << std::endl;
}
/// Get layer ID
unsigned int GblDetectorLayer::getLayerID() const {
return layer;
}
/// Get radiation length.
double GblDetectorLayer::getRadiationLength() const {
return xbyx0;
......@@ -391,8 +435,10 @@ const Matrix<double, 3, 6> GblDetectorLayer::getRigidBodyDerGlobal(
return global2meas * drdm * dmdg;
}
/// Get rigid body derivatives in local (measurement) frame.
/// Get rigid body derivatives in local (alignment) frame (with N=(0,0,1)).
/**
* Normal to measurement plane has to be (0,0,1) in local frame.
*
* \param[in] position position (of prediction or measurement)
* \param[in] direction track direction
*/
......
cpp/examples/exampleUtil.h
View file @ 135f53a2
......@@ -80,6 +80,8 @@ public:
GblSimpleHelix(double aRinv, double aPhi0, double aDca, double aDzds,
double aZ0);
virtual ~GblSimpleHelix();
double getPhi(double aRadius) const;
double getArcLengthR(double aRadius) const;
double getArcLengthXY(double xPos, double yPos) const;
void moveToXY(double xPos, double yPos, double& newPhi0, double& newDca,
double& newZ0) const;
......@@ -104,12 +106,13 @@ private:
*/
class GblDetectorLayer {
public:
GblDetectorLayer(const std::string aName, const int aDim,
const double thickness, Eigen::Vector3d& aCenter,
GblDetectorLayer(const std::string aName, const unsigned int aLayer,
const int aDim, const double thickness, Eigen::Vector3d& aCenter,
Eigen::Vector2d& aResolution, Eigen::Vector2d& aPrecision,
Eigen::Matrix3d& measTrafo, Eigen::Matrix3d& alignTrafo);
virtual ~GblDetectorLayer();
void print() const;
unsigned int getLayerID() const;
double getRadiationLength() const;
Eigen::Vector2d getResolution() const;
Eigen::Vector2d getPrecision() const;
......@@ -122,7 +125,8 @@ public:
private:
std::string name; ///< name
int measDim; ///< measurement dimension (1 or 2)
unsigned int layer; ///< layer ID
unsigned int measDim; ///< measurement dimension (1 or 2)
double xbyx0; ///< normalized material thickness
Eigen::Vector3d center; ///< center
Eigen::Vector2d resolution; ///< measurements resolution
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment

imprint | privacy | Conditions of Use