diff --git a/DDSegmentation/include/DDSegmentation/MegatileLayerGridXY.h b/DDSegmentation/include/DDSegmentation/MegatileLayerGridXY.h
index 01516591bd7f7702f9f3cd7622c8983007a7e865..fbba4557b1b611c88b7d19d3137eed68e3a6e6f9 100644
--- a/DDSegmentation/include/DDSegmentation/MegatileLayerGridXY.h
+++ b/DDSegmentation/include/DDSegmentation/MegatileLayerGridXY.h
@@ -10,200 +10,135 @@
 
 #include "DDSegmentation/CartesianGrid.h"
 
+#include <cassert>
+
 #define MAX_LAYERS 100
-#define MAX_WAFERS 100
+
+/*
+
+a megatile is a rectangule in x-y, split into a grid along x and y, with an exactly integer number of cells in x and y.
+
+this class assumes a mostly-common megatile size, with possibility for a number of "special" megatiles of non-standard size / segmentation
+
+the segmentation of standard megatiles is defined layer-by-layer.
+
+some changes wrt previous version from Kotera et al.
+- significantly simplified. 
+- complications due to end-of-slab moved to higher level detector drivers.
+
+D. Jeans - Nov 2016
+
+ */
+
 
 namespace DD4hep {
-namespace DDSegmentation {
-
-class MegatileLayerGridXY: public CartesianGrid {
-public:
-	/// Default constructor passing the encoding string
-	MegatileLayerGridXY(const std::string& cellEncoding = "");
-	/// Default constructor used by derived classes passing an existing decoder
-	MegatileLayerGridXY(BitField64* decoder);
-	/// destructor
-	virtual ~MegatileLayerGridXY();
-
-	/// determine the position based on the cell ID
-	virtual Vector3D position(const CellID& cellID) const;
-	/// determine the cell ID based on the position
-	virtual CellID cellID(const Vector3D& localPosition, const Vector3D& globalPosition, const VolumeID& volumeID) const;
-
-	/// access the field name used for X
-	const std::string& fieldNameX() const {
-		return _xId;
-	}
-	/// access the field name used for Y
-	const std::string& fieldNameY() const {
-		return _yId;
-	}
-
-	// total size of surface in X,Y
-//	double totalSizeX(int inLayer)  {
-//		return _totalSizeX[ inLayer ];
-//	}
-//	double totalSizeY() const  {
-//		return _totalSizeY;
-//	}
-
-	// number of megatiles along X
-//	int nMegaX(int inLayer=0) const {
-//		return _nMegaY[inLayer];
-//	}
-
-	// the dead region at edge of megatiles (e.g. guard ring width)
-//	double deadWidth() {
-//		return _deadWidth;
-//	}
-	
-		int NStripsY() {
-			return _nStripsY;
-		}
-
-		int NStripsX() {
-			return _nStripsX;
-		}
-
-		int NMegaY() {
-			return _nMegaY;
-		}
-	/// set the coordinate offset in X
-//	void setSlabOffsetX( int layer, double offset) {
-//	  _slabOffsetX[layer] = offset;
-//	}
-	/// set the coordinate offset in Y
-//	void setSlabOffsetY( int layer, double offset) {
-//	  _slabOffsetY[layer] = offset;
-//	}
-
-	void setIsRegulatingEBU( int layer, int wafer, bool isRegulating ) {
-		_isRegulatingEBU[layer][wafer] = isRegulating;
-	}
-
-	/// set the coordinate offset in X
-	void setWaferOffsetX( int layer, int wafer, double offset) {
-	  _waferOffsetX[layer][wafer] = offset;
-	}
-	/// set the coordinate offset in Y
-	void setWaferOffsetY( int layer, int wafer, double offset) {
-	  _waferOffsetY[layer][wafer] = offset;
-	}
-
-	// total size of surface in X,Y
-  void setTotalSizeX(int inLayer, double x) {
-    _totalSizeX[inLayer] = x;
-	  _calculated=false;
-	}
-        
-	void setTotalSizeY( double y) {
-    _totalSizeY = y;
-	  _calculated=false;
-	}
-
-	// number of megatiles along X
-//	void setNMegaY( int n) {
-//		_nMegaY = n;
-//		_calculated=false;
-//	}
-
-	// the dead region at edge of megatiles (e.g. guard ring width)
-//	void setDeadWidth(double x) {
-//		_deadWidth = x;
-//	  _calculated=false;
-//	}
-
-
-	/// set the field name used for X
-	void setFieldNameX(const std::string& fieldName) {
-		_xId = fieldName;
-	}
-	/// set the field name used for Y
-	void setFieldNameY(const std::string& fieldName) {
-		_yId = fieldName;
-	}
-
-//	void setWaferIndexX(double waferIndexX) {
-//			_waferIndexX = waferIndexX;
-//			_calculated=false;
-//	}
-
-//	void setWaferIndexY(double waferIndexY) {
-//			_waferIndexY = waferIndexY;
-//			_calculated=false;
-//	}
-
-	/** \brief Returns a vector<double> of the cellDimensions of the given cell ID
-	    in natural order of dimensions, e.g., dx/dy/dz, or dr/r*dPhi
-
-	    Returns a vector of the cellDimensions of the given cell ID
-	    \param cellID is ignored as all cells have the same dimension
-	    \return std::vector<double> size 2:
-	    -# size in x
-	    -# size in y
-	*/
-	virtual std::vector<double> cellDimensions(const CellID& cellID) const;
-	virtual std::vector<double> layercellDimensions(const int ilayer) const;
-
-protected:
-
-	// total size of surface in X,Y
-	double  _totalSizeX[MAX_LAYERS];
-	double  _totalSizeY;
-  // number of megatiles along X : together with _totalSizeX, this defines the size of the square megatile
-	int  _nMegaY;
-  // number of cells per megatile in X, Y
-	int  _nStripsX;
-	int  _nStripsY;
-	int  _nCells;
-
-	// the dead region at edge of megatiles (e.g. guard ring width): assumed constant in each layer
-	double _deadWidth;
-
-  // size of megatile (including edge region)
-	mutable double  _megaSize;
-	// the grid size in X
-	mutable double  _gridSizeT;
-	// the grid size in Y
-	mutable double  _gridSizeL;
-	mutable double  _gridSizeS;
-
-	// calculate the derived quantities
-	void calculateDerivedQuantities() const;
-	mutable bool   _calculated;
-
-	/// the coordinate offset in X
-	double _offsetX;
-	/// the coordinate offset in Y
-	double _offsetY;
-	double _slabOffsetX;
-	double _slabOffsetY;
-	double _waferOffsetX[MAX_LAYERS][MAX_WAFERS];
-	double _waferOffsetY[MAX_LAYERS][MAX_WAFERS];
-
-	unsigned int _waferIndexX;
-	unsigned int _waferIndexY;
-
-	bool _isRegulatingEBU[MAX_LAYERS][MAX_WAFERS];
-
-	/// the field name used for X
-	std::string _xId;
-	/// the field name used for Y
-	std::string _yId;
-	/// encoding field used for the layer
-	std::string _identifierLayer; 
-	/// encoding field used for the wafer
-	std::string _identifierWafer; 
-//  std::string _identifierWaferY; 
-
-	//	std::string _identifierGR;
-
-	std::string _layerConfig;
-
-	std::string _identifierModule;
-
-};
-
-} /* namespace DDSegmentation */
+  namespace DDSegmentation {
+
+    class MegatileLayerGridXY: public CartesianGrid {
+    public:
+      /// Default constructor passing the encoding string
+      MegatileLayerGridXY(const std::string& cellEncoding = "");
+      /// destructor
+      virtual ~MegatileLayerGridXY();
+
+      /// determine the position based on the cell ID
+      virtual Vector3D position(const CellID& cellID) const;
+      /// determine the cell ID based on the position
+      virtual CellID cellID(const Vector3D& localPosition, const Vector3D& globalPosition, const VolumeID& volumeID) const;
+
+      // set size of megatile in X,Y
+      void setMegaTileSizeXY(double x, double y) {
+        _megaTileSizeX = x;
+        _megaTileSizeY = y;
+      }
+
+      /// set the coordinate offset in X, Y
+      void setMegaTileOffsetXY( double x, double y) {
+        _megaTileOffsetX = x;
+        _megaTileOffsetY = y;
+      }
+
+      void setMegaTileCellsXY( unsigned int layer, int ix, int iy ) {
+	assert ( layer < MAX_LAYERS );
+	_nCellsX[layer] = ix;
+	_nCellsY[layer] = iy;
+      }
+
+      void setSpecialMegaTile( unsigned int layer, unsigned int tile, 
+			       double sizex, double sizey,
+			       double offsetx, double offsety,
+			       unsigned int ncellsx, unsigned int ncellsy );
+
+      /// access the field name used for X
+      const std::string& fieldNameX() const {
+        return _xId;
+      }
+      /// access the field name used for Y
+      const std::string& fieldNameY() const {
+        return _yId;
+      }
+      
+      /// set the field name used for X
+      void setFieldNameX(const std::string& fieldName) {
+        _xId = fieldName;
+      }
+      /// set the field name used for Y
+      void setFieldNameY(const std::string& fieldName) {
+        _yId = fieldName;
+      }
+
+      virtual std::vector<double> cellDimensions(const CellID& cellID) const;
+      virtual std::vector<double> cellDimensions(const unsigned int ilayer, const unsigned int iwafer) const;
+
+      
+    protected:
+
+      struct segInfo {
+	double megaTileSizeX;
+	double megaTileSizeY;
+	double megaTileOffsetX;
+	double megaTileOffsetY;
+	unsigned int nCellsX;
+	unsigned int nCellsY;
+      };
+
+      mutable segInfo _currentSegInfo;
+
+      void getSegInfo( unsigned int layerIndex, unsigned int waferIndex) const;
+
+      // the "usual" megatiles
+      //  megatile size and offset is constant in all layers
+      //  the segmentation may change layer-to-layer (e.g. orthogonal strips)
+
+      // total size of surface in X,Y
+      double  _megaTileSizeX; // [MAX_LAYERS][MAX_WAFERS];
+      double  _megaTileSizeY; //[MAX_LAYERS][MAX_WAFERS];
+
+      double  _megaTileOffsetX; 
+      double  _megaTileOffsetY; 
+
+      // number of cells per megatile in X, Y
+      unsigned int _nCellsX[MAX_LAYERS];
+      unsigned int _nCellsY[MAX_LAYERS];
+
+      std::map < std::pair < unsigned int, unsigned int > , segInfo > specialMegaTiles_layerWafer;
+
+      
+      /// the field name used for X
+      std::string _xId;
+      /// the field name used for Y
+      std::string _yId;
+      /// encoding field used for the layer
+      std::string _identifierLayer;
+      /// encoding field used for the wafer
+      std::string _identifierWafer;
+
+      std::string _layerConfig;
+
+      std::string _identifierModule;
+
+    };
+
+  } /* namespace DDSegmentation */
 } /* namespace DD4hep */
