Skip to content

Commit

Permalink
Merge pull request #173 from myliu-hub/master
Browse files Browse the repository at this point in the history 
Shorten initialization time and Solving digital problems
  • Loading branch information
@mirguest
mirguest authored Aug 18, 2021
2 parents 0e915f6 + ca24e1d commit 85375ea
Show file tree
Hide file tree
Showing 7 changed files with 223 additions and 61 deletions.
3 changes: 1 addition & 2 deletions Detector/DetDriftChamber/compact/det.xml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -74,7 +74,7 @@

<limits>
<limitset name="DC_limits">
<limit name="step_length_max" particles="*" value="0.5" unit="mm" />
<limit name="step_length_max" particles="*" value="0.1" unit="mm" />
</limitset>
</limits>

Expand Down Expand Up @@ -113,7 +113,6 @@
<segmentation type="GridDriftChamber" cell_size="SDT_chamber_cell_width" epsilon0="Epsilon" detector_length="DC_length" identifier_phi="cellID" DC_rbegin="DC_chamber_layer_rbegin" DC_rend="DC_chamber_layer_rend" DC_rmin="SDT_chamber_radius_min" DC_rmax="SDT_chamber_radius_max" safe_distance="DC_safe_distance" layerID="layer" layer_width="SDT_chamber_layer_width"/>


<!-- <id>system:8,chamber:1,layer:8,cellID:16</id> -->
<id>system:5,layer:7:9,chamber:8,cellID:32:16</id>
</readout>
</readouts>
Expand Down
30 changes: 19 additions & 11 deletions Detector/DetDriftChamber/src/driftchamber/DriftChamber.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -59,6 +59,7 @@ static dd4hep::Ref_t create_detector(dd4hep::Detector& theDetector,
double chamber_layer_rend = theDetector.constant<double>("DC_chamber_layer_rend");
int chamber_layer_number = floor((chamber_layer_rend-chamber_layer_rbegin)/chamber_layer_width);

double safe_distance = theDetector.constant<double>("DC_safe_distance");
double epsilon = theDetector.constant<double>("Epsilon");

// =======================================================================
Expand All @@ -83,12 +84,6 @@ static dd4hep::Ref_t create_detector(dd4hep::Detector& theDetector,
// - chamber volume
dd4hep::Tube det_chamber_solid(chamber_radius_min, chamber_radius_max, chamber_half_length);
dd4hep::Volume det_chamber_vol(det_name+"_chamber_vol", det_chamber_solid, chamber_mat);
if ( x_det.isSensitive() ) {
det_chamber_vol.setRegion(theDetector,x_det.regionStr());
det_chamber_vol.setLimitSet(theDetector,x_det.limitsStr());
det_chamber_vol.setSensitiveDetector(sens);
sd.setType("tracker");
}

// - wall
double chamber_inner_wall_rmin = theDetector.constant<double>("SDT_chamber_inner_wall_radius_min");
Expand Down Expand Up @@ -177,6 +172,18 @@ static dd4hep::Ref_t create_detector(dd4hep::Detector& theDetector,
DCHseg->setGeomParams(chamber_id, layerIndex, layer_Phi, rmid, epsilon, offset);
DCHseg->setWiresInLayer(chamber_id, layerIndex, numWire);


dd4hep::Tube layer_vol_solid(rmin,rmax,chamber_half_length);
dd4hep::Volume layer_vol(det_name+"_layer_vol",layer_vol_solid,det_mat);
current_vol_ptr = &layer_vol;

if ( x_det.isSensitive() ) {
layer_vol.setRegion(theDetector,x_det.regionStr());
layer_vol.setLimitSet(theDetector,x_det.limitsStr());
layer_vol.setSensitiveDetector(sens);
sd.setType("tracker");
}

// - wire vol
//phi <-------------------> -phi
// | F8 F7 F6 F5| Only on the outermost layer.
Expand All @@ -185,15 +192,14 @@ static dd4hep::Ref_t create_detector(dd4hep::Detector& theDetector,
// | |
// | F0 F1 F2 F3|
// -----------------------
// if(layer_id == 0 || layer_id == 1 || layer_id == 2 || layer_id == 99) {
for(int icell=0; icell< numWire; icell++) {
double wire_phi = (icell+0.5)*layer_Phi + offset;
// - signal wire
dd4hep::Transform3D transform_module(dd4hep::Rotation3D(),dd4hep::Position(rmid*std::cos(wire_phi),rmid*std::sin(wire_phi),0.));
dd4hep::PlacedVolume module_phy = (*current_vol_ptr).placeVolume(module_vol,transform_module);
// - Field wire
dd4hep::PlacedVolume Module_phy;
double radius[9] = {rmid-chamber_layer_width*0.5,rmid-chamber_layer_width*0.5,rmid-chamber_layer_width*0.5,rmid-chamber_layer_width*0.5,rmid,rmid+chamber_layer_width*0.5,rmid+chamber_layer_width*0.5,rmid+chamber_layer_width*0.5,rmid+chamber_layer_width*0.5};
double radius[9] = {rmid-chamber_layer_width*0.5+safe_distance,rmid-chamber_layer_width*0.5+safe_distance,rmid-chamber_layer_width*0.5+safe_distance,rmid-chamber_layer_width*0.5+safe_distance,rmid,rmid+chamber_layer_width*0.5-safe_distance,rmid+chamber_layer_width*0.5-safe_distance,rmid+chamber_layer_width*0.5-safe_distance,rmid+chamber_layer_width*0.5-safe_distance};
double phi[9] = {wire_phi+layer_Phi*0.25,wire_phi,wire_phi-layer_Phi*0.25,wire_phi-layer_Phi*0.5,wire_phi-layer_Phi*0.5,wire_phi-layer_Phi*0.5,wire_phi-layer_Phi*0.25,wire_phi,wire_phi+layer_Phi*0.25};
int num = 5;
if(layer_id==(chamber_layer_number-1)) {
Expand All @@ -206,9 +212,11 @@ static dd4hep::Ref_t create_detector(dd4hep::Detector& theDetector,
Module_phy = (*current_vol_ptr).placeVolume(Module_vol,transform_Module);
}
}
}

// }
dd4hep::Transform3D transform_layer(dd4hep::Rotation3D(),
dd4hep::Position(0,0,0));
dd4hep::PlacedVolume layer_phy = det_chamber_vol.placeVolume(layer_vol ,transform_layer);
layer_phy.addPhysVolID("layer", layer_id);
}

// - place in det
// - chamber
Expand Down
13 changes: 5 additions & 8 deletions Detector/DetSegmentation/include/DetSegmentation/GridDriftChamber.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -61,6 +61,11 @@ class GridDriftChamber : public Segmentation {
const VolumeID& aVolumeID) const;
virtual double distanceTrackWire(const CellID& cID, const TVector3& hit_start, const TVector3& hit_end) const;
virtual void cellposition(const CellID& cID, TVector3& Wstart, TVector3& Wend) const;
TVector3 LineLineIntersect(TVector3& p1, TVector3& p2, TVector3& p3, TVector3& p4) const;
virtual TVector3 distanceClosestApproach(const CellID& cID, const TVector3& hitPos) const;
virtual TVector3 Line_TrackWire(const CellID& cID, const TVector3& hit_start, const TVector3& hit_end) const;
virtual TVector3 IntersectionTrackWire(const CellID& cID, const TVector3& hit_start, const TVector3& hit_end) const;
virtual TVector3 wirePos_vs_z(const CellID& cID, const double& zpos) const;

// double phi(const CellID& cID) const;
inline double cell_Size() const { return m_cellSize; }
Expand Down Expand Up @@ -109,14 +114,6 @@ class GridDriftChamber : public Segmentation {

inline auto returnAllWires() const { return m_wiresPositions; }

// inline TVector3 returnWirePosition(double angle, int sign) const {
// TVector3 w(0, 0, 0);
// w.SetX(_currentRadius * std::cos(angle));
// w.SetY(_currentRadius * std::sin(angle));
// w.SetZ(sign * m_detectorLength / 2.0);
// return w;
// }

void updateParams(int chamber, int layer) const{
auto it_end = layer_params.cend();
--it_end;
Expand Down
136 changes: 132 additions & 4 deletions Detector/DetSegmentation/src/GridDriftChamber.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -102,29 +102,157 @@ void GridDriftChamber::cellposition(const CellID& cID, TVector3& Wstart,
Wend = returnWirePosition(phi_end, 1);
}

TVector3 GridDriftChamber::LineLineIntersect(TVector3& p1, TVector3& p2, TVector3& p3, TVector3& p4) const {

TVector3 p13, p43, p21;
double d1343, d4321, d1321, d4343, d2121;
double numer, denom;
double mua, mub;
TVector3 pa, pb;

p13.SetX(p1.X() - p3.X());
p13.SetY(p1.Y() - p3.Y());
p13.SetZ(p1.Z() - p3.Z());
p43.SetX(p4.X() - p3.X());
p43.SetY(p4.Y() - p3.Y());
p43.SetZ(p4.Z() - p3.Z());
/* if (ABS(p43.X()) < EPS && ABS(p43.Y()) < EPS && ABS(p43.Z()) < EPS) */
/* return(FALSE); */
p21.SetX(p2.X() - p1.X());
p21.SetY(p2.Y() - p1.Y());
p21.SetZ(p2.Z() - p1.Z());
/* if (ABS(p21.X()) < EPS && ABS(p21.Y()) < EPS && ABS(p21.Z()) < EPS) */
/* return(FALSE); */

d1343 = p13.X() * p43.X() + p13.Y() * p43.Y() + p13.Z() * p43.Z();
d4321 = p43.X() * p21.X() + p43.Y() * p21.Y() + p43.Z() * p21.Z();
d1321 = p13.X() * p21.X() + p13.Y() * p21.Y() + p13.Z() * p21.Z();
d4343 = p43.X() * p43.X() + p43.Y() * p43.Y() + p43.Z() * p43.Z();
d2121 = p21.X() * p21.X() + p21.Y() * p21.Y() + p21.Z() * p21.Z();

denom = d2121 * d4343 - d4321 * d4321;
/* if (ABS(denom) < EPS) */
/* return(FALSE); */
numer = d1343 * d4321 - d1321 * d4343;

mua = numer / denom;
mub = (d1343 + d4321 * (mua)) / d4343;

pa.SetX(p1.X() + mua * p21.X());
pa.SetY(p1.Y() + mua * p21.Y());
pa.SetZ(p1.Z() + mua * p21.Z());
pb.SetX(p3.X() + mub * p43.X());
pb.SetY(p3.Y() + mub * p43.Y());
pb.SetZ(p3.Z() + mub * p43.Z());

return pb - pa;
}

double GridDriftChamber::distanceTrackWire(const CellID& cID, const TVector3& hit_start,
const TVector3& hit_end) const {

TVector3 Wstart = {0,0,0};
TVector3 Wend = {0,0,0};
cellposition(cID,Wstart,Wend);

TVector3 a = hit_end - hit_start;
TVector3 b = Wend - Wstart;
TVector3 a = (hit_end - hit_start).Unit();
TVector3 b = (Wend - Wstart).Unit();
TVector3 c = Wstart - hit_start;

double num = std::abs(c.Dot(a.Cross(b)));
double denum = (a.Cross(b)).Mag();

double DCA = 0;

if (denum) {
DCA = num / denum;
if (denum) {
DCA = num / denum;
}

return DCA;
}

TVector3 GridDriftChamber::distanceClosestApproach(const CellID& cID, const TVector3& hitPos) const {
// Distance of the closest approach between a single hit (point) and the closest wire

TVector3 Wstart = {0,0,0};
TVector3 Wend = {0,0,0};
cellposition(cID,Wstart,Wend);

TVector3 temp = (Wend + Wstart);
TVector3 Wmid(temp.X() / 2.0, temp.Y() / 2.0, temp.Z() / 2.0);

double hitPhi = hitPos.Phi();
if (hitPhi < 0) {
hitPhi = hitPhi + 2 * M_PI;
}

TVector3 PCA = Wstart + ((Wend - Wstart).Unit()).Dot((hitPos - Wstart)) * ((Wend - Wstart).Unit());
TVector3 dca = hitPos - PCA;

return dca;
}

TVector3 GridDriftChamber::Line_TrackWire(const CellID& cID, const TVector3& hit_start, const TVector3& hit_end) const {
// The line connecting a particle track to the closest wire
// Returns the vector connecting the both
TVector3 Wstart = {0,0,0};
TVector3 Wend = {0,0,0};
cellposition(cID,Wstart,Wend);

TVector3 P1 = hit_start;
TVector3 P2 = hit_end;
TVector3 P3 = Wstart;
TVector3 P4 = Wend;

TVector3 intersect = LineLineIntersect(P1, P2, P3, P4);
return intersect;
}

// Get the wire position for a z
TVector3 GridDriftChamber::wirePos_vs_z(const CellID& cID, const double& zpos) const {

TVector3 Wstart = {0,0,0};
TVector3 Wend = {0,0,0};
cellposition(cID,Wstart,Wend);

double t = (zpos - Wstart.Z())/(Wend.Z()-Wstart.Z());
double x = Wstart.X()+t*(Wend.X()-Wstart.X());
double y = Wstart.Y()+t*(Wend.Y()-Wstart.Y());

TVector3 wireCoord(x, y, zpos);
return wireCoord;
}

TVector3 GridDriftChamber::IntersectionTrackWire(const CellID& cID, const TVector3& hit_start, const TVector3& hit_end) const {
// Intersection between the particle track and the wire assuming that the track between hit_start and hit_end is linear

TVector3 Wstart = {0,0,0};
TVector3 Wend = {0,0,0};
cellposition(cID,Wstart,Wend);

TVector3 P1 = hit_start;
TVector3 V1 = hit_end-hit_start;

TVector3 P2 = Wstart;
TVector3 V2 = Wend - Wstart;

TVector3 denom = V1.Cross(V2);
double mag_denom = denom.Mag();

TVector3 intersect(0, 0, 0);

if (mag_denom !=0)
{
TVector3 num = ((P2-P1)).Cross(V2);
double mag_num = num.Mag();
double a = mag_num / mag_denom;

intersect = P1 + a * V1;

}
return intersect;
}


}
}
Loading

0 comments on commit 85375ea

Please sign in to comment.