6 added + 1 removed + 2 modified, total 9 files
lcsim/src/org/lcsim/digisim
diff -N AbsValueDiscrimination.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ AbsValueDiscrimination.java 27 Jun 2005 21:34:31 -0000 1.1
@@ -0,0 +1,123 @@
+package org.lcsim.digisim;
+
+import org.lcsim.digisim.AbstractCalHitModifier;
+import java.util.Vector;
+import java.util.Map;
+import java.util.Set;
+
+/**
+ * A modifier for basic double-threshold transformations.
+ *
+ * Double here means that the threshold is taken as an absolute value,
+ * but two thresholds are used, so all hits with energies in the range
+ * (-trheshold, +threshold) are discarded. The reason for this is
+ * that in some cases (e.g. gaussian noise), negative noise can cancel
+ * out positive noise.
+ *
+ * @author Guilherme Lima
+ * @version $Id: AbsValueDiscrimination.java,v 1.1 2005/06/27 21:34:31 lima Exp $
+ */
+class AbsValueDiscrimination extends AbstractCalHitModifier {
+
+ /** Returns a reference to a new instance */
+ public AbsValueDiscrimination newInstance(Digitizer digitizer,
+ ModifierParameters modpars) {
+ return new AbsValueDiscrimination(digitizer, modpars);
+ }
+
+ /** private constructor, call newInstance() instead */
+ private AbsValueDiscrimination(Digitizer digi, ModifierParameters modpars) {
+ super(digi, modpars);
+ }
+
+ /** default constructor is only called for global instance */
+ private AbsValueDiscrimination() {
+ super(null, new ModifierParameters( ModifierParameters.Type.AbsValueDiscrimination.toString(), ModifierParameters.Type.AbsValueDiscrimination, new Vector<Double>()) );
+ }
+
+ public static boolean registerMe() {
+ String type = ModifierParameters.Type.AbsValueDiscrimination.toString();
+ AbstractCalHitModifier.registerModifier(type, _me);
+ return true;
+ }
+
+ //*** Defining modifier interface ***
+
+ /** initialize parameters */
+ public void init(Vector<Double> params) {
+ _par = params;
+ if(_debug>0) {
+ int i=0;
+ for( Double x : params ) {
+ System.out.println(_name+".init(): _par["+(i++)+"] = "+x);
+ }
+ }
+ }
+
+ /**
+ * This is where real work is performed.
+ * @param hitmap map with transient hits. Both input and output.
+ */
+ public void processHits(Map<Long,TempCalHit> hitmap) {
+ if(_debug>0) {
+ System.out.println(" name="+_name+", debug="+_debug);
+ }
+
+ // apply threshold
+ Vector<Long> toBeErased = new Vector<Long>();
+ for( Long id : hitmap.keySet() ) {
+ TempCalHit ihit = hitmap.get(id);
+ double energy = ihit.getTotalEnergy();
+
+ if( isBetweenThresholds( energy ) ) {
+ // Tried to erase here, but got problems in C++:
+ // Docs say that iterators
+ // are not safe to use after the erase, and it is used in ++it...
+ // hitmap.erase(it); //
+ // Thus, for now, just save the keys to erase the elements later
+ toBeErased.add( id );
+ }
+ }
+
+ // Remove elements between thresholds from input collection
+ for( Long key : toBeErased ) {
+ hitmap.remove( key );
+ }
+ }
+
+ /** debugging printout */
+ public void print() {
+ System.out.println( this.toString() );
+ System.out.println( _name + ": a gain+threshold modifier" );
+ String ret = new String();
+ if(_par!=null) {
+ ret += " Parameters:";
+ for(int i=0; i<_par.size(); ++i) {
+ ret += " "+ _par.get(i);
+ }
+ }
+ System.out.println(ret);
+ }
+
+ /** transformation functions: discrimination */
+ private boolean isBetweenThresholds(double energy) {
+ // assign roles to the parameters
+ final double threshNominal = _par.get(0);
+ final double threshWidth = _par.get(1);
+
+ double threshold = threshNominal;
+ if(threshWidth>0.0) {
+ threshold = threshNominal + threshWidth * _random.nextGaussian();
+ }
+
+ if( Math.abs(energy)<threshold ) return true;
+ else return false;
+ }
+
+ //=== FIELDS ===
+
+ /** Parameters: gain and threshold */
+ private Vector<Double> _par;
+ /** Global instance, used only for registration */
+ private static AbsValueDiscrimination _me = new AbsValueDiscrimination();
+}
lcsim/src/org/lcsim/digisim
diff -N Crosstalk.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ Crosstalk.java 27 Jun 2005 21:34:31 -0000 1.1
@@ -0,0 +1,204 @@
+package org.lcsim.digisim;
+
+import java.util.Vector;
+import java.util.Map;
+import java.util.HashMap;
+import org.lcsim.geometry.segmentation.SegmentationImpl;
+// import org.lcsim.geometry.util.IDEncoder;
+
+/**
+ * A modifier for simulation of crosstalk in calorimeter cells
+ *
+ * @author Guilherme Lima
+ * @version $Id: Crosstalk.java,v 1.1 2005/06/27 21:34:31 lima Exp $
+ */
+class Crosstalk extends AbstractCalHitModifier {
+
+ //*** Defining modifier interface ***
+
+ /** initialize parameters */
+ public void init(Vector<Double> params) {
+ _par = params;
+ if(_debug>0) {
+ int i=0;
+ for( Double x : params ) {
+ System.out.println(_name+".init(): _par["+(i++)+"] = "+x);
+ }
+ }
+
+ // decode parameters
+ _fracNom = params.get(0);
+ _fracSigma = params.get(1);
+// _timeNom = params.get(2);
+// _timeSigma = params.get(3);
+// _meanNumNoisyCells = params.get(4);
+// _system = params.get(5).intValue();
+// _barend = params.get(6).intValue();
+
+ // @todo we need to load here the number of cells per layer,
+ // for quick determination of neighbors
+
+ // initialize segmentation and encoding
+ CellSelector selector = _digitizer.getCellSelector();
+// _encoder = selector.getIDEncoder();
+ }
+
+ /**
+ * This is where real work is performed.
+ * @param hitmap map with transient hits. Both input and output.
+ */
+ public void processHits(Map<Long,TempCalHit> hitmap) {
+ if(_debug>0) {
+ System.out.println(" name="+_name+", debug="+_debug);
+ }
+
+ _segm = _digitizer.getCellSelector().getSegmentation();
+ assert _segm!=null : "No segmentation available?!";
+
+ if( !_segm.supportsNeighbours() ) return;
+
+ // save input map to be used as reference
+ Map<Long,TempCalHit> origMap = new HashMap<Long,TempCalHit>();
+ // copy map elements, avoiding shallow copies!!
+ for( Long id : hitmap.keySet() ) {
+ TempCalHit oldhit = hitmap.get(id);
+ TempCalHit newhit = new TempCalHit(oldhit);
+ origMap.put(id,newhit);
+ }
+
+ // loop over hit cells
+ for( Long sourceID : origMap.keySet() ) {
+ TempCalHit srchit = origMap.get( sourceID );
+
+ _segm.setID( sourceID.longValue() );
+ long neigs[] = _segm.getNeighbourIDs(0,1,1);
+
+ // loop over neighbours
+ for(int i=0; i<neigs.length; ++i) {
+ long destID = neigs[i];
+
+ // ** temporary crosscheck **
+ assert (sourceID&0xffc00000ffffffffL) == (destID&0xffc00000ffffffffL)
+ : "Error: neighbor with incompatible cellID";
+
+ // Only first neighbors get any crosstalk
+ if( !firstNeighbors(sourceID, destID) ) continue;
+
+ // crosstalk amplitude and timing
+ double frac = _fracNom + _fracSigma * _random.nextGaussian();
+ double ampl = frac * srchit.getTotalEnergy();
+ double time = srchit.getPrimaryTime();
+
+ if(_debug>1) {
+ System.out.println("i="+i+"/"+neigs.length
+ +", ampl="+ampl+", time="+time);
+ }
+
+ // check if neighbour already in map
+ // Note that changes are applied to hitmap, not to origMap
+ TempCalHit desthit = hitmap.get(destID);
+ if( desthit!=null ) {
+ // add crosstalk to existing cell
+ if(_debug>0) {
+ System.out.print("Crosstalk:"
+ +" from sourceID="+Long.toHexString(sourceID)
+ +" (E="+srchit.getTotalEnergy()+")"
+ +" to existing destID="+Long.toHexString(destID)
+ +" (bef="+desthit.getTotalEnergy()
+ +", add="+ampl );
+ }
+
+ desthit.addContribution( sourceID, ampl, time );
+ if(_debug>0) {
+ System.out.println(", after="+desthit.getTotalEnergy()+")");
+ }
+ }
+ else {
+ // create new hit and add it to hitmap
+ desthit = new TempCalHit( destID, sourceID, ampl, time );
+ if(_debug>0) {
+ System.out.println("Crosstalk:"
+ +" from sourceID="+Long.toHexString(sourceID)
+ +" (E="+srchit.getTotalEnergy()+")"
+ +" to new destID="+Long.toHexString(destID)
+ +" (E="+ampl+")" );
+ }
+ hitmap.put( destID, desthit );
+ }
+ }
+ }
+ }
+
+ /**
+ * Return true for first neighbors, false otherwise.\
+ * @param cellid1
+ */
+ public boolean firstNeighbors(long cellid1, long cellid2) {
+ int ithe1 = (int)(cellid1 >> 32) & 0x7ff;
+ int ithe2 = (int)(cellid2 >> 32) & 0x7ff;
+ int iphi1 = (int)(cellid1 >> 43) & 0x7ff;
+ int iphi2 = (int)(cellid2 >> 43) & 0x7ff;
+ int deltheta = Math.abs(ithe2-ithe1);
+ int delphi = Math.abs(iphi2-iphi1);
+ if( ithe1==ithe2 && delphi==1 ) return true;
+ if( iphi1==iphi2 && deltheta==1 ) return true;
+ return false;
+ }
+
+ /** debugging printout */
+ public void print() {
+ System.out.println( this.toString() );
+ System.out.println( _name + ": a crosstalk modifier" );
+ String ret = new String();
+ if(_par!=null) {
+ ret += " Parameters:";
+ for(int i=0; i<_par.size(); ++i) {
+ ret += " "+ _par.get(i);
+ }
+ }
+ System.out.println(ret);
+ }
+
+ /** Returns a reference to a new instance */
+ public Crosstalk newInstance(Digitizer digitizer,
+ ModifierParameters modpars) {
+ return new Crosstalk(digitizer, modpars);
+ }
+
+ /** private constructor, call newInstance() instead */
+ private Crosstalk(Digitizer digitizer, ModifierParameters modpars) {
+ super(digitizer, modpars);
+ }
+
+ /** default constructor is only called for global instance */
+ private Crosstalk() {
+ super( null,
+ new ModifierParameters( ModifierParameters.Type.Crosstalk.toString(),
+ ModifierParameters.Type.Crosstalk,
+ new Vector<Double>()
+ )
+ );
+ }
+
+ public static boolean registerMe() {
+ String type = ModifierParameters.Type.Crosstalk.toString();
+ AbstractCalHitModifier.registerModifier(type, _me);
+ return true;
+ }
+
+ //=== FIELDS ===
+
+ /** Parameters */
+ Vector<Double> _par;
+ /** mean amplitude */
+ private double _fracNom;
+ /** width used for smearing of amplitude */
+ private double _fracSigma;
+
+ /** Global instance, used only for registration */
+ private static Crosstalk _me = new Crosstalk();
+ /** Geometry-aware class to find neighbours */
+ private SegmentationImpl _segm;
+// /** cellID encoder */
+// private IDEncoder _encoder;
+}
lcsim/src/org/lcsim/digisim
diff -N ExponentialNoise.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ ExponentialNoise.java 27 Jun 2005 21:34:31 -0000 1.1
@@ -0,0 +1,123 @@
+package org.lcsim.digisim;
+
+import java.util.Vector;
+import org.apache.commons.math.special.Erf;
+import org.apache.commons.math.distribution.DistributionFactory;
+import org.apache.commons.math.distribution.ExponentialDistribution;
+import org.apache.commons.math.MathException;
+
+/**
+ * Noise parametrization based on an exponential distribution. This
+ * class extends RandomNoise, which internally uses the first 5
+ * parameters. Hence the sixth parameter is the mean of the
+ * exponential distribution.
+ *
+ * @author Guilherme Lima
+ * @version $Id: ExponentialNoise.java,v 1.1 2005/06/27 21:34:31 lima Exp $
+ */
+class ExponentialNoise extends RandomNoise {
+
+ /** Initialize parameters */
+ public void init(Vector<Double> params) {
+ super.init(params);
+
+ if(_debug>0) {
+ int i=0;
+ for( Double x : params ) {
+ System.out.println(_name+".init(): _par["+(i++)+"] = "+x);
+ }
+ }
+
+ // decode parameters
+ _beta = params.get(5);
+
+ // create distribution
+ System.out.println("ExpoNoise.init(): beta="+_beta);
+ _expo = DistributionFactory.newInstance().createExponentialDistribution(_beta);
+ }
+
+ /** Probability to get noise above the threshold */
+ public double getProbabilityAboveThreshold() {
+ _probAboveThreshold = Math.exp( - _threshold / _beta );
+ System.out.println("ExponentialNoise: probAboveThreshold = "
+ +_probAboveThreshold);
+ return _probAboveThreshold;
+ }
+
+ /**
+ * Return a random variable following this distribution.
+ */
+ public double drawRandomNoise() {
+ // gaussian distribution
+ try {
+ double flat = _random.nextDouble();
+ double expo = _expo.inverseCumulativeProbability( flat );
+ org.lcsim.util.aida.AIDA.defaultInstance().cloud1D("allExpo").fill(expo);
+ return expo;
+ }
+ catch(MathException x) {
+ assert false
+ : "ExponentialNoise.drawRandomNoise(): Exception caught: " + x;
+ }
+ return 0;
+ }
+
+ /**
+ * Return a random variable following this distribution, but above
+ * the threshold provided during initialization.
+ */
+ public double drawRandomNoiseAboveThreshold() {
+ double draw = _probAboveThreshold*_random.nextDouble();
+ double cumulProb = 1.0 - _probAboveThreshold + draw;
+ assert cumulProb < 1.0 : "cumulProb="+ cumulProb+", draw="+ draw
+ + ", probAboveThreshold="+ _probAboveThreshold ;
+ assert cumulProb >= 0.0 : "cumulProb="+ cumulProb+", draw="+ draw
+ + ", probAboveThreshold="+ _probAboveThreshold ;
+
+ double noise = 0;
+ try {
+ noise = _expo.inverseCumulativeProbability( cumulProb );
+ }
+ catch(MathException x) {
+ assert false : "MathException caught: " + x;
+ }
+ org.lcsim.util.aida.AIDA.defaultInstance().cloud1D("noiseOnlyExpo").fill(noise);
+ return noise;
+ }
+
+ /** Returns a reference to a new ExponentialNoise modifier */
+ public ExponentialNoise newInstance(Digitizer digitizer,
+ final ModifierParameters modpars) {
+ return new ExponentialNoise(digitizer, modpars);
+ }
+
+ /** Register this modifier, making it available for any digitizer */
+ public static boolean registerMe() {
+ String type = ModifierParameters.Type.ExponentialNoise.toString();
+ AbstractCalHitModifier.registerModifier(type, _me);
+ return true;
+ }
+
+ //*** Interface to be defined in modifiers ***
+
+ /** Constructor */
+ private ExponentialNoise(Digitizer digi, final ModifierParameters modpars) {
+ super(digi,modpars);
+ }
+
+ /** This default constructor is only called for global instance */
+ private ExponentialNoise() {
+ super(null, new ModifierParameters( ModifierParameters.Type.ExponentialNoise.toString(), ModifierParameters.Type.ExponentialNoise, new Vector<Double>()) );
+ }
+
+ //=== FIELDS ===
+
+ /** slope of the exponential */
+ private double _beta;
+ /** probability above threshold */
+ private double _probAboveThreshold;
+ /** gaussian distribution */
+ private ExponentialDistribution _expo;
+ /** Global instance, used only for registration */
+ private static ExponentialNoise _me = new ExponentialNoise();
+}
lcsim/src/org/lcsim/digisim
diff -N GaussianNoise.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ GaussianNoise.java 27 Jun 2005 21:34:31 -0000 1.1
@@ -0,0 +1,170 @@
+package org.lcsim.digisim;
+
+import java.util.Vector;
+import org.apache.commons.math.special.Erf;
+import org.apache.commons.math.distribution.DistributionFactory;
+import org.apache.commons.math.distribution.NormalDistribution;
+import org.apache.commons.math.MathException;
+import org.lcsim.util.aida.AIDA;
+
+/**
+ * Noise parametrization based on a gaussian distribution. This class
+ * extends RandomNoise, which internally uses the first 5 parameters.
+ * Hence the next two parameters are the gaussian mean and width.
+ *
+ * @author Guilherme Lima
+ * @version $Id: GaussianNoise.java,v 1.1 2005/06/27 21:34:31 lima Exp $
+ */
+class GaussianNoise extends RandomNoise {
+
+ /** Initialize parameters */
+ public void init(Vector<Double> params) {
+ super.init(params);
+
+ if(_debug>0) {
+ int i=0;
+ for( Double x : params ) {
+ System.out.println(_name+".init(): _par["+(i++)+"] = "+x);
+ }
+ }
+
+ // decode parameters
+ _mean = params.get(5);
+ _sigma = params.get(6);
+ if(_sigma<0) {
+ _keepNegative = true;
+ _sigma = Math.abs(_sigma);
+ }
+
+ // create distribution
+ _gauss = DistributionFactory.newInstance().createNormalDistribution(_mean,_sigma);
+ }
+
+ //*** Defining modifier interface ***
+
+ /** Returns probability above threshold */
+ public double getProbabilityAboveThreshold() {
+ double z = (_threshold - _mean) / ( Math.sqrt(2.0)*_sigma );
+ // analytically correct, but gives large roundoff errors. GL050627
+ // _probAboveThreshold = 0.5 * ( 1 - Erf.erf(z) );
+
+ _probAboveThreshold = 0.5 * tail( z );
+ System.out.println("GaussianNoise: probAboveThreshold = "
+ +_probAboveThreshold);
+
+ if(_keepNegative) return 2*_probAboveThreshold;
+ return _probAboveThreshold;
+ }
+
+ /**
+ * Return a random variable following this distribution.
+ */
+ public double drawRandomNoise() {
+ // gaussian distribution
+ try {
+ double flat = _random.nextDouble();
+ double gauss = _gauss.inverseCumulativeProbability( flat );
+ AIDA.defaultInstance().cloud1D("drawGauss").fill(gauss);
+ return gauss;
+ }
+ catch(MathException x) {
+ assert false: "GaussianNoise.drawRandomNoise: exception caught: "+x;
+ }
+ return 0;
+ }
+
+ /**
+ * Return a random variable following this distribution, but beyond
+ * threshold provided during initialization.
+ */
+ public double drawRandomNoiseAboveThreshold() {
+ double draw, cumulProb;
+ if(_keepNegative) {
+ draw = 2 * _probAboveThreshold*_random.nextDouble();
+ if(draw>_probAboveThreshold) cumulProb = draw - _probAboveThreshold;
+ else cumulProb = 1.0 - _probAboveThreshold + draw;
+ }
+ else {
+ draw = _probAboveThreshold * _random.nextDouble();
+ cumulProb = 1.0 - _probAboveThreshold + draw;
+ }
+
+ assert cumulProb < 1.0 : "cumulProb="+ cumulProb+", draw="+ draw
+ + ", probAboveThreshold="+ _probAboveThreshold ;
+ assert cumulProb >= 0.0 : "cumulProb="+ cumulProb+", draw="+ draw
+ + ", probAboveThreshold="+ _probAboveThreshold ;
+
+ double noise = 0;
+ try {
+ noise = _gauss.inverseCumulativeProbability( cumulProb );
+ }
+ catch(MathException e) {
+ System.out.println("MathException caught: "+e);
+ }
+ AIDA.defaultInstance().cloud1D("noiseOnlyGauss").fill(noise);
+ return noise;
+ }
+
+ /** Returns a reference to a new GaussianNoise modifier */
+ public GaussianNoise newInstance(Digitizer digitizer,
+ final ModifierParameters modpars) {
+ return new GaussianNoise(digitizer, modpars);
+ }
+
+ /** Register this modifier, making it available for any digitizer */
+ public static boolean registerMe() {
+ String type = ModifierParameters.Type.GaussianNoise.toString();
+ AbstractCalHitModifier.registerModifier(type, _me);
+ return true;
+ }
+
+ //*** Interface to be defined in modifiers ***
+
+ /** Constructor */
+ private GaussianNoise(Digitizer digi, final ModifierParameters modpars) {
+ super(digi,modpars);
+ }
+
+ /** This default constructor is only called for global instance */
+ private GaussianNoise() {
+ super(null, new ModifierParameters( ModifierParameters.Type.GaussianNoise.toString(), ModifierParameters.Type.GaussianNoise, new Vector<Double>()) );
+ }
+
+ /**
+ * This implementation, based on z>>1, avoids roundoff errors on 1-erf(z),
+ * where erf(z) = 1 - tail(z)
+ * from MathWorld at http://mathworld.wolfram.com/Erf.html
+ */
+ private double tail(double x) {
+ double factor = Math.exp(-x*x)/Math.sqrt(Math.PI);
+ double z = 1 / x;
+ double base = z*z / 2;
+
+ double series = z;
+ series += - z*base;
+ series += 3*z*base*base;
+ series += - 15*z*base*base*base;
+ series += 105*z*base*base*base*base;
+ if(_debug>0) {
+ double last = 105*z*base*base*base*base / series;
+ System.out.println(" GaussianNoise.tail(z): estimated error < "
+ +(100*last)+"%");
+ }
+ return factor*series;
+ }
+
+ //=== FIELDS ===
+
+ /** mean amplitude */
+ private double _mean;
+ /** width used for smearing of amplitude */
+ private double _sigma;
+ /** probability above threshold */
+ private double _probAboveThreshold;
+ /** gaussian distribution */
+ private NormalDistribution _gauss;
+ /** flag to control negative values */
+ private boolean _keepNegative;
+ /** Global instance, used only for registration */
+ private static GaussianNoise _me = new GaussianNoise();
+}
lcsim/src/org/lcsim/digisim
diff -N SmearedGain.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ SmearedGain.java 27 Jun 2005 21:34:32 -0000 1.1
@@ -0,0 +1,60 @@
+package org.lcsim.digisim;
+
+import java.util.Vector;
+import java.util.Map;
+import java.util.Set;
+
+/**
+ * A simple modifier for function-based smeared linear transformations
+ * on energy.
+ *
+ * @author Guilherme Lima
+ * @version $Id: SmearedGain.java,v 1.1 2005/06/27 21:34:32 lima Exp $
+ */
+class SmearedGain extends FunctionModifier {
+
+ /** Returns a reference to a new SmearedGain modifier */
+ public SmearedGain newInstance(Digitizer digitizer,
+ final ModifierParameters modpars) {
+ return new SmearedGain(digitizer, modpars);
+ }
+
+ /** Register this modifier, making it available for any digitizer */
+ public static boolean registerMe() {
+ String type = ModifierParameters.Type.SmearedGain.toString();
+ AbstractCalHitModifier.registerModifier(type, _me);
+ return true;
+ }
+
+ //*** Interface to be defined in modifiers ***
+
+ /** Constructor */
+ private SmearedGain(Digitizer digi, final ModifierParameters modpars) {
+ super(digi,modpars);
+ }
+
+ /** This default constructor is only called for global instance */
+ private SmearedGain() {
+ super(null, new ModifierParameters( ModifierParameters.Type.SmearedGain.toString(), ModifierParameters.Type.SmearedGain, new Vector<Double>()) );
+ }
+
+
+ /** Smeared linear transformations on energy */
+ public double transformEnergy(final TempCalHit hit) {
+ // assign roles to the parameters
+ final double gainNominal = _par.get(0);
+ final double gainWidth = _par.get(1);
+
+ double gain = gainNominal;
+ if(gainWidth>0.0) {
+ gain = gainNominal + gainWidth * _random.nextGaussian();
+ }
+ if(_debug>10) System.out.println("gain = " + gain );
+
+ double adc = hit.getTotalEnergy() * gain;
+ return adc;
+ }
+
+ /** Global instance, used only for registration */
+ private static SmearedGain _me = new SmearedGain();
+}
lcsim/src/org/lcsim/digisim
diff -N SmearedTiming.java
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ SmearedTiming.java 27 Jun 2005 21:34:32 -0000 1.1
@@ -0,0 +1,105 @@
+package org.lcsim.digisim;
+
+import org.lcsim.digisim.AbstractCalHitModifier;
+import java.util.Vector;
+import java.util.Map;
+import java.util.Set;
+
+/**
+ * A modifier for basic gain transformations
+ *
+ * @author Guilherme Lima
+ * @version $Id: SmearedTiming.java,v 1.1 2005/06/27 21:34:32 lima Exp $
+ */
+class SmearedTiming extends AbstractCalHitModifier {
+
+ /** Returns a reference to a new instance */
+ public SmearedTiming newInstance(Digitizer digitizer,
+ ModifierParameters modpars) {
+ return new SmearedTiming(digitizer, modpars);
+ }
+
+ /** private constructor, call newInstance() instead */
+ private SmearedTiming(Digitizer digi, ModifierParameters modpars) {
+ super(digi, modpars);
+ }
+
+ /** default constructor is only called for global instance */
+ private SmearedTiming() {
+ super(null, new ModifierParameters( ModifierParameters.Type.SmearedTiming.toString(), ModifierParameters.Type.SmearedTiming, new Vector<Double>()) );
+ }
+
+ public static boolean registerMe() {
+ String type = ModifierParameters.Type.SmearedTiming.toString();
+ AbstractCalHitModifier.registerModifier(type, _me);
+ return true;
+ }
+
+ //*** Defining modifier interface ***
+
+ /** initialize parameters */
+ public void init(Vector<Double> params) {
+ _par = params;
+ if(_debug>0) {
+ int i=0;
+ for( Double x : params ) {
+ System.out.println(_name+".init(): _par["+(i++)+"] = "+x);
+ }
+ }
+ }
+
+ /**
+ * This is where real work is performed.
+ * @param hitmap map with transient hits. Both input and output.
+ */
+ public void processHits(Map<Long,TempCalHit> hitmap) {
+ if(_debug>0) {
+ System.out.println(" name="+_name+", debug="+_debug);
+ }
+
+ // apply gain
+ for( Long id : hitmap.keySet() ) {
+ TempCalHit ihit = hitmap.get(id);
+ double timing = transformTime( ihit );
+// cout<<"hit id="<< hex << it->first << dec
+// << ", Ein="<< ihit.getTotalEnergy()
+// <<", Eout="<< energy << endl;
+ ihit.setTime( timing );
+ }
+ }
+
+ /** Smeared linear transformations on time */
+ public double transformTime(final TempCalHit hit) {
+ // assign roles to the parameters
+ final double smearNominal = _par.get(2);
+ final double smearWidth = _par.get(3);
+
+ double smear = smearNominal;
+ if(smearWidth>0.0) {
+ smear = smearNominal + smearWidth*_random.nextGaussian();
+ }
+ double newTime = hit.getPrimaryTime() * smear;
+ return newTime;
+ }
+
+ /** debugging printout */
+ public void print() {
+ System.out.println( this.toString() );
+ System.out.println( _name + ": a gain+threshold modifier" );
+ String ret = new String();
+ if(_par!=null) {
+ ret += " Parameters:";
+ for(int i=0; i<_par.size(); ++i) {
+ ret += " "+ _par.get(i);
+ }
+ }
+ System.out.println(ret);
+ }
+
+ //=== FIELDS ===
+
+ /** Parameters: gain and threshold */
+ private Vector<Double> _par;
+ /** Global instance, used only for registration */
+ private static SmearedTiming _me = new SmearedTiming();
+}
lcsim/src/org/lcsim/digisim
diff -u -r1.3 -r1.4
--- ModifierParameters.java 25 May 2005 17:50:22 -0000 1.3
+++ ModifierParameters.java 27 Jun 2005 21:34:32 -0000 1.4
@@ -6,7 +6,7 @@
* Encapsulates the configuration parameters for any DigiSim modifier.
*
* @author Guilherme Lima
- * @version $Id: ModifierParameters.java,v 1.3 2005/05/25 17:50:22 lima Exp $
+ * @version $Id: ModifierParameters.java,v 1.4 2005/06/27 21:34:32 lima Exp $
*/
public class ModifierParameters {
/**
@@ -26,8 +26,12 @@
/** Enumerates all possible modifier types */
public enum Type {
GainDiscrimination,
- SmearedLinear,
- RandomNoise,
+ AbsValueDiscrimination,
+ SmearedGain,
+ SmearedTiming,
+ Crosstalk,
+ GaussianNoise,
+ ExponentialNoise,
HotCell,
DeadCell,
SiPMSaturation
lcsim/src/org/lcsim/digisim
diff -u -r1.3 -r1.4
--- DigiSimMain.java 25 May 2005 17:50:21 -0000 1.3
+++ DigiSimMain.java 27 Jun 2005 21:34:32 -0000 1.4
@@ -3,6 +3,7 @@
import java.util.List;
import java.util.Vector;
import java.util.Map;
+import java.util.HashMap;
import java.io.File;
import org.lcsim.util.Driver;
@@ -20,7 +21,7 @@
* The main driver for digitization simulation.
*
* @author Guilherme Lima
- * @version $Id: DigiSimMain.java,v 1.3 2005/05/25 17:50:21 lima Exp $
+ * @version $Id: DigiSimMain.java,v 1.4 2005/06/27 21:34:32 lima Exp $
*/
public class DigiSimMain extends Driver {
@@ -46,7 +47,7 @@
_hitmgr.process(event);
_digi.process(event);
- // test calhitmapmgr
+ //.. test calhitmapmgr
CalHitMapMgr hitmgr = CalHitMapMgr.getInstance();
Map<Long,SimCalorimeterHit> emhits, hadhits;
emhits = hitmgr.getCollHitMap("EcalBarrHits");
@@ -54,41 +55,127 @@
// System.out.println("# hits retrieved:"
// +" EM="+emhits.size()+", HAD="+hadhits.size());
- // make sure digisim works
+ //.. test digisim works
// using fixed factors, compare sim and raw hits
- List<LCRelation> raw2simLinks = event.get( LCRelation.class,
- "EcalBarrRaw2sim" );
- for( LCRelation rel : raw2simLinks ) {
- RawCalorimeterHit rawhit = (RawCalorimeterHit)rel.getFrom();
- SimCalorimeterHit simhit = (SimCalorimeterHit)rel.getTo();
- int simE = (int)(simhit.getEnergy()*1000000.0);
- if(simE!=rawhit.getAmplitude()) {
- System.out.println("Discrepancy: simE=" + simE
- + ", rawE=" + rawhit.getAmplitude() );
+ try {
+ List<LCRelation> raw2simLinks
+ = event.get( LCRelation.class, "EcalBarrRaw2sim" );
+ for( LCRelation rel : raw2simLinks ) {
+ RawCalorimeterHit rawhit = (RawCalorimeterHit)rel.getFrom();
+ SimCalorimeterHit simhit = (SimCalorimeterHit)rel.getTo();
+ int simE = (int)(simhit.getEnergy()*1000000.0);
+ assert simE == rawhit.getAmplitude()
+ : "Discrepancy on EcalBarr: simE=" + simE
+ + ", rawE=" + rawhit.getAmplitude() ;
}
}
+ catch(IllegalArgumentException e) {
+ // no problem if data is not there
+ }
+
+ // Test crosstalk modifier: Use fixed 10% values, check neighbours,
+ // weights, etc.
+ try {
+ List<SimCalorimeterHit> simhits
+ = event.get(SimCalorimeterHit.class, "EcalEndcapHits");
+ List<RawCalorimeterHit> rawhits
+ = event.get(RawCalorimeterHit.class, "EcalEndcapRawHits");
+ List<LCRelation> raw2simLinks
+ = event.get(LCRelation.class, "EcalEndcapRaw2sim");
+
+ System.out.println("Checking crosstalk for EcalEndcaps...");
+ System.out.println("# hits: sim="+( simhits==null? 0 : simhits.size() )
+ +", raw=" +( rawhits==null? 0 : rawhits.size() )
+ +" - #links="+( raw2simLinks==null? 0: raw2simLinks.size()));
+
+ Map<Long,Double> weightMap = new HashMap<Long,Double>();
+ int ilink = 0;
+ for( LCRelation rel : raw2simLinks ) {
+ ++ilink;
+ RawCalorimeterHit rawhit = (RawCalorimeterHit)rel.getFrom();
+ SimCalorimeterHit simhit = (SimCalorimeterHit)rel.getTo();
+
+ Long rawid = new Long(rawhit.getCellID());
+ long simid = (simhit!=null ? simhit.getCellID() : 0);
+ Double temp = weightMap.get(rawid);
+ if(temp==null) {
+ temp = new Double( rel.getWeight() );
+ weightMap.put( rawid, temp );
+ }
+ else {
+ temp += rel.getWeight();
+// double newval = temp.doubleValue() + rel.getWeight();
+// weightMap.put( new Long(rawid), new Double(rel.getWeight()) );
+ }
+
+ // discard links where simid==rawid or simid==zero
+ if(rawid==simid) {
+ System.out.println("link "+ilink+" not crosstalk.");
+ continue;
+ }
+ if(simid==0.0) {
+ System.out.println("link "+ilink+", simid=zero.");
+ continue;
+ }
+
+ // make sure sim,raw cells are good neighbors
+// IDDecoder decoder = new IDDecoder();
+ assert ((simid & 0x7f) == (rawid & 0x7f))
+ : "link "+ilink+", layer discrepancy: "+ simid +" "+rawid;
+ assert ((simid & 0x1f80) == (rawid & 0x1f80))
+ : "link "+ilink+", system+barrel discrepancy: "+ simid +" "+rawid;
+
+ // check theta,phi distance
+ int rawthe = (int)(rawid >> 32) & 0x7ff;
+ int simthe = (int)(simid >> 32) & 0x7ff;
+ int rawphi = (int)(rawid >> 43) & 0x7ff;
+ int simphi = (int)(simid >> 43) & 0x7ff;
+ int delthe = rawthe - simthe + 1;
+ int delphi = rawphi - simphi + 1;
+ assert ((delthe>=0) && (delthe<=2))
+ : "link "+ilink+", theta discrepancy: "+ simid +" "+rawid;
+ assert ((delphi>=0) && (delphi<=2))
+ : "link "+ilink+", phi discrepancy: "+ simid +" "+rawid;
+
+ int simE = (int)(simhit.getEnergy()*1000000.0);
+ int rawE = (int)(rawhit.getAmplitude());
+ System.out.println("X-check: simid="+Long.toHexString(simid)
+ +", simE=" + simE
+ + " - RawID="+Long.toHexString(rawid.longValue())
+ +", rawE=" + rawE );
+ }
+ }
+ catch(IllegalArgumentException e) {
+ // a collection requested does not exist in current event
+ }
}
/** Debugging printout control */
public void setDebug() {
// format: "Digitizer:Modifier"
-// _digi.setDebug("EcalEndcapDigitizer:EMECFixedGain", 1);
-// _digi.setDebug("HcalBarrDigitizer:HBSiPMSaturat", 1);
+// _digi.setDebug("EcalEndcapDigitizer:EMECcrosstalk", 1);
+// _digi.setDebug("HcalBarrDigitizer:HBcrosstalk", 1);
}
//***** FIELDS *****
// event counter
private int _nevt;
+
// subdrivers
private CalHitMapDriver _hitmgr;
private DigiSimDriver _digi;
+
// static registration globals
private static boolean _gaindisc = GainDiscrimination.registerMe();
- private static boolean _smearlin = SmearedLinear.registerMe();
- private static boolean _rndnoise = RandomNoise.registerMe();
+ private static boolean _gaussnoise = GaussianNoise.registerMe();
+ private static boolean _exponoise = ExponentialNoise.registerMe();
private static boolean _hotcell = HotCell.registerMe();
private static boolean _deadcell = DeadCell.registerMe();
private static boolean _sipmsatur = SiPMSaturation.registerMe();
+ private static boolean _xtalk = Crosstalk.registerMe();
+ private static boolean _absdisc = AbsValueDiscrimination.registerMe();
+ private static boolean _gain = SmearedGain.registerMe();
+ private static boolean _timing = SmearedTiming.registerMe();
//***** MAIN ENTRY POINT *****
@@ -114,6 +201,11 @@
loop.loop(-1); // -1 for all
loop.dispose();
+ org.lcsim.util.aida.AIDA.defaultInstance().saveAs("myhistos.aida");
System.out.println("Analyzed "+digidriver._nevt+" events");
}
+
+ private void log(String text) {
+ System.out.println(text);
+ }
}
lcsim/src/org/lcsim/digisim
diff -N SmearedLinear.java
--- SmearedLinear.java 25 May 2005 17:50:27 -0000 1.4
+++ /dev/null 1 Jan 1970 00:00:00 -0000
@@ -1,74 +0,0 @@
-package org.lcsim.digisim;
-
-import java.util.Vector;
-import java.util.Map;
-import java.util.Set;
-
-/**
- * A simple modifier for function-based smeared linear transformations
- * on energy and time.
- *
- * @author Guilherme Lima
- * @version $Id: SmearedLinear.java,v 1.4 2005/05/25 17:50:27 lima Exp $
- */
-class SmearedLinear extends FunctionModifier {
-
- /** Returns a reference to a new SmearedLinear modifier */
- public SmearedLinear newInstance(Digitizer digitizer,
- final ModifierParameters modpars) {
- return new SmearedLinear(digitizer, modpars);
- }
-
- /** Register this modifier, making it available for any digitizer */
- public static boolean registerMe() {
- String type = ModifierParameters.Type.SmearedLinear.toString();
- AbstractCalHitModifier.registerModifier(type, _me);
- return true;
- }
-
- //*** Interface to be defined in modifiers ***
-
- /** Constructor */
- private SmearedLinear(Digitizer digi, final ModifierParameters modpars) {
- super(digi,modpars);
- }
-
- /** This default constructor is only called for global instance */
- private SmearedLinear() {
- super(null, new ModifierParameters( ModifierParameters.Type.SmearedLinear.toString(), ModifierParameters.Type.SmearedLinear, new Vector<Double>()) );
- }
-
-
- /** Smeared linear transformations on energy */
- public double transformEnergy(final TempCalHit hit) {
- // assign roles to the parameters
- final double gainNominal = _par.get(0);
- final double gainWidth = _par.get(1);
-
- double gain = gainNominal;
- if(gainWidth>0.0) {
- gain = gainNominal + gainWidth * _random.nextGaussian();
- }
- if(_debug>10) System.out.println("gain = " + gain );
-
- double adc = hit.getTotalEnergy() * gain;
- return adc;
- }
-
- /** Smeared linear transformations on time */
- public double transformTime(final TempCalHit hit) {
- // assign roles to the parameters
- final double smearNominal = _par.get(2);
- final double smearWidth = _par.get(3);
-
- double smear = smearNominal;
- if(smearWidth>0.0) {
- smear = smearNominal + smearWidth*_random.nextGaussian();
- }
- double newTime = hit.getPrimaryTime() * smear;
- return newTime;
- }
-
- /** Global instance, used only for registration */
- private static SmearedLinear _me = new SmearedLinear();
-}
CVSspam 0.2.8