GeomConverter/src/org/lcsim/geometry/segmentation
diff -u -r1.24 -r1.25
--- GridXYZ.java 20 May 2008 18:19:40 -0000 1.24
+++ GridXYZ.java 23 May 2008 00:27:33 -0000 1.25
@@ -101,6 +101,7 @@
double dy = gridSizeY;
int zBins = detector.getLayering().getLayerCount();
double Rmax = ((CylindricalCalorimeter) detector).getOuterRadius();
+ double Rmin = this.getRMin();
int size = 0;
int refLayer = values[layerIndex];
@@ -115,6 +116,8 @@
int xBins = (int) Math.ceil(2 * Rmax / dx);
int yBins = xBins;
+ int rBinMin = (int) Math.floor(Rmin / dx); // Inner radius, in units of cell width
+ int rBinMax = (int) Math.ceil(Rmax / dx); // Outer radius, in units of cell width
// project theta,phi to sensitive middle of current layer
double z = getDepthSensitiveMid(ilay);
@@ -155,7 +158,26 @@
// System.out.println("Adding neighbor: "+ilay+" "+ix+" "
// +iy+", total="+ (size+1));
- result[size++] = encoder.setValue(yIndex, iy);
+ long candidate = encoder.setValue(yIndex, iy);
+ if (ix*ix + iy*iy >= (rBinMax-1)*(rBinMax-1) || ix*ix + iy*iy <= (rBinMin+1)*(rBinMin+1)) {
+ // We're near the edge. Verify that this is actually a legal value:
+ this.setID(candidate);
+ double candAbsZ = Math.abs(this.getZ());
+ double candRhoSq = (this.getX()*this.getX() + this.getY()*this.getY());
+ boolean candValid = (candAbsZ >= getZMin() && candAbsZ <= getZMax() && candRhoSq >= (getRMin()*getRMin()) && candRhoSq <= (getRMax()*getRMax()));
+ if (candValid) {
+ // OK -- neighbour is at a legal spacepoint
+ result[size++] = candidate;
+ } else {
+ // Failure -- this can happen if the "neighbour" was just us wandering off the grid
+ // Not serious -- don't add the candidate to the output, but don't panic either.
+ }
+ // Reset ID to previous value
+ this.setID(id);
+ } else {
+ // We're far from the edge. Go ahead and add it without spending CPU time on a check
+ result[size++] = candidate;
+ }
}
}
}
GeomConverter/src/org/lcsim/geometry/segmentation
diff -u -r1.25 -r1.26
--- NonprojectiveCylinder.java 20 May 2008 18:19:40 -0000 1.25
+++ NonprojectiveCylinder.java 23 May 2008 00:27:33 -0000 1.26
@@ -138,6 +138,7 @@
IDEncoder gnEncoder = new IDEncoder(descriptor);
BaseIDDecoder gnDecoder = new BaseIDDecoder(descriptor);
gnEncoder.setValues(values);
+ long origID = gnEncoder.getID();
gnDecoder.setID(gnEncoder.getID());
int nMax = (2*layerRange + 1)*(2*zRange + 1)*(2*phiRange + 1) - 1;
@@ -147,7 +148,40 @@
double theta = getTheta();
double phi = getPhi();
- int zBins = (int) Math.floor((getZMax() - getZMin()) / getGridSizeZ());
+ // Be careful with # of Z bins.
+ // The rules are:
+ // * A cell is valid if its center has (z >= getZMin() && z <= getZMax())
+ // * The center of a bin is always at z = (integer + 0.5) * gridSizeZ
+ // => Bins are arranged about z=0, with a bin edge at z=0
+ // => If getZMin()==-getZMax() then there is always an even number of bins
+ //
+ // If things are nicely symmetric, i.e. getZMin()==-getZMax(), then this turns
+ // into a very easy problem: if we have n bins ON EACH SIDE (for integer n) then
+ // n - 0.5 <= zMax/gridSizeZ < n + 0.5
+ // and so we can work backwards (and add in the other half) to obtain the total
+ // number of bins 2n as 2*Math.round((zMax-zMin)/(2.0*gridSizeZ)).
+ //
+ // If things are not symmetric but span the z=0 point then we're OK as long as
+ // we're careful -- the easiest way to look at it is to do each half of the
+ // detector separately:
+ // Number of bins in z>0: Math.round((zMax - 0.0)/gridSizeZ)
+ // Number of bins in z<0: Math.round((0.0 - zMin)/gridSizeZ)
+ //
+ // If things are asymmetric and do not include z=0 then things are more tricky again.
+ // We'd have to work like this:
+ // a) find the first bin-center z1 after GetZMin()
+ // b) find the last bin-center z2 before GetZMax()
+ // c) count number of bins as Math.round((z2-z1)/gridSizeZ)
+ // where
+ // z1 = gridSizeZ * (Math.ceil ((zMin/gridSizeZ)+0.5) - 0.5)
+ // z2 = gridSizeZ * (Math.floor((zMax/gridSizeZ)-0.5) + 0.5)
+ //
+ // So the most general form for the number of bins should be
+ // Math.round( (gridSizeZ*(Math.floor((zMax/gridSizeZ)-0.5)+0.5) - gridSizeZ*(Math.floor((zMax/gridSizeZ)-0.5)+0.5) ) / gridSizeZ )
+ //
+ // ... but we will assume that the detector is symmetric here -- that is true
+ // for any barrel calorimeter generated from a compact.xml geometry file.
+ int zBins = 2 * (int) Math.round( (getZMax()-getZMin())/(2.0*getGridSizeZ()) );
int size = 0;
for (int i = -layerRange; i <= layerRange; ++i)
@@ -159,24 +193,34 @@
gnEncoder.setValue(layerIndex, ilay);
double dphi = this.computeDeltaPhiForLayer(ilay);
- int phiBins = (int) Math.floor(2 * Math.PI / dphi);
+ int phiBins = (int) Math.round(2 * Math.PI / dphi); // Use round() not floor() since a floor() was already applied in the definition of dphi
- double cylR = getRadiusSensitiveMid(ilay);
- double x = cylR * Math.cos(phi);
- double y = cylR * Math.sin(phi);
- double z = cylR / Math.tan(theta);
-
- long id = this.findCellContainingXYZ(x,y,z);
- if(id==0) continue;
+ if (i != 0) {
+ double cylR = getRadiusSensitiveMid(ilay);
+ double x = cylR * Math.cos(phi);
+ double y = cylR * Math.sin(phi);
+ double z = cylR / Math.tan(theta);
+
+ long id = this.findCellContainingXYZ(x,y,z);
+ if(id==0) continue;
+ // save indices in a new array, as values[] keeps the original ref
+ gnDecoder.setID(id);
+ } else {
+ // This is the original layer => center cell is just the original cell
+ gnDecoder.setID(origID);
+ }
- // save indices in a new array, as values[] keeps the original ref
- gnDecoder.setID(id);
for (int j = -zRange; j <= zRange; ++j)
{
int iz = gnDecoder.getValue(zIndex) + j;
- if (iz < -zBins / 2 || iz >= zBins / 2)
+ if (iz < -zBins / 2 || iz >= zBins / 2) {
+ // We make the implicit assumption that detector is symmetric, and so number
+ // of bins is even and looks like [-n, +n-1]. For example, for the
+ // very simple case of two bins, the range of bin indices is {-1, 0}.
+ // In this instance we're outside the valid range, so skip this bin.
continue;
+ }
gnEncoder.setValue(zIndex, iz);
for (int k = -phiRange; k <= phiRange; ++k)
