lcsim/src/org/lcsim/recon/muon
diff -u -r1.2 -r1.3
--- MuonCalSegmentFinder.java 31 Oct 2005 22:31:25 -0000 1.2
+++ MuonCalSegmentFinder.java 6 Dec 2005 16:16:39 -0000 1.3
@@ -1,11 +1,19 @@
package org.lcsim.recon.muon;
+import java.util.Vector;
+import java.util.List;
+import java.util.Map;
+import java.util.HashMap;
import org.lcsim.geometry.*;
import org.lcsim.util.step.TrackStepper; // C.M. Stepper
import org.lcsim.event.EventHeader;
import org.lcsim.geometry.field.Solenoid;
import org.lcsim.geometry.FieldMap;
-
+import org.lcsim.geometry.subdetector.MultiLayerTracker;
+import org.lcsim.geometry.layer.LayerStack;
+import org.lcsim.geometry.layer.LayerSlice;
+import org.lcsim.event.SimCalorimeterHit;
+import org.lcsim.recon.cluster.density.MyTools;
/**
* Muon Segment Finder for muon calorimeter
@@ -18,35 +26,42 @@
final public class MuonCalSegmentFinder extends BarrelCalSegmentFinder
{
- public MuonCalSegmentFinder(Detector det, String name)
+ public MuonCalSegmentFinder(Detector det, String subdetectorName)
{
- super(det, name);
+ super(det, subdetectorName);
// Set the nearest neighbor cuts to default values
setPhiNNCut(3);
setThetaNNCut(2);
+ // getting info from the coil
+ MultiLayerTracker coilsub = (MultiLayerTracker)det.getSubdetectors().get("SolenoidCoilBarrel");
- // The following kludge is required to get inner radius of the
- // solenoid because the Coil class does not implement the
- // DetectorComponent interface. Hence we must look through all
- // components to find the coil
- coilRadius = 2500.0;
- // FIXME: need to get coilRadius from some detector class...
-// if(detName.contains("SDJan03")) coilRadius = 271.0;
-
-// for (Enumeration e = det.getComponents(); e.hasMoreElements();)
-// {
-// Object detComp = (Object)e.nextElement();
-// if (detComp.getClass().getName().equals(
-// "org.lcsim.geometry.component.Coil"))
-// {
-// Coil coil = (Coil)detComp;
-// coilRadius = coil.getDouble("r_inner");
-// breacoilk;
-// }
-// }
-
+ /**
+ * Get materials from MultiLayerTracker class.
+ *
+ * WARNING: This approach is prone to error, because the slices are
+ * not guaranteed to be the same thickness in every layer.
+ */
+ LayerStack layers = coilsub.getLayering().getLayers();
+ double[] rmin = coilsub.getInnerR();
+ double[] zmax = coilsub.getOuterZ();
+
+ int nsublayers = layers.getNumberOfLayers();
+ coilSlices = new Vector<CoilSubLayer>();
+ for(int i=0; i<nsublayers; ++i) {
+ List<LayerSlice> slices
+ = coilsub.getLayering().getLayers().getLayer(i).getSlices();
+ double raux = rmin[i];
+ for ( Object o : slices ) {
+ LayerSlice ls = (LayerSlice) o;
+ String matName = ls.getMaterial().getName();
+ double thickness = ls.getThickness();
+ double matdedx = dedx.getDeDx(matName)/10; // Eloss this material
+ coilSlices.add( new CoilSubLayer(matName, thickness, raux, zmax[i], matdedx) );
+ raux += thickness;
+ }
+ }
}
@@ -65,7 +80,7 @@
numLayersHit = 0;
double rr=0.;
double zmxCoil= 0.;
- FieldMap fm = det.getFieldMap();
+ FieldMap fm = det.getFieldMap();
if( fm instanceof Solenoid)
{
Solenoid sol = (Solenoid) fm;
@@ -73,16 +88,16 @@
rr = Math.sqrt(sol.getOuterRadius2());
System.out.println("MCSF: Coil outer_radius="+rr+" zmxCoil="+zmxCoil);
}
- //coilRadius = rr;
+
// Forms a muon segment by looking for calorimeter hits that match
// cells which the track passes through
rpVect = stpr.getNewRp();
- //First in the Air Before the COIL
-
+
+ //First in the Air Before the COIL
rNow = Math.sqrt(rpVect[0]*rpVect[0]+rpVect[1]*rpVect[1]);
- r=coilRadius;
- xNSteps=((r-rNow)>1.)? ((r-rNow)*10.):10.; //C.M. step 1/10 of the range
+ r = coilSlices.get(0).getRmin();
+ xNSteps=((r-rNow)>1.)? ((r-rNow)*10.):10.; //C.M. step 1/10 of the range
pos = new double[]{rpVect[0], rpVect[1], rpVect[2]};
field = det.getFieldMap().getField(pos); // C.M.-June-03
@@ -108,11 +123,38 @@
}
// In the COIL
-
- r = r + 50.; // Coil Thickness 50cm
+ int islice = 0;
+
+ for(int ii=0; ii<coilSlices.size(); ++ii) {
+ CoilSubLayer slice = coilSlices.get(ii);
+ if(stpr.getStopTkELow()) break; // C.M.
+ double zmax = slice.getZmax();
+ double rmin = slice.getRmin();
+ dr = slice.getThickness();
+ double totalDedx = slice.getdEdx() * dr;
+
+ if( Math.abs(rpVect[2])>zmax ) break;
+ // Go on to the next layer
+ r = r + dr;
+ xNSteps=(r-Math.sqrt(rpVect[0]*rpVect[0]+rpVect[1]*rpVect[1]))*10.;
+
+ field = det.getFieldMap().getField( rpVect ); // C.M.-June-03
+ double radius = Math.sqrt(rpVect[0]*rpVect[0]+rpVect[1]*rpVect[1]);
+// System.out.println("slice "+islice+", mat="+slice.getMaterial()
+// +", r="+r
+// +", zmax="+zmax+", rmin="+rmin
+// +", dr="+dr
+// +", totalDedx="+totalDedx
+// +", field="+field[2]
+// +", rad="+radius+", z="+rpVect[2] );
+
+ stepConditions = steprConditions(xNSteps, field[2], totalDedx);
+ stpr.tkSteps(r,zmax,stepConditions);
+ rpVect = stpr.getNewRp();
+ ++islice;
+ }
- double[] posmu = new double[]{rpVect[0], rpVect[1], rpVect[2]};
- field = det.getFieldMap().getField(posmu); // C.M.-June-03
+ field = det.getFieldMap().getField(rpVect); // C.M.-June-03
System.out.println("MuonCalSegFinder: muonRadius="+rmin+", field="+muField);
muField = field[2];
@@ -145,7 +187,7 @@
rNow=Math.sqrt(rpVect[0]*rpVect[0]+rpVect[1]*rpVect[1]);
xNSteps = (r-rNow)*10.;
double muFieldNorm = 0.;
- stepConditions = steprConditions(xNSteps,muFieldNorm,meanDEdxPerLay(aname));
+ stepConditions = steprConditions(xNSteps,muFieldNorm,meanDEdxPerLay(subdetName));
stpr.tkSteps(r,zmax,stepConditions);
rpVect = stpr.getNewRp();
// Get intersections with succesive layers as long as track is
@@ -163,6 +205,7 @@
encoder.setValue("theta", itheta);
encoder.setValue("phi", iphi);
Long cellid = encoder.getID();
+ System.out.println("adding to cellList: cellid="+MyTools.printID(cellid));
cellList.add(cellid);
if(stpr.getStopTkELow())break; // C.M.
@@ -170,7 +213,7 @@
if(nr!=(activLayers-1)) r = r + dr;
else {r=r+0.5*dr;}
xNSteps=(r-Math.sqrt(rpVect[0]*rpVect[0]+rpVect[1]*rpVect[1]))*10.;
- stepConditions = steprConditions(xNSteps,muFieldNorm,meanDEdxPerLay(aname));
+ stepConditions = steprConditions(xNSteps,muFieldNorm,meanDEdxPerLay(subdetName));
stpr.tkSteps(r,zmax,stepConditions);
rpVect = stpr.getNewRp();
nr++;
@@ -184,13 +227,14 @@
// HistogramFolder.setDefaultFolder("/MUDet"); // C.M.-27Jan03- Folder
-
- matchHitsFast(segm, dataMgr.getCollHitMap(muonName) );
+ Map<Long,SimCalorimeterHit> muhits = dataMgr.getCollHitMap(muonName);
+ System.out.println("# muon hits: "+muhits.size());
+ matchHitsFast(segm, dataMgr.getCollHitMap(muonName), true );
// HistogramFolder.setDefaultFolder(".."); // C.M.-27Jan03- Folder
}
private double coilRadius;
private static double materDedxMudet; // C.M. get the dEdx
private String muonName = "MuonBarrHits";
+ private Vector<CoilSubLayer> coilSlices;
}
-