-#endif /* DDSegmentation_MEGATILELAYERGRIDXY_H_ */
+#endif /* DDSegmentation_WAFERGRIDXY_H_ */
diff --git a/DDSegmentation/src/MegatileLayerGridXY.cpp b/DDSegmentation/src/MegatileLayerGridXY.cpp
index 93b77bd2cba828dd70b9faa7c8cd6477fe4975bb..5446eee02364f04518e6d9911e297a14666912a4 100644
--- a/DDSegmentation/src/MegatileLayerGridXY.cpp
+++ b/DDSegmentation/src/MegatileLayerGridXY.cpp
@@ -6,13 +6,13 @@
  *              D Jeans UTokyo
  */
 
-//#define VERBOSEcoterra 1
-
-
 #include "DDSegmentation/MegatileLayerGridXY.h"
 
+#undef NDEBUG
 #include <cassert>
 
+#include <algorithm>
+
 namespace DD4hep {
   namespace DDSegmentation {
 
@@ -21,62 +21,7 @@ namespace DD4hep {
       CartesianGrid(cellEncoding) {
       // define type and description
       _type = "MegatileLayerGridXY";
-      // _description = "Cartesian segmentation in the local XY-plane for both Normal wafer and Magic wafer(depending on the layer dimensions)";
-      _description = "Cartesian segmentation in the local XY-plane: megatiles with dead areas; integer number of megatiles and cells";
-
-//////// setup in the compact-steering file (i.e., ILD_o3_v05.xml)//// 
-// for example,
-// <segmentation type="MegatileLayerGridXY" yer_nCells="36" layer_nStripsX="2" layer_nStripsY="36" deadWidth="0" layer_configuration="TL"/>
-////
-//
-      registerParameter("nMegaY", "number of megatiles along Z", _nMegaY , 1 , SegmentationParameter::NoUnit, true);
-//      _nMegaY is given by Ecal_n_wafers_per_tower from "model_parameters_ILD_o3_v05.xml"
-
-      registerParameter("layer_nCells", "division of megatile into square tiles", _nCells, 36 , SegmentationParameter::NoUnit, true);
-      registerParameter("layer_nStripsX", "division of megatile into strips (1)", _nStripsX, 4 , SegmentationParameter::NoUnit, true);
-      registerParameter("layer_nStripsY", "division of megatile into strips (2)", _nStripsY, 36 , SegmentationParameter::NoUnit, true);
-
-      registerParameter("deadWidth", "width of dead region at edge of megatile", _deadWidth, 0., SegmentationParameter::LengthUnit, true);
-
-      registerIdentifier("identifier_x", "Cell ID identifier for X", _xId, "cellX");
-      registerIdentifier("identifier_y", "Cell ID identifier for Y", _yId, "cellY");
-
-      registerParameter("identifier_wafer", "Cell encoding identifier for wafer", _identifierWafer, std::string("wafer"),
-                        SegmentationParameter::NoUnit, true);
-
-      registerParameter("identifier_layer", "Cell encoding identifier for layer", _identifierLayer, std::string("layer"),
-                        SegmentationParameter::NoUnit, true);
-
-      registerParameter("identifier_module", "Cell encoding identifier for module", _identifierModule, std::string("module"),
-                        SegmentationParameter::NoUnit, true);
-
-      registerParameter("layer_configuration", "layer configuration (S, T, L)", _layerConfig, std::string("TLS"), SegmentationParameter::NoUnit, true);
-
-      _calculated=false;
-
-    }
-
-    /// default constructor using an encoding string
-    MegatileLayerGridXY::MegatileLayerGridXY(BitField64* decode) :
-      CartesianGrid(decode) {
-      // define type and description
-      _type = "MegatileLayerGridXY";
-      // _description = "Cartesian segmentation in the local XY-plane for both Normal wafer and Magic wafer(depending on the layer dimensions)";
-      _description = "Cartesian segmentation in the local XY-plane: megatiles with dead areas; integer number of megatiles and cells";
-
-//////// setup in the compact-steering file (i.e., ILD_o3_v05.xml)//// 
-// for example,
-// <segmentation type="MegatileLayerGridXY" yer_nCells="36" layer_nStripsX="2" layer_nStripsY="36" deadWidth="0" layer_configuration="TL"/>
-////
-//
-      registerParameter("nMegaY", "number of megatiles along Z", _nMegaY , 1 , SegmentationParameter::NoUnit, true);
-//      _nMegaY is given by Ecal_n_wafers_per_tower from "model_parameters_ILD_o3_v05.xml"
-
-      registerParameter("layer_nCells", "division of megatile into square tiles", _nCells, 36 , SegmentationParameter::NoUnit, true);
-      registerParameter("layer_nStripsX", "division of megatile into strips (1)", _nStripsX, 4 , SegmentationParameter::NoUnit, true);
-      registerParameter("layer_nStripsY", "division of megatile into strips (2)", _nStripsY, 36 , SegmentationParameter::NoUnit, true);
-
-      registerParameter("deadWidth", "width of dead region at edge of megatile", _deadWidth, 0., SegmentationParameter::LengthUnit, true);
+      _description = "Cartesian segmentation in the local XY-plane: megatiles, containing integer number of tiles/strips/cells";
 
       registerIdentifier("identifier_x", "Cell ID identifier for X", _xId, "cellX");
       registerIdentifier("identifier_y", "Cell ID identifier for Y", _yId, "cellY");
@@ -90,9 +35,10 @@ namespace DD4hep {
       registerParameter("identifier_module", "Cell encoding identifier for module", _identifierModule, std::string("module"),
                         SegmentationParameter::NoUnit, true);
 
-      registerParameter("layer_configuration", "layer configuration (S, T, L)", _layerConfig, std::string("TLS"), SegmentationParameter::NoUnit, true);
-
-      _calculated=false;
+      for (int i=0; i<MAX_LAYERS; i++) {
+	_nCellsX[i]=0;
+	_nCellsY[i]=0;
+      }
 
     }
 
@@ -101,291 +47,116 @@ namespace DD4hep {
 
     }
 
-
-    void MegatileLayerGridXY::calculateDerivedQuantities() const {
-
-			if (_calculated) return; // already up to date
-
-      _megaSize  =   _totalSizeY / _nMegaY ; // size of square megatile in this layer
-      _gridSizeL =  (_megaSize - 2.*_deadWidth) / _nStripsX ; // size of cell within megatile (x)
-      _gridSizeT =  (_megaSize - 2.*_deadWidth) / _nStripsY ;
-      _gridSizeS = (_megaSize - 2.*_deadWidth) / _nCells;
-
-      std::cout << "calculateDerivedQuantities: total X width: " << _totalSizeX 
-				<< " / nMega: " << _nMegaY << " = " << _megaSize 
-				<< " ; dead width = " << _deadWidth << " ; nstripsX,Y, nCells : " << _nStripsX 
-				<< " " << _nStripsY << " " << _nCells 
-				<< " ; stripX,Y/cell size: " << _gridSizeL 
-				<< " " << _gridSizeT << " " << _gridSizeS << std::endl;
-
-      _calculated=true;
-
-      return;
-    }
-
     /// determine the position based on the cell ID
     Vector3D MegatileLayerGridXY::position(const CellID& cID) const {
-
-      calculateDerivedQuantities();
+      // this is local position within the megatile
 
       _decoder->setValue(cID);
-      unsigned int LayerIndex;
-      unsigned int WaferIndex;
-      unsigned int CellIndexX;
-      unsigned int CellIndexY;
-			unsigned int ModuleIndex; 
-
-      Vector3D cellPosition;
-
-      LayerIndex = (*_decoder)[_identifierLayer];
-      WaferIndex = (*_decoder)[_identifierWafer];
-      CellIndexX = (*_decoder)[_xId];
-      CellIndexY = (*_decoder)[_yId];
-			ModuleIndex = (*_decoder)[_identifierModule];
-			
-			double regulatingEBUSizeX = _waferOffsetX[LayerIndex][WaferIndex] * 2;
-			double regulatingEBUSizeY = _waferOffsetY[LayerIndex][WaferIndex] * 2;
-
-      std::vector<double> celldims = layercellDimensions( LayerIndex );
-			cellPosition.X = -regulatingEBUSizeX / 2.;
-			cellPosition.Y = -regulatingEBUSizeY / 2.;   // _megaSize includes two _deadWidth.
-      cellPosition.Y += (CellIndexY + 0.5) * celldims[1];
-			double XfromEdge = 0;
-
-#if VERBOSEcoterra
-	if ( cellPosition.X > 100000. ) std::cout << "HERE too large cellPosition.X place 1" << std::endl;
-	if ( _isRegulatingEBU[LayerIndex][WaferIndex] ) {
-		std::cout << "regulatingEBUSizeX, regulatingEBUSizeY(layer) = " << regulatingEBUSizeX << ", " << regulatingEBUSizeY 
-		<< "(" << LayerIndex << ")" << std::endl;
-		 std::cout << "CellIndexX, Y = " << CellIndexX << ", " << CellIndexY << ", celldims[0], [1] = " << celldims[0] 
-				<< ", " << celldims[1] << std::endl;
-		std::cout << "cellPosition.X, Y = " << cellPosition.X << ", " << cellPosition.Y << std::endl;
-	}
-#endif
+      unsigned int layerIndex = (*_decoder)[_identifierLayer];
+      unsigned int waferIndex = (*_decoder)[_identifierWafer];
+      int cellIndexX = (*_decoder)[_xId];
+      int cellIndexY = (*_decoder)[_yId];
 
+      // segmentation info for this megatile ("wafer")
+      getSegInfo(layerIndex, waferIndex);
 
-			int n_cellX = regulatingEBUSizeX / celldims[0];
-			if  (  _isRegulatingEBU[LayerIndex][WaferIndex] && celldims[0] > celldims[1]*1.5 
-						&&CellIndexX == (unsigned int)(n_cellX - 1) ) {
-					XfromEdge = ( n_cellX -1) * celldims[0] + ( regulatingEBUSizeX-( n_cellX -1) * celldims[0] )/2;
-#if VERBOSEcoterra
-						std::cout << "This hit is out of regulation EBU=> merged into neighbour strip; XfromEdge = " << XfromEdge << std::endl;
-#endif
-			} else {
-				 XfromEdge = (CellIndexX + 0.5) * celldims[0];
-			}
-	
-      ////////////////cellPosition.X += _deadWidth + XfromEdge;
-      cellPosition.X += XfromEdge;
-
-#if VERBOSEcoterra
-	std::cout << "XfromEdge = " << XfromEdge << ", CellIndexX = " << CellIndexX << std::endl;
-  if ( cellPosition.X > 100000. ) std::cout << "HERE too large cellPosition.X place 2" << std::endl;
-#endif
-	
-#if VERBOSEcoterra 
-		if ( cellPosition.X  > 50000 || cellPosition.Y  > 50000 || cellPosition.Z  > 50000 
-				|| cellPosition.X  < -50000 || cellPosition.Y < -50000 || cellPosition.Z < -50000 ) {
-			std::cout << "+ _deadWidth + ( CellIndexX + 0.5 )*celldims[0]: cellPosition.X = " << cellPosition.X 
-							<< ", cellPosition.Y = " << cellPosition.Y << std::endl;
-			std::cout << "CellIndexX * celldims[0] = " << CellIndexX << " x " << celldims[0] 
-							<< ", CellIndexY * celldims[1] = " << CellIndexY << " x " << celldims[1] 
-							<< ", _isRegulatingEBU[" << LayerIndex << "][" << WaferIndex << "] = " << _isRegulatingEBU[LayerIndex][WaferIndex]
-							<< std::endl;
-			std::cout << "in MegatileLayerGridXY: regulatingEBUSizeX = " << regulatingEBUSizeX << std::endl;
-			std::cout << "x,y,z = " << cellPosition.X << ", " << cellPosition.Y  << ", " << cellPosition.Z << std::endl;
-	 		std::cout << "LayerIndex = " << LayerIndex << ", WaferIndex = " << WaferIndex << std::endl;
-		}
-#endif
+      Vector3D cellPosition(0,0,0);
+      cellPosition.X = cellIndexX * (_currentSegInfo.megaTileSizeX / _currentSegInfo.nCellsX ) + _currentSegInfo.megaTileOffsetX;
+      cellPosition.Y = cellIndexY * (_currentSegInfo.megaTileSizeY / _currentSegInfo.nCellsY ) + _currentSegInfo.megaTileOffsetY;
 
-/// for Endcaps coterra g728
-			if ( ModuleIndex > 4 ) {
-				double endcapx = cellPosition.Y;
-				cellPosition.Y = cellPosition.X;
-				cellPosition.X = endcapx;
-#if VERBOSEcoterra
-				if ( 192 < WaferIndex && WaferIndex < 197 ) {
-					std::cout << "including BAD POSITION"<< std::endl;
-					std::cout << "x,y,z = " << cellPosition.Y << ", " << cellPosition.X << ", " << cellPosition.Z << std::endl;
-					std::cout << "LayerIndex = " << LayerIndex << ", WaferIndex = " << WaferIndex << ", CellIndexX = " << CellIndexX
-										<< ", CellIndexY = " << CellIndexY << ", ModuleIndex = " << ModuleIndex 
-										<< ", _isRegulatingEBU[LayerIndex][WaferIndex] = " << _isRegulatingEBU[LayerIndex][WaferIndex]  
-										<< std::endl;
-				}
-#endif
-			}
+      if ( abs( cellPosition.X )>10000 || abs( cellPosition.Y )>10000 ) {
+	std::cout << "crazy cell position: " << cellPosition.X << " " << cellPosition.Y << std::endl;
+	std::cout << "layer, wafer, cellx,y indices: " << layerIndex << " " << waferIndex << " " << cellIndexX << " " << cellIndexY << std::endl;
+	assert(0 && "crazy cell position?");
+      }
 
       return cellPosition;
     }
 
+
     /// determine the cell ID based on the position
     CellID MegatileLayerGridXY::cellID(const Vector3D& localPosition, const Vector3D& /* globalPosition */, const VolumeID& vID) const {
-      calculateDerivedQuantities();
+      // this is the local position within a megatile, local coordinates
+
+      // get the layer, wafer, module indices from the volumeID
       _decoder->setValue(vID);
       unsigned int layerIndex = (*_decoder)[_identifierLayer];
-			unsigned int waferIndex = (*_decoder)[_identifierWafer];
-			unsigned int moduleIndex = (*_decoder)[_identifierModule];
-
-      double _x = localPosition.X; 
-      double _y = localPosition.Y;
-
-//// coterra g728.1857 for Endcaps
-			if ( moduleIndex > 4 ) {
-				_x = localPosition.Y;
-				_y = localPosition.X;
-			}
-////
-
-      int _cellIndexX(-1);
-      int _cellIndexY(-1);
-
-	//  std::cout << "hello3 " << _x << " " << _y << std::endl;
-			std::vector<double> celldims = layercellDimensions( layerIndex );
-
-#if VERBOSEcoterra
-			std::cout << "in cellID before area check" << std::endl;
-			std::cout << "_x = localPosition.X: " << _x << " = " << localPosition.X << ", _waferOffsetX[" 
-								<< layerIndex << "][" << waferIndex << "] = " << _waferOffsetX[layerIndex][waferIndex] << std::endl;
-			std::cout << "_y = localPosition.Y: " << _y << " = " << localPosition.Y << ", _waferOffsetY["
-								<< layerIndex << "][" << waferIndex << "] = " << _waferOffsetY[layerIndex][waferIndex] << std::endl;
-			std::cout << "_megaSize = " << _megaSize << std::endl;
-#endif
-			double WaferOffsetX = _waferOffsetX[layerIndex][waferIndex];
-			double WaferOffsetY = _waferOffsetY[layerIndex][waferIndex];
-#if VERBOSEcoterra
-			std::cout << "WaferOffsetX = " << _waferOffsetX[layerIndex][waferIndex] 
-			<< " layer, wafer = " << layerIndex << ", " << waferIndex
-			<< std::endl;
-#endif
-      // _waferOffsetX has negative sign. <<<---??? not true!!
-			if ( _x < -WaferOffsetX || _x > WaferOffsetX ||
-	  	  	 _y < -WaferOffsetY || _y > WaferOffsetY ) {
-#if VERBOSEcoterra
-			std::cout << "This hit falls down in to dead width." << std::endl;
-			std::cout << "localPosition.X/Y = " << localPosition.X << ", " << localPosition.Y << std::endl;
-			std::cout << "layerIndex = " << layerIndex << ", waferIndex = " << waferIndex 
-								<< ", _isRegulatingEBU[" << layerIndex << "][" << waferIndex << "] = " << _isRegulatingEBU[layerIndex][waferIndex] 
-								<< std::endl; 
-#endif
+      unsigned int waferIndex = (*_decoder)[_identifierWafer];
 
-//      	 _x += WaferOffsetX;
-//      	 _y += WaferOffsetY;
+      // segmentation info for this megatile ("wafer")
+      getSegInfo(layerIndex, waferIndex);
 
-				} else {
+      double localX = localPosition.X;
+      double localY = localPosition.Y;
 
-			// localPosition.X is a point in a coordination which center is on the center of wafer.
-      // _waferOffsetX is 1/2 of wafer_dim_x, 
-      // and wafer_dim_x is already subtracted with the guard ring (deadWidth);
-      // _x should be interpreted in a coodination who's origin is on the edge of wafre.
+      // correct for offset : move origin to corner of megatile
+      localX -= _currentSegInfo.megaTileOffsetX;
+      localY -= _currentSegInfo.megaTileOffsetY;
 
-      	 _x += WaferOffsetX;
-      	 _y += WaferOffsetY;
+      // the cell index (counting from the corner)
+      int _cellIndexX = int ( localX / ( _currentSegInfo.megaTileSizeX / _currentSegInfo.nCellsX ) );
+      int _cellIndexY = int ( localY / ( _currentSegInfo.megaTileSizeY / _currentSegInfo.nCellsY ) );
 
-#if VERBOSEcoterra
-  if ( _isRegulatingEBU[layerIndex][waferIndex] ) {
-		if ( localPosition.X < -10000 || localPosition.X > 10000 || localPosition.Y < -10000 || localPosition.Y > 10000 ) {
-			std::cout << "in cellID on hit " << std::endl;
-			std::cout << "_x = localPosition.X + _waferOffsetX; => " << _x << " = " << localPosition.X << " + " << WaferOffsetX << std::endl;
-			std::cout << "_y = localPosition.Y + _waferOffsetY; => " << _y << " = " << localPosition.Y << " + " << WaferOffsetY << std::endl;
-			std::cout << "_megaSize = " << _megaSize << std::endl;
-		}
-	}
-#endif
+      (*_decoder)[_xId] = _cellIndexX;
+      (*_decoder)[_yId] = _cellIndexY;
 
-	  // edge of active area of megatile
-	  			_cellIndexX = int ( _x / celldims[0] );
-	  			_cellIndexY = int ( _y / celldims[1] );
-#if VERBOSEcoterra
-      std::cout << "in cellID after subtract other  wafers" << std::endl;
-			std::cout << "_x = " << _x << " and _y = " << _y << " _deadWidth = " << _deadWidth << "; in cellID" << std::endl;
-				if ( _isRegulatingEBU[layerIndex][waferIndex] && celldims[0] < celldims[1]*1.5 ) {
-					int temp_n_cellX = WaferOffsetX * 2 / celldims[0];
-			//		if ( _cellIndexX > temp_n_cellX - 1 ) {
-						std::cout << "EXESS hit in transeverse layer WaferOffsetX = " << WaferOffsetX 
-						<< ", celldims[0] = " << celldims[0]
-						<< ", _cellIndexX = " << _cellIndexX
-						<< ", _x = " << _x
-						<< ", _megaSize = " << _megaSize
-						<< ", temp_n_cellX = " << temp_n_cellX
-						<< std::endl;
-			//		}
-				}
-				//if ( _x > 10000 || _x < -10000 ) {
-					std::cout << "HERE checking _x and celldims " << std::endl;
-					std::cout << " _cellIndexX, Y = " << _cellIndexX << ", " << _cellIndexY << std::endl;
-					std::cout << " celldims[0], [1] = " << celldims[0] << ", " << celldims[1] << std::endl;
-				//}
-#endif 
-					if ( _isRegulatingEBU[layerIndex][waferIndex] && celldims[0] > celldims[1]*1.5 ) {
-						int n_cellX = WaferOffsetX * 2 / celldims[0]; // How many strip-row can stay on this EBU.
-						if ( _cellIndexX > n_cellX - 1 ) { // hit position is out of this EBU.
-							_cellIndexX -= 1;
-#if VERBOSEcoterra
-						std::cout << "This hit is out of regulation EBU=> merged into neighbour strip; _cellIndexX = " << _cellIndexX << std::endl;
-#endif
-						}
-					}
-				}
-	if ( _isRegulatingEBU[layerIndex][waferIndex] ) {
-      	 //_cellIndexX = 0;
-      	 //_cellIndexY = 0;
-	}
-
-
-      	(*_decoder)[_xId] = _cellIndexX;
-      	(*_decoder)[_yId] = _cellIndexY;
-     
-#if VERBOSEcoterra
-				int CellIndexX = (*_decoder)[_xId]; 
-				if ( CellIndexX > 10 ) std::cout << "HERE check just before return _decoder: _cellIndexX = " << _cellIndexX << std::endl;
-#endif
-      	return _decoder->getValue();
+      return _decoder->getValue();
     }
 
-    std::vector<double> MegatileLayerGridXY::cellDimensions(const CellID& cID) const {
 
+    std::vector<double> MegatileLayerGridXY::cellDimensions(const CellID& cID) const {
       _decoder->setValue( cID );
-      unsigned int _layerIndex = (*_decoder)[_identifierLayer];
-
-      return layercellDimensions(_layerIndex);
+      unsigned int layerIndex = (*_decoder)[_identifierLayer];
+      unsigned int waferIndex = (*_decoder)[_identifierWafer];
+      return cellDimensions(layerIndex, waferIndex);
+    }
 
+    void MegatileLayerGridXY::setSpecialMegaTile( unsigned int layer, unsigned int tile, 
+						  double sizex, double sizey,
+						  double offsetx, double offsety,
+						  unsigned int ncellsx, unsigned int ncellsy ) {
+
+      std::pair <int, int> tileid(layer, tile);
+      segInfo sinf;
+      sinf.megaTileSizeX = sizex;
+      sinf.megaTileSizeY = sizey;
+      sinf.megaTileOffsetX = offsetx;
+      sinf.megaTileOffsetY = offsety;
+      sinf.nCellsX = ncellsx;
+      sinf.nCellsY = ncellsy;
+      specialMegaTiles_layerWafer[tileid] = sinf;
     }
 
-    std::vector<double> MegatileLayerGridXY::layercellDimensions(const int layerIndex) const {
 
-      calculateDerivedQuantities();
+    void MegatileLayerGridXY::getSegInfo( unsigned int layerIndex, unsigned int waferIndex) const {
 
-      unsigned int mm = layerIndex % _layerConfig.length();
-      
-      const char laytype = _layerConfig.at(mm);
+      assert ( layerIndex < MAX_LAYERS && "layer index too high" );
 
-      double xsize(0);
-      double ysize(0);
+      std::pair < unsigned int, unsigned int > tileid(layerIndex, waferIndex);
+      if ( specialMegaTiles_layerWafer.find( tileid ) == specialMegaTiles_layerWafer.end() ) { // standard megatile
+	_currentSegInfo.megaTileSizeX   = _megaTileSizeX;
+        _currentSegInfo.megaTileSizeY   = _megaTileSizeY;
+        _currentSegInfo.megaTileOffsetX = _megaTileOffsetX;
+        _currentSegInfo.megaTileOffsetY = _megaTileOffsetY;
+        _currentSegInfo.nCellsX         = _nCellsX[layerIndex];
+        _currentSegInfo.nCellsY         = _nCellsY[layerIndex];
+      } else { // special megatile
+	_currentSegInfo = specialMegaTiles_layerWafer.find( tileid )->second;
+      }	
+    }
 
-      if ( laytype=='S' ) { // square cell
-				xsize = ysize = _gridSizeS;
-      } else if ( laytype=='T' ) { // transverse
-				xsize = _gridSizeT;
-				ysize = _gridSizeL;
-      } else if ( laytype=='L' ) { // transverse
-				xsize = _gridSizeL;
-				ysize = _gridSizeT;
-      } else {
-				assert(0);
-      }
-      
-      //      std::cout << laytype << " " << xsize << " " << ysize << std::endl;
+    std::vector<double> MegatileLayerGridXY::cellDimensions(const unsigned int layerIndex, const unsigned int waferIndex) const {
+      // calculate the cell size for a given wafer in a given layer
+
+      getSegInfo(layerIndex, waferIndex);
+
+      double xsize = _currentSegInfo.megaTileSizeX/_currentSegInfo.nCellsX;
+      double ysize = _currentSegInfo.megaTileSizeY/_currentSegInfo.nCellsY;
 
 #if __cplusplus >= 201103L
-      //      return {_gridSizeX[_layerIndex], _gridSizeY[_layerIndex]};
-      //      return {_gridSizeX, _gridSizeY};
       return {xsize, ysize};
 #else
       std::vector<double> cellDims(2,0.0);
-      //cellDims[0] = _gridSizeX[_layerIndex];
-      //cellDims[1] = _gridSizeY[_layerIndex];
-      //cellDims[0] = _gridSizeX;
-      //cellDims[1] = _gridSizeY;
       cellDims[0] = xsize;
       cellDims[1] = ysize;
       return cellDims;