LISTSERV 16.5 - LCD-CVS Archives

+/* + * BarrelCylinderSegmentationBase.java + * + * 20051124 - G.Lima - Created + */ + +package org.lcsim.geometry.segmentation; + +import org.jdom.Element; +import org.lcsim.geometry.layer.Layer; +import org.lcsim.geometry.layer.LayerStack; + +/** + * Base implementation of a barrel segmentation, it provides common + * code implementation. + * + * @author Guilherme Lima + * @versoin $Id: BarrelCylinderSegmentationBase.java,v 1.1 2005/12/07 10:39:37 lima Exp $ + */ +public abstract class BarrelCylinderSegmentationBase extends SegmentationBase +{ + public BarrelCylinderSegmentationBase(Element e) + { + super(e); + } + + /** + * Return the layer number based on the radial distance to cylinder axis. + * + * @param r radial distance to cylinder axis + * @return layer number of layer corresponding to that distance + * (may be either in absorber or live material) + */ + public int getLayerBin(double r) { + // In order to be general, we should not assume that all + // layers have the same thickness. Therefore, one has to + // guess the starting layer (based on average thickness), and + // then navigate through layers until one finds the right one + double depth = r - getRMin(); + double mean_t = (getRMax()-getRMin()) / getNumberOfLayers(); + + int ilay = (int)Math.floor( depth / mean_t ); + LayerStack stack = getLayering().getLayers(); + Layer layer = stack.getLayer(ilay); + double depHi = stack.getThicknessToLayerBack(ilay); + double depLo = depHi - layer.getThickness(); + for(;;) { + if( depth>depLo && depth<=depHi ) return ilay; + if( depth<=depLo ) { + --ilay; + depHi = depLo; + layer = stack.getLayer(ilay); + depLo -= layer.getThickness(); + } + if( depth>depHi ) { + ++ilay; + depLo = depHi; + layer = stack.getLayer(ilay); + depHi += layer.getThickness(); + } + } + } +}

GeomConverter/src/org/lcsim/geometry/segmentation

GridXYZ.java 1.16 -> 1.17

diff -u -r1.16 -r1.17
--- GridXYZ.java	28 Oct 2005 00:34:23 -0000	1.16
+++ GridXYZ.java	7 Dec 2005 10:39:37 -0000	1.17
@@ -14,6 +14,10 @@

 import org.jdom.Element;
 
 import org.lcsim.geometry.util.IDDescriptor;

+import org.lcsim.geometry.util.IDEncoder;
+import org.lcsim.geometry.layer.Layer;
+import org.lcsim.geometry.layer.LayerStack;
+import org.lcsim.geometry.subdetector.CylindricalCalorimeter;

 
 /**
  * @author jeremym

@@ -27,104 +31,228 @@

     private double gridSizeX = 0;
     private double gridSizeY = 0;
     private double gridSizeZ = 0;

     private int xIndex = -1;
     private int yIndex = -1;
     private int zIndex = -1;

     private double[] _localPos = { 0, 0, 0};
     private double[] _globalPos = { 0, 0, 0};

-    
-    public double getGridSizeX()
-    {
-        return gridSizeX;
-    }
-    
-    public double getGridSizeY()
-    {
-        return gridSizeY;
-    }
-    
-    public double getGridSizeZ()
-    {
-        return gridSizeZ;
-    }
-

     /** Creates a new instance of GridXYZ */
     public GridXYZ(Element node) throws DataConversionException
     {
         super(node);

         if ( node.getAttribute("gridSizeX") != null )
         {
             gridSizeX = node.getAttribute("gridSizeX").getDoubleValue();
         }

         if ( node.getAttribute("gridSizeY") != null )
         {
             gridSizeY = node.getAttribute("gridSizeY").getDoubleValue();
         }

         if ( node.getAttribute("gridSizeZ") != null )
         {
             gridSizeZ = node.getAttribute("gridSizeZ").getDoubleValue();
         }
     }

-    
-    /* FIXME: Implementation requires very specific information about the subdetector topology. --JM */
-    public long[] getNeighbourIDs(int a, int b, int c)

+
+    /** Find neighbouring cells to the current cell.
+     *
+     * Cell neighbors are found based on the direction (theta,phi) of
+     * the reference cell w.r.t the origin.
+     *
+     * @see org.lcsim.geometry.segmentation.SegmentationBase#setID()
+     * @return array of cellIDs for the neighbouring cells
+     */
+    public long[] getNeighbourIDs(int layerRange, int xRange, int yRange)

-        return new long[3];

+//       System.out.println("Nonproj neighbs: "+layerRange+" "+xRange+" "+yRange);
+      IDEncoder encoder = new IDEncoder( descriptor );
+      encoder.setValues(values);
+      long saveID = encoder.getID();
+
+      int klay = values[layerIndex];
+      int kx   = values[xIndex];
+      int ky   = values[yIndex];
+//       System.out.println("Ref.ID: ref="+klay+" "+kx+" "+ky
+// 			 +" (hex "+Long.toHexString(saveID)+")");
+//       long system = (saveID>>7) & 0x1ff;
+//       System.out.println("hit pos: x="+getX()+", y="+getY()+", z="+getZ()+", system="+system);
+
+      int nMax = (2*layerRange + 1)*(2*xRange + 1)*(2*yRange + 1) - 1;
+      long[] result = new long[nMax];
+
+      // theta and phi are used to find central neighbors in other layers
+      double theta = getTheta();
+      double phi = getPhi();
+
+      // this is not general, but specific to endcaps
+      double dx = gridSizeX;
+      double dy = gridSizeY;
+      int zBins = detector.getLayering().getLayerCount();
+      double Rmax = ((CylindricalCalorimeter)detector).getOuterRadius();
+
+      int size = 0;
+      int refLayer = values[layerIndex];
+      boolean isNegativeZ = getZ()<0 ? true : false;
+      for (int i=-layerRange; i<=layerRange; ++i) {
+	int ilay = refLayer + i;
+
+	if (ilay<0 || ilay>=zBins ) continue;
+	encoder.setValue(layerIndex,ilay);
+
+	int xBins = (int)Math.ceil( 2*Rmax / dx );
+	int yBins = xBins;
+
+	// project theta,phi to sensitive middle of current layer
+	double z = getDepthSensitiveMid( ilay );
+	if( isNegativeZ ) z = -z;
+	double rho = z * Math.tan(theta);
+	double x = rho * Math.cos(phi);
+	double y = rho * Math.sin(phi);
+//  	System.out.println("layer "+ilay+", cylR="+rho
+//  			   +", theta="+(theta*180./Math.PI)
+// 			   +", phi="+(phi*180./Math.PI)
+//  			   +", rvec=( "+x+", "+y+", "+z+")");
+
+ 	long id = this.findCellContainingXYZ(x,y,z);
+ 	if(id==0) continue;
+
+	this.setID(id);
+// 	System.out.println("from findXYZ: x="+getX()+", y="+getY()+", z="+getZ());
+ 	for (int j=-xRange; j<=xRange; ++j) {
+	  int ix = values[xIndex] + j;
+
+ 	  if (ix<-xBins/2 || ix>=xBins/2) continue;
+ 	  encoder.setValue(xIndex,ix);
+
+ 	  for (int k=-yRange; k<=yRange; ++k) {
+	    // skip reference cell
+	    if (i==0 && j==0 && k==0) continue;
+
+ 	    int iy = values[yIndex] + k;
+	    if (iy<-yBins/2 || iy>=yBins/2) continue;
+
+//    	    System.out.println("Adding neighbor: "+ilay+" "+ix+" "+iy
+//  			       +", total="+ (size+1));
+ 	    result[size++] = encoder.setValue(yIndex,iy);
+	  }
+ 	}
+      }
+
+      // reset encoder and decoder
+      this.setID(saveID);
+      if(size < result.length) {
+// 	  int kklay = values[layerIndex];
+// 	  int kkx = values[xIndex];
+// 	  int kky = values[yIndex];
+//  	  System.out.println("In a border? Neighbors="+size+"/"+result.length
+//  			     +", ref="+kklay+" "+kkx+" "+kky);
+
+	  long[] temp = new long[size];
+	  System.arraycopy(result,0,temp,0,size);
+	  result = temp;
+      }
+
+//       int kklay = values[layerIndex];
+//       int kkx = values[xIndex];
+//       int kky = values[yIndex];
+//       System.out.println("Neighbors="+size+"/"+result.length
+// 			 +", ref="+kklay+" "+kkx+" "+kky
+// 			 +" (hex "+Long.toHexString(encoder.getID())+")");
+
+//       System.out.println("Leaving getNeighbours: current id="
+// 			 +values[layerIndex]+" "+values[xIndex]
+// 			 +" "+values[yIndex]);
+      return result;
+    }
+
+    // Not for public use, this is needed to calculate positions
+    // and number of cells, etc.
+    public double getZMin() {
+	return ((CylindricalCalorimeter)detector).getZMin();

+
+    public double getZMax() {
+	return ((CylindricalCalorimeter)detector).getZMax();
+    }
+
+    public double getRMin() {
+	return ((CylindricalCalorimeter)detector).getInnerRadius();
+    }
+
+    public double getRMax() {
+	return ((CylindricalCalorimeter)detector).getOuterRadius();
+    }
+
+    /** FIXME: Doesn't belong here, as it is always generally derivable from x and y. */

     public double getPhi()
     {

-        return SegmentationUtil.getPhi(getX(), getY() );

+        double phi = atan2(getY(), getX() );
+
+        if ( phi < 0 )
+        {
+            phi += 2 * PI;
+        }
+
+        return phi;

+
+    /** FIXME: Doesn't belong here, as it is always computable given x and y. */

     public double getTheta()
     {

-        return SegmentationUtil.getTheta(getX(), getY(), getZ() );

+        double theta = atan(getCylindricalRadiusFromPosition() / getZ() );
+
+        if ( theta < 0 )
+        {
+            theta += PI;
+        }
+
+        return theta;

     public double getX()
     {
         return getPosition()[0];
     }

     public double getY()
     {
         return getPosition()[1];
     }

     public double getZ()
     {
         return getPosition()[2];
     }

     public double[] getPosition()
     {
         return _globalPos;
     }

     private void computePosition()
     {
         computeLocalPosition();
         computeGlobalPosition();
     }

     private void computeLocalPosition()
     {
         computeLocalX();
         computeLocalY();
         computeLocalZ();
     }

     private void computeGlobalPosition()
     {
         _globalPos = transformLocalToGlobal(_localPos);
     }

     private void computeLocalX()
     {
         if ( xIndex != -1 )

@@ -132,7 +260,7 @@

             _localPos[0] = (((double)getValue(xIndex)) + 0.5) * gridSizeX;
         }
     }

     private void computeLocalY()
     {
         if ( yIndex != -1 )

@@ -140,7 +268,7 @@

             _localPos[1] = (((double)getValue(yIndex)) + 0.5) * gridSizeY;
         }
     }

     private void computeLocalZ()
     {
         if ( zIndex != -1 )

@@ -148,19 +276,19 @@

             _localPos[2] = (((double)getValue(zIndex)) + 0.5) * gridSizeZ;
         }
     }

     /** Overridden to cache the global position. */
     public void setID(long id)
     {
         super.setID(id);
         computePosition();
     }

     public double getCylindricalRadiusFromPosition()
     {
         return sqrt(getX() * getX() + getY() * getY() );
     }

     /**
      * Check for a valid grid size for each of XYZ.
      * If the grid is not valid, then leave those

@@ -172,20 +300,120 @@

     public void setIDDescription(IDDescriptor id)
     {
         super.setIDDescription(id);

         if ( gridSizeX != 0 )
         {
             xIndex = id.indexOf("x");
         }

         if ( gridSizeY != 0 )
         {
             yIndex = id.indexOf("y");
         }

         if ( gridSizeZ != 0 )
         {
             zIndex = id.indexOf("z");
         }
     }

-}

\ No newline at end of file

+
+    /**
+     * Returns positive distance from IR to center of sensitive slice
+     * of any layer
+     * @param layer layer index
+     */
+    private double getDepthSensitiveMid(int ilay) {
+	LayerStack stack = detector.getLayering().getLayers();
+	Layer layer = stack.getLayer(ilay);
+
+	double preLayers = 0;
+	if(ilay>0) preLayers = stack.getSectionThickness(0,ilay-1);
+
+	return this.getZMin() + preLayers + layer.getThicknessToSensitiveMid();
+    }
+
+    /**
+     * Return the cell which contains a given point (x,y,z).  If point is
+     * not contained in this component, zero is returned.
+     *
+     * @param x,y,z cartesian coordinates of the point
+     * @return ID of cell containing the point (maybe either in
+     * absorber or live material)
+     */
+    public long findCellContainingXYZ(double x, double y, double z) {
+
+	// validate point
+	double absz = Math.abs(z);
+	if(absz<getZMin()) return 0;
+	if(absz>getZMax()) return 0;
+	double rho = Math.sqrt(x*x+y*y);
+	if(rho<getRMin()) return 0;
+	if(rho>getRMax()) return 0;
+
+	// ok, point is valid, so a valid ID should be returned
+	int ix = getXBin(x);
+	int iy = getYBin(y);
+	int ilay = getLayerBin(absz);
+
+	IDEncoder enc = new IDEncoder(descriptor);
+	enc.setValues(values);
+	enc.setValue(layerIndex, ilay);
+	enc.setValue(xIndex, ix);
+	enc.setValue(yIndex, iy);
+	long resultID = enc.getID();
+	return resultID;
+    }
+
+    public int getLayerBin(double z) {
+      // In order to be general, we should not assume that all
+      // layers have the same thickness.  Therefore, one has to
+      // guess the starting layer (based on average thickness), and
+      // then navigate through layers until one finds the right one
+      double depth = z - getZMin();
+      double mean_t = (getZMax()-getZMin()) / getNumberOfLayers();
+
+      int ilay = (int)Math.floor( depth / mean_t );
+      LayerStack stack = getLayering().getLayers();
+      Layer layer = stack.getLayer(ilay);
+      double depHi = stack.getThicknessToLayerBack(ilay);
+      double depLo = depHi - layer.getThickness();
+      for(;;) {
+	if( depth>depLo && depth<=depHi ) return ilay;
+	if( depth<=depLo ) {
+	  --ilay;
+	  depHi = depLo;
+	  layer = stack.getLayer(ilay);
+	  depLo -= layer.getThickness();
+	}
+	if( depth>depHi ) {
+	  ++ilay;
+	  depLo = depHi;
+	  layer = stack.getLayer(ilay);
+	  depHi += layer.getThickness();
+	}
+      }
+    }
+
+    public int getXBin(double x) {
+	return getBin( x, gridSizeX );
+    }
+    public int getYBin(double y) {
+	return getBin( y, gridSizeY );
+    }
+
+    public int getBin(double u, double gridSizeU) {
+      int numBins = (int)Math.floor( 2*getRMax() / gridSizeU );
+      double u0 = gridSizeU / 2;
+      int iu = (int)Math.floor((u-u0)/gridSizeU + 0.5);
+      return iu;
+    }
+
+    public double getGridSizeX() { return gridSizeX; }
+    public double getGridSizeY() { return gridSizeY; }
+    public double getGridSizeZ() { return gridSizeZ; }
+
+    private int flag1 = 0;
+    private int flag2 = 0;
+    private int flag3 = 0;
+    private int flag4 = 0;
+}

- - return result; - } - - // Not for public use, this is needed to calculate positions - // and number of cells, etc. - /* - * FIXME: Are these next four functions really needed? --JM - */ - private double getZMin() { - return ((CylindricalSubdetector)detector).getZMin(); - }

- private double getZMax() { - return ((CylindricalSubdetector)detector).getZMax(); - } - - private double getRMin() { - return ((CylindricalSubdetector)detector).getInnerRadius(); - } - - private double getRMax() { - return ((CylindricalSubdetector)detector).getOuterRadius();

- if(rho<getRMin()) { - System.out.println("find: fail rho<Rmin: rho="+rho); - return 0; - } - if(rho>getRMax()) { - System.out.println("find: fail rho>Rmax: rho="+rho); - return 0; - }

- public int getLayerBin(double r) { - // In order to be general, we should not assume that all - // layers have the same thickness. Therefore, one has to - // guess the starting layer (based on average thickness), and - // then navigate trhough layers until one finds the right one - double depth = r - getRMin(); - double mean_t = (getRMax()-getRMin()) / getNumberOfLayers(); - - int ilay = (int)Math.floor( depth / mean_t ); - LayerStack stack = getLayering().getLayers(); - Layer layer = stack.getLayer(ilay); - double depHi = stack.getThicknessToLayerBack(ilay); - double depLo = depHi - layer.getThickness(); - for(;;) { - if( depth>depLo && depth<=depHi ) return ilay; - if( depth<=depLo ) { - --ilay; - depHi = depLo; - layer = stack.getLayer(ilay); - depLo -= layer.getThickness(); - } - if( depth>depHi ) { - ++ilay; - depLo = depHi; - layer = stack.getLayer(ilay); - depHi += layer.getThickness(); - } - } - } -

ProjectiveZPlane(Element node) throws DataConversionException { super(node); thetaBins = node.getAttribute("thetaBins").getIntValue(); phiBins = node.getAttribute("phiBins").getIntValue(); }

int nMax = (2*deltaLayer + 1)*(2*deltaTheta + 1)*(2*deltaPhi + 1) - 1; int size = 0; long[] result = new long[nMax]; for (int i=-deltaLayer;i<=deltaLayer;i++) { int l = values[layerIndex] + i;

+ + /** + * Return the cell which contains a given point (x,y,z), or zero. + * + * @param x,y,z cartesian coordinates of the point + * @return ID of cell containing the point (maybe either in + * absorber or live material) + */ + public long findCellContainingXYZ(double x, double y, double z) { + + // validate point + if(z<getZMin()) return 0; + if(z>getZMax()) return 0; + double rho = Math.sqrt(x*x+y*y); + if(rho<getRMin()) return 0; + if(rho>getRMax()) return 0; + + // ok, point is valid, so a valid ID should be returned + int ilay = getLayerBin(rho); + + double phi = Math.atan2( y, x ); + if(phi<0) phi += 2*Math.PI; + int iphi = (int)( phi / phiBins ); + + double theta = Math.atan2( rho, z ); + if(theta<0) { + // If never prints out, the whole if can be dropped + System.out.println("ProjCylinder: Is this really needed?!?"); + theta += 2.*Math.PI; + } + int itheta = (int)(theta/thetaBins); + + IDEncoder enc = new IDEncoder(descriptor); + enc.setValue("layer",ilay); + enc.setValue(thetaIndex,itheta); + enc.setValue(phiIndex,iphi); + long resultID = enc.getID(); + + return resultID; + } + + /** + * Return the layer number based on the z-coordinate + * + * @param z z-coordinate + * @return layer number of layer corresponding to that distance + * (may be either in absorber or live material) + * @throws RuntimeException if abs(z)<zMin + */ + public int getLayerBin(double z) { + // In order to be general, we should not assume that all + // layers have the same thickness. Therefore, one has to + // guess the starting layer (based on average thickness), and + // then navigate through layers until one finds the right one + double depth = Math.abs(z) - getZMin(); + if(depth<0) throw new RuntimeException("ProjectiveZPlane: Error: z < zMin, z="+z+", zMin="+this.getZMin()); + + double mean_t = (getZMax()-getZMin()) / getNumberOfLayers(); + + int ilay = (int)Math.floor( depth / mean_t ); + LayerStack stack = getLayering().getLayers(); + Layer layer = stack.getLayer(ilay); + double depHi = stack.getThicknessToLayerBack(ilay); + double depLo = depHi - layer.getThickness(); + for(;;) { + if( depth>depLo && depth<=depHi ) return ilay; + if( depth<=depLo ) { + --ilay; + depHi = depLo; + layer = stack.getLayer(ilay); + depLo -= layer.getThickness(); + } + if( depth>depHi ) { + ++ilay; + depLo = depHi; + layer = stack.getLayer(ilay); + depHi += layer.getThickness(); + } + } + }

+ + /** + * Returns the cell which contains a given point (x,y,z). + * + * @param x,y,z cartesian coordinates of the point + * @return ID of cell containing the point (maybe either in + * absorber or live material), or zero if the point is + * not inside this component. + */ + abstract public long findCellContainingXYZ(double x, double y, double z); + + /** + * Returns the cell which contains a given point (x,y,z). + * + * @param pos 3-dim array with cartesian coordinates of the point + * @return ID of cell containing the point (maybe either in + * absorber or live material), or zero if the point is + * not inside this component. + */ + public long findCellContainingXYZ(double[] pos) { + return findCellContainingXYZ( pos[0], pos[1], pos[2] ); + } + + /** + * Returns the cell which contains a given point (x,y,z). + * + * @param pos Hep3Vector with cartesian coordinates of the point + * @return ID of cell containing the point (maybe either in + * absorber or live material), or zero if the point is + * not inside this component. + */ + public long findCellContainingXYZ(Hep3Vector pos) { + return findCellContainingXYZ( pos.x(), pos.y(), pos.z() ); + } +

Commit in `GeomConverter/src/org/lcsim/geometry/segmentation` on MAIN
`BarrelCylinderSegmentationBase.java`	+63		added 1.1
`GridXYZ.java`	+274	-46	1.16 -> 1.17
`NonprojectiveCylinder.java`	+14	-78	1.17 -> 1.18
`ProjectiveCylinder.java`	+41	-1	1.19 -> 1.20
`ProjectiveZPlane.java`	+105	-21	1.17 -> 1.18
`SegmentationBase.java`	+63	-6	1.7 -> 1.8
`SegmentationUtil.java`	+8	-8	1.1 -> 1.2
	+568	-160