lcsim/src/org/lcsim/recon/tracking/seedtracker

diff -u -r1.9 -r1.10
--- MaterialManager.java	3 Mar 2011 21:37:58 -0000	1.9
+++ MaterialManager.java	11 Apr 2011 17:21:07 -0000	1.10
@@ -38,7 +38,7 @@

  * This class should now group together SiTrackerEndcap2 layers correctly.
  * 
  * @author Jeremy McCormick

- * @version $Id: MaterialManager.java,v 1.9 2011/03/03 21:37:58 jeremy Exp $

+ * @version $Id: MaterialManager.java,v 1.10 2011/04/11 17:21:07 mgraham Exp $

  */
 public class MaterialManager
 {

@@ -49,6 +49,7 @@

     private List<MaterialPolyconeSegment> _matpc = new ArrayList<MaterialPolyconeSegment>();
     private List<MaterialCylinder> _matcyl = new ArrayList<MaterialCylinder>();
     private List<MaterialDisk> _matdsk = new ArrayList<MaterialDisk>();

+    private List<MaterialXPlane> _matxpl = new ArrayList<MaterialXPlane>();

     private HashMap<ISolid, Double> solid_vol_map = new HashMap<ISolid, Double>(400);
     private static double _rmax;
     private static double _zmax = 1800.;

@@ -135,7 +136,40 @@

             return pathGroups;
         }
     }

-    

+
+        /**
+     * Get the path groups for SiTrackerFixedTarget2
+     */
+    static private class SiTrackerFixedTarget2VolumeGrouper implements SubdetectorVolumeGrouper
+    {
+        public List<List<String>> getPathGroups(Subdetector subdet, IPhysicalVolume topVol)
+        {
+            List<List<String>> pathGroups = new ArrayList<List<String>>();
+            // Positive and negative endcap loop.
+            for (IDetectorElement endcaps : subdet.getDetectorElement().getChildren())
+            {
+                // Layer loop.
+                for (IDetectorElement layer : endcaps.getChildren())
+                {
+                    List<String> modulePaths = new ArrayList<String>();
+                    // System.out.println(layer.getName());
+
+                    // Module loop.
+                    for (IDetectorElement module : layer.getChildren())
+                    {
+                        String path = "";
+                        PhysicalVolumeNavigator.getLeafPaths(modulePaths, module.getGeometry().getPhysicalVolume(),
+                                path);
+                    }
+                    // Add module paths to this layer.
+                    pathGroups.add(modulePaths);
+                }
+            }
+            return pathGroups;
+        }
+    }
+
+

     /**
      * Default VolumeGroup for endcaps with physical layers.
      */

@@ -182,6 +216,8 @@

         SubdetectorVolumeGrouper endcapGrouper2 = new SiTrackerEndap2VolumeGrouper();
         subdetGroups.put(SiTrackerEndcap2.class, endcapGrouper2);
         subdetGroups.put(SiTrackerFixedTarget2.class, endcapGrouper2);

+        SubdetectorVolumeGrouper fixedtargetGrouper = new SiTrackerFixedTarget2VolumeGrouper();
+        subdetGroups.put(SiTrackerFixedTarget2.class, fixedtargetGrouper);

 
         // Add VolumeGrouper for PolyconeSupport.
         subdetGroups.put(PolyconeSupport.class, new PolyconeSupportVolumeGrouper());        

@@ -285,11 +321,16 @@

                             System.out.println("VolumeGroupInfo ...");
                             System.out.println("    rmax = " + vgi.rmax);
                             System.out.println("    rmin = " + vgi.rmin);

+                            System.out.println("    xmin = " + vgi.xmin);
+                            System.out.println("    xmax = " + vgi.xmax);
+                            System.out.println("    ymin = " + vgi.ymin);
+                            System.out.println("    ymax = " + vgi.ymax);

                             System.out.println("    zmin = " + vgi.zmin);
                             System.out.println("    zmax = " + vgi.zmax);
                             System.out.println("    X0 = " + vgi.X0);
                             System.out.println("    weighted_r = " + vgi.weighted_r);
                             System.out.println("    weighted_z = " + vgi.weighted_z);

+                            System.out.println("    weighted_z = " + vgi.weighted_y);

                             System.out.println("    vtot_tube_only = " + vgi.vtot_tube_only);
                             System.out.println("    vtot = " + vgi.vtot);
                         }

@@ -335,14 +376,32 @@

             // Calculate the average radiation length for this volume
 
             // Determine if this volume should be modeled as barrel or disk

-            if (isCylinder(vgi.rmin, vgi.rmax, vgi.zmin, vgi.zmax))

+            if (isXPlane(vgi.xmin, vgi.xmax))
+            {
+                // Calculate the weighted radius of the elements
+                double zlen = vgi.zmax - vgi.zmin;
+                double ylen = vgi.ymax - vgi.ymin;
+                double thickness = vtot / (ylen * zlen * vgi.X0);
+                double x=(vgi.xmax-vgi.xmin)/2;
+
+                if (DEBUG) {
+                    System.out.println("Treating as a XPlane...x0: " + vgi.X0 + "| zmin: " + vgi.zmin +
+                            "| zmax: " + vgi.zmax + "| vtot: " + vtot +
+                            "| thickness: " + thickness + "| rmin: " + vgi.rmin +
+                            "| rmax: " + vgi.rmax);
+                    System.out.println();
+                }
+
+                _matxpl.add(new MaterialXPlane(vgi.ymin,vgi.ymax, vgi.zmin, vgi.zmax, x, thickness));
+            }
+            else if(isCylinder(vgi.rmin, vgi.rmax, vgi.zmin, vgi.zmax))

             {
                 // Calculate the weighted radius of the elements
                 double zlen = vgi.zmax - vgi.zmin;
                 double thickness = vtot / (2. * Math.PI * vgi.weighted_r * zlen * vgi.X0);
 
                 if (DEBUG) {

-                    System.out.println("x0: " + vgi.X0 + "| zmin: " + vgi.zmin +

+                    System.out.println("Treating as a Cylinder...x0: " + vgi.X0 + "| zmin: " + vgi.zmin +

                             "| zmax: " + vgi.zmax + "| vtot: " + vtot +
                             "| thickness: " + thickness + "| rmin: " + vgi.rmin +
                             "| rmax: " + vgi.rmax);

@@ -397,6 +456,10 @@

     {
         return _matdsk;
     }

+    public List<MaterialXPlane> getMaterialXPlanes()
+    {
+        return _matxpl;
+    }

 
     public List<MaterialPolyconeSegment> getMaterialPolyconeSegments()
     {

@@ -418,6 +481,11 @@

         return (rmax - rmin) * Math.abs(zmax + zmin) < (zmax - zmin) * (rmax + rmin);
     }
 

+    private boolean isXPlane(double xmin, double xmax)
+    {
+        return (xmax-xmin)<1.0;//1mm
+    }
+

     // special handling for Polycone...
     private void handlePolycone(IPhysicalVolume pv)
     {

@@ -474,6 +542,11 @@

         double rmin = 1.e10;
         double zmin = 1.e10;
         double zmax = -1.e10;

+        //mg 3/14/11   MaterialXPlane info
+        double ymin = 1.e10;
+        double ymax=-1.e10;
+        double xmin = 1.e10;
+        double xmax=-1.e10;

     }
 
     /**

@@ -490,6 +563,13 @@

         double weighted_z = 0.0;
         double vtot_tube_only = 0.;
         double vtot = 0.0;

+        //mg 3/14/11   MaterialXPlane info
+        double ymin = 1.e10;
+        double ymax=-1.e10;
+        double xmin = 1.e10;
+        double xmax=-1.e10;
+        double weighted_y = 0.0;
+

     }
 
     // This function performs all the calculations on lists of physical volumes

@@ -536,7 +616,13 @@

             vgi.zmax = Math.max(vi.zmax, vgi.zmax);
             vgi.rmin = Math.min(vi.rmin, vgi.rmin);
             vgi.rmax = Math.max(vi.rmax, vgi.rmax);

-

+            //mg 3/14/11  also store y,x information
+            vgi.ymin = Math.min(vi.ymin, vgi.ymin);
+            vgi.ymax = Math.max(vi.ymax, vgi.ymax);
+            double y0 = pv.getLtoGTransform().getTranslation().y();
+            vgi.weighted_y += y0 * wgt;
+            vgi.xmin = Math.min(vi.xmin, vgi.xmin);
+            vgi.xmax = Math.max(vi.xmax, vgi.xmax);

         }
 
         // finish weighted R/Z calculations + perform X0 calculation

@@ -545,6 +631,8 @@

             vgi.weighted_r /= totwgt;
             vgi.weighted_z /= totwgt;
             vgi.X0 = vgi.vtot / totwgt;

+            //mg 3/14/11  also y info
+            vgi.weighted_y/=totwgt;

         }
 
         return vgi;

@@ -615,6 +703,11 @@

                 double r = Math.sqrt(transformed.x() * transformed.x() + transformed.y() * transformed.y());
                 vi.rmin = Math.min(vi.rmin, r);
                 vi.rmax = Math.max(vi.rmax, r);

+                //mg 3/14/11  also store ymin,ymax
+                vi.ymin = Math.min(transformed.y(), vi.ymin);
+                vi.ymax = Math.max(transformed.y(), vi.ymax);
+                vi.xmin = Math.min(transformed.x(), vi.xmin);
+                vi.xmax = Math.max(transformed.x(), vi.xmax);

             }
         } // Note: this information will NOT be used most of the time...
         // Polycones that are top-level elements (e.g. the beampipe) are

lcsim/src/org/lcsim/recon/tracking/seedtracker

diff -u -r1.1 -r1.2
--- MaterialXPlane.java	14 Mar 2011 20:10:17 -0000	1.1
+++ MaterialXPlane.java	11 Apr 2011 17:21:07 -0000	1.2
@@ -16,14 +16,16 @@

     private double _ymax;
     private double _zmin;
     private double _zmax;

+    private double _x;

     private double _t_RL;
     
     /** Creates a new instance of MaterialCylinder */

-    public MaterialXPlane(double ymin, double ymax, double zmin, double zmax, double t_RL) {

+    public MaterialXPlane(double ymin, double ymax, double zmin, double zmax, double x, double t_RL) {

         _ymin = ymin;
         _ymax = ymax;
         _zmin = zmin;
         _zmax = zmax;

+        _x = x;

         _t_RL = t_RL;
     }
     

@@ -43,6 +45,10 @@

         return _zmax;
     }
     

+    public double x() {
+        return _x;
+    }
+

     public double ThicknessInRL() {
         return _t_RL;
     }

lcsim/src/org/lcsim/recon/tracking/seedtracker

diff -u -r1.5 -r1.6
--- MultipleScattering.java	23 Nov 2010 01:30:11 -0000	1.5
+++ MultipleScattering.java	11 Apr 2011 17:21:07 -0000	1.6
@@ -57,6 +57,7 @@

         //  Retrieve the cylinder and disk material models from the material manager
         List<MaterialCylinder> matcyl = _materialmanager.getMaterialCylinders();
         List<MaterialDisk> matdsk = _materialmanager.getMaterialDisks();

+        List<MaterialXPlane> matxpl = _materialmanager.getMaterialXPlanes();

         
         //  Find the largest path length to a hit
         double smax = 9999.;

@@ -105,7 +106,31 @@

                 }
             }
         }

-       

+        //mg 3/14/11  add in XPlanes
+        for (MaterialXPlane xpl : matxpl) {
+               double scmin = HelixUtils.PathToZPlane(helix, xpl.zmin());
+            double scmax = HelixUtils.PathToZPlane(helix, xpl.zmax());
+              if (scmin > scmax) {
+                double temp = scmin;
+                scmin = scmax;
+                scmax = temp;
+            }
+            List<Double> pathlist = HelixUtils.PathToXPlane(helix, xpl.x(), smax, _mxint);
+            for (Double s : pathlist) {
+                Hep3Vector dir = HelixUtils.Direction(helix, s);
+                Hep3Vector pos = HelixUtils.PointOnHelix(helix, s);
+                double y=pos.y();
+                double z=pos.z();
+                if (y >= xpl.ymin() && y <= xpl.ymax()&&z >= xpl.zmin() && z <= xpl.zmax()) {
+                    Hep3Vector xhat = VecOp.unit(new BasicHep3Vector(pos.x(),0. , 0.));
+                    double cth = Math.abs(VecOp.dot(dir, xhat));
+                    double radlen = xpl.ThicknessInRL() / Math.max(cth, .001);
+                    double angle = msangle(helix.p(_bfield), radlen);
+//                    System.out.println("MultipleScatter:  angle="+angle);
+                    scatters.add(new ScatterAngle(s, angle));
+                }
+            }
+        }

         //  Sort the multiple scatters by their path length
         Collections.sort(scatters);
         return scatters;

lcsim/src/org/lcsim/fit/helicaltrack

diff -u -r1.12 -r1.13
--- HelixUtils.java	1 Feb 2011 22:45:05 -0000	1.12
+++ HelixUtils.java	11 Apr 2011 17:21:07 -0000	1.13
@@ -120,7 +120,43 @@

         if (Math.abs(safeslope) < _minslope) safeslope = _minslope * Math.signum(safeslope);
         return zdist / safeslope;
     }

-    

+
+        /**
+     * Return the x-y path length to an x-plane.
+     * @param helix HelicalTrackFit to be used in calculating the path length
+     * @param x location of x-plane
+     * @return path length from the DCA to the x-plane
+     */
+    public static List<Double>  PathToXPlane(HelicalTrackFit helix, double x,double smax, int mxint) {
+        //  Create a list to hold the path lengths
+        List<Double> pathlist = new ArrayList<Double>();
+        //  Retrieve helix dca and RC
+        double x0 = helix.x0();
+        double y0 = helix.y0();
+
+        double xc=helix.xc();
+        double yc=helix.yc();
+        double RC = helix.R();
+        double safeslope = helix.slope();
+        if (Math.abs(safeslope) < _minslope) safeslope = _minslope * Math.signum(safeslope);
+        double y=yc+Math.sqrt(RC*RC-Math.pow(x-xc,2));
+        double dphi=Math.atan2(y-y0, x-x0);
+        Double s = dphi * Math.abs(RC);
+        //  Loop over crossings until we exceed one of the limits
+        while (s < smax && pathlist.size() < mxint) {
+            //  Add this odd-numbered crossing to the list
+            pathlist.add(s);
+                //  Advance to the next even-numbered crossing
+            s += 2. * (Math.PI - dphi) * Math.abs(RC);
+            //  Check to see if we should add it
+            if (s < smax && pathlist.size() < mxint) pathlist.add(s);
+            //  Add this even-numbered crossing to the list
+            s += 2. * dphi * Math.abs(RC);
+        }
+        return pathlist;
+    }
+
+

     /**
      * Return a list of x-y path lengths to a cylinder centered on the z axis.
      * @param helix HelicalTrackFit to be used in calculating the path length