diff --git a/Documentation/Examples/stack_example.cc b/Documentation/Examples/stack_example.cc
index 6f6dd03824da770253d2a2557e2c467c0ce7752c..3d28e5566b2f53df940bf2769e85c4331d276e87 100644
--- a/Documentation/Examples/stack_example.cc
+++ b/Documentation/Examples/stack_example.cc
@@ -21,6 +21,7 @@
 
 using namespace corsika::units::si;
 using namespace corsika::stack;
+using namespace corsika;
 using namespace std;
 
 void fill(corsika::stack::super_stupid::SuperStupidStack& s) {
@@ -44,7 +45,6 @@ void read(corsika::stack::super_stupid::SuperStupidStack& s) {
     assert(p.GetPID() == corsika::particles::Code::Electron);
     assert(p.GetEnergy() == 1.5_GeV * (i++));
   }
-  // assert(total_energy == 82.5_GeV);
 }
 
 int main() {
diff --git a/Framework/Geometry/BaseTrajectory.h b/Framework/Geometry/BaseTrajectory.h
index 8d13f97a3a96c46549f8a4b8a5baad8951e33353..d8cf7ff5633de7847db2ede94634346048af08bb 100644
--- a/Framework/Geometry/BaseTrajectory.h
+++ b/Framework/Geometry/BaseTrajectory.h
@@ -1,4 +1,14 @@
 
+/*
+ * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+
 /*
  * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
  *
diff --git a/Framework/Geometry/BaseVector.h b/Framework/Geometry/BaseVector.h
index c41d8edf205f84b2066f83a61a253130225ab207..106fb632446f476913bb93b01e085912182ebb1e 100644
--- a/Framework/Geometry/BaseVector.h
+++ b/Framework/Geometry/BaseVector.h
@@ -1,4 +1,14 @@
 
+/*
+ * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+
 /*
  * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
  *
diff --git a/Framework/Geometry/CoordinateSystem.cc b/Framework/Geometry/CoordinateSystem.cc
index 35e2941678d970b1868626a252785055d4413be4..6f818462b0399f7026864d7de3fc8da6aba877c9 100644
--- a/Framework/Geometry/CoordinateSystem.cc
+++ b/Framework/Geometry/CoordinateSystem.cc
@@ -1,4 +1,14 @@
 
+/*
+ * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+
 /*
  * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
  *
diff --git a/Framework/Geometry/FourVector.h b/Framework/Geometry/FourVector.h
index 0a917a153d23a623f8d34881300a2ebccdc3ed6d..1cd15c6afc2f0fa4e21bb1168f2fdb8211078b25 100644
--- a/Framework/Geometry/FourVector.h
+++ b/Framework/Geometry/FourVector.h
@@ -1,4 +1,14 @@
 
+/*
+ * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+
 /*
  * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
  *
diff --git a/Framework/Geometry/Helix.h b/Framework/Geometry/Helix.h
index caf2c34d0009bb2e33ba0bcbf0babcdf7136ba41..4d9badf8c757721507613a5abb1c4257dbe3ae0b 100644
--- a/Framework/Geometry/Helix.h
+++ b/Framework/Geometry/Helix.h
@@ -1,4 +1,14 @@
 
+/*
+ * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+
 /*
  * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
  *
diff --git a/Framework/Geometry/Line.h b/Framework/Geometry/Line.h
index ce795d8d9ffc918c9b6559278ee6a8cf58bdd282..8672aa178ecbc2bcc9dffd199d2348db0f9e86d4 100644
--- a/Framework/Geometry/Line.h
+++ b/Framework/Geometry/Line.h
@@ -1,4 +1,14 @@
 
+/*
+ * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+
 /*
  * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
  *
diff --git a/Framework/Geometry/RootCoordinateSystem.h b/Framework/Geometry/RootCoordinateSystem.h
index 0e641b3a349e22173fd894e7bd005a8f71530a2b..b5f496aa48e2e5ee0c89aba0c84e39b4df1fbd67 100644
--- a/Framework/Geometry/RootCoordinateSystem.h
+++ b/Framework/Geometry/RootCoordinateSystem.h
@@ -1,4 +1,14 @@
 
+/*
+ * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+
 /*
  * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
  *
diff --git a/Framework/Geometry/Sphere.h b/Framework/Geometry/Sphere.h
index 5a38820a3fecc600d6c2b55177ababc777b7c13c..01d74543ec4407dc5a3acc12edee73c92c8d1504 100644
--- a/Framework/Geometry/Sphere.h
+++ b/Framework/Geometry/Sphere.h
@@ -1,4 +1,14 @@
 
+/*
+ * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+
 /*
  * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
  *
diff --git a/Framework/Geometry/Trajectory.h b/Framework/Geometry/Trajectory.h
index 5be99e363ebe2b734510e36b2b8a21a5d86ecb38..69bcfd6f8ad6a8ea4314531c48865e06179bbe91 100644
--- a/Framework/Geometry/Trajectory.h
+++ b/Framework/Geometry/Trajectory.h
@@ -1,4 +1,14 @@
 
+/*
+ * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+
 /*
  * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
  *
diff --git a/Framework/Geometry/Volume.h b/Framework/Geometry/Volume.h
index c6c953d9511132e8663a25a3f05cbbd82b98916a..95306edca53ff546b54fe362babbdcf17a73334c 100644
--- a/Framework/Geometry/Volume.h
+++ b/Framework/Geometry/Volume.h
@@ -1,4 +1,14 @@
 
+/*
+ * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+
 /*
  * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
  *
diff --git a/Framework/Geometry/testFourVector.cc b/Framework/Geometry/testFourVector.cc
index a4063f0ef9a9b7d3528142c2d06094e9e70c7755..824426aad0bdaa4173fff93b26486d1eb423fea6 100644
--- a/Framework/Geometry/testFourVector.cc
+++ b/Framework/Geometry/testFourVector.cc
@@ -1,4 +1,14 @@
 
+/*
+ * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+
 /*
  * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
  *
diff --git a/Framework/Particles/NuclearData.xml b/Framework/Particles/NuclearData.xml
index a1c35c85eafa071fefa62ad21b5b6c4adfeea82f..0fff50c33839093fc8bfacb9a7632f36110a2492 100644
--- a/Framework/Particles/NuclearData.xml
+++ b/Framework/Particles/NuclearData.xml
@@ -1,5 +1,8 @@
 <chapter name="Nuclear Data">
 
+<particle id="1000000000" name="nucleus" A="0" Z="0" >
+</particle> 
+
 <particle id="1000010010" name="hydrogen" A="1" Z="1" >
 </particle> 
 
diff --git a/Framework/Particles/ParticleData.xml b/Framework/Particles/ParticleData.xml
index 16e7878abe42147c0d2ec3a285465c975ecec93c..6c8c3d0a06fec27481d508c030ddad6a7df8f925 100644
--- a/Framework/Particles/ParticleData.xml
+++ b/Framework/Particles/ParticleData.xml
@@ -4,7 +4,7 @@
 <particle id="0" name="void" spinType="0" chargeType="0" colType="0" 
           m0="0.00000"> 
 </particle> 
- 
+
 <!--
 <particle id="1" name="d" antiName="dbar" spinType="2" chargeType="-1" colType="1" 
           m0="0.33000"> 
diff --git a/Framework/Particles/ParticleProperties.cc b/Framework/Particles/ParticleProperties.cc
index 7df29e065beb35100afddb77a1e22ff70480eb01..4932492ddf33028b456670e9ca2da0a5b8d0efe2 100644
--- a/Framework/Particles/ParticleProperties.cc
+++ b/Framework/Particles/ParticleProperties.cc
@@ -1,4 +1,14 @@
 
+/*
+ * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+
 /*
  * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
  *
diff --git a/Framework/Particles/pdxml_reader.py b/Framework/Particles/pdxml_reader.py
index 2e463a0e1d7c9a227d9f49160d4916ecf3390761..99544cf5029c933ef9a43478dae73d28c2efbb56 100755
--- a/Framework/Particles/pdxml_reader.py
+++ b/Framework/Particles/pdxml_reader.py
@@ -251,7 +251,7 @@ def read_nuclei_db(filename, particle_db, classnames):
         }
     
     return particle_db
-    
+
 
 
 
diff --git a/Framework/Particles/testParticles.cc b/Framework/Particles/testParticles.cc
index 30cbbca868befbb4ae3ea6c5db558ab1038d5547..1ca8b2188b5004c22c3168b5ad90290aa10a8ffb 100644
--- a/Framework/Particles/testParticles.cc
+++ b/Framework/Particles/testParticles.cc
@@ -45,6 +45,8 @@ TEST_CASE("ParticleProperties", "[Particles]") {
   SECTION("Names") {
     REQUIRE(Electron::GetName() == "e-");
     REQUIRE(PiMinus::GetName() == "pi-");
+    REQUIRE(Nucleus::GetName() == "nucleus");
+    REQUIRE(Iron::GetName() == "iron");
   }
 
   SECTION("PDG") {
diff --git a/Framework/StackInterface/ParticleBase.h b/Framework/StackInterface/ParticleBase.h
index 5bf2ece2b1e26af2ccfad2b34e5c3a9534250d86..59b133a96dbf81071b1c821dc79a0c5da2d59dbd 100644
--- a/Framework/StackInterface/ParticleBase.h
+++ b/Framework/StackInterface/ParticleBase.h
@@ -22,6 +22,26 @@ namespace corsika::stack {
    The base class to define the readout of particle properties from a
    particle stack. Every stack must implement this readout via the
    ParticleBase class.
+
+   The StackIterator template argument is derived from StackIteratorInterface, which is of
+   type <code> template <typename StackData, template <typename> typename
+   ParticleInterface> class StackIteratorInterface : public
+   ParticleInterface<StackIteratorInterface<StackData, ParticleInterface>>
+   </code>
+
+   where StackData must refer to a Stack type, and
+   ParticleInterface<StackIteratorInterface> is the corresponding particle readout class.
+
+   Thus, StackIteratorInterface is a CRTP class, injecting the full StackIteratorInterface
+   machinery into the ParticleInterface (aka ParticleBase) type!
+
+   The declartion of a StackIteratorInterface type simultaneously declares the
+   corresponding ParticleInterface type.
+
+   Furthermore, the operator* of the StackIteratorInterface returns a
+   static_cast to the ParticleInterface type, allowing a direct
+   readout of the particle data from the iterator.
+
   */
 
   template <typename StackIterator>
@@ -31,6 +51,7 @@ namespace corsika::stack {
     ParticleBase() = default;
 
   private:
+    /*
     // those copy constructors and assigments should never be implemented
     ParticleBase(ParticleBase&) = delete;
     ParticleBase operator=(ParticleBase&) = delete;
@@ -38,10 +59,11 @@ namespace corsika::stack {
     ParticleBase operator=(ParticleBase&&) = delete;
     ParticleBase(const ParticleBase&) = delete;
     ParticleBase operator=(const ParticleBase&) = delete;
-
+    */
   public:
-    /** delete this particle on the stack. The corresponding iterator
-     *  will be invalidated by this operation
+    /**
+     * Delete this particle on the stack. The corresponding iterator
+     * will be invalidated by this operation
      */
     void Delete() { GetIterator().GetStack().Delete(GetIterator()); }
 
@@ -55,7 +77,7 @@ namespace corsika::stack {
       return GetStack().AddSecondary(GetIterator(), args...);
     }
 
-    //  protected: // todo should [MAY]be proteced, but don't now how to 'friend Stack'
+    // protected: // todo should [MAY]be proteced, but don't now how to 'friend Stack'
     // Function to provide CRTP access to inheriting class (type)
     /**
      * return the corresponding StackIterator for this particle
@@ -66,7 +88,8 @@ namespace corsika::stack {
     }
 
   protected:
-    /** @name Access to underlying stack data
+    /**
+        @name Access to underlying stack data
         @{
     */
     auto& GetStackData() { return GetIterator().GetStackData(); }
@@ -78,7 +101,7 @@ namespace corsika::stack {
     /**
      * return the index number of the underlying iterator object
      */
-    int GetIndex() const { return GetIterator().GetIndex(); }
+    unsigned int GetIndex() const { return GetIterator().GetIndex(); }
   };
 
 } // namespace corsika::stack
diff --git a/Framework/StackInterface/Stack.h b/Framework/StackInterface/Stack.h
index bc9aa2e5c0ee0388df07a6ea4a90b4c3865c4485..ff92fa2383e2614b2f6fa4438d6ca492d77608e6 100644
--- a/Framework/StackInterface/Stack.h
+++ b/Framework/StackInterface/Stack.h
@@ -27,36 +27,67 @@ namespace corsika::stack {
      ParticleInterface, which is one of the essential template
      parameters for the Stack.
 
-     Important: ParticleInterface must inherit from ParticleBase !
+     <b>Important:</b> ParticleInterface must inherit from ParticleBase !
    */
 
   template <typename>
   class ParticleInterface; // forward decl
 
   /**
-     Interface definition of a Stack object. The Stack implements the
+     The Stack class provides (and connects) the main particle data storage machinery.
+
+     The StackData type is the user-provided bare data storage
+     object. This can be of any complexity, from a simple struct
+     (fortran common block), to a combination of different and
+     distributed data sources.
+
+     The user-provided ParticleInterface template type is the base
+     class type of the StackIteratorInterface class (CRTP) and must
+     provide all functions to read single particle data from the
+     StackData, given an 'unsigned int' index.
+
+     The Stack implements the
      std-type begin/end function to allow integration in normal for
-     loops etc.
+     loops, ranges, etc.
    */
 
   template <typename StackData, template <typename> typename ParticleInterface>
   class Stack : public StackData {
 
   public:
-    typedef Stack<StackData, ParticleInterface> StackType;
+    typedef StackData StackImpl; ///< this is the type of the user-provided data structure
+    template <typename SI>
+    using PIType = ParticleInterface<SI>;
+    // typedef ParticleInterface<StackIteratorInterface> StackParticleInterface;  ///<
+    // this is the type of the user-provided ParticleInterface typedef Stack<StackData,
+    // ParticleInterface> StackType;
+
+    /**
+     * Via the StackIteratorInterface and ConstStackIteratorInterface
+     * specialization, the type of the StackIterator
+     * template class is declared for a particular stack data
+     * object. Using CRTP, this also determines the type of
+     * ParticleInterface template class simultaneously.
+     */
     typedef StackIteratorInterface<StackData, ParticleInterface> StackIterator;
     typedef ConstStackIteratorInterface<StackData, ParticleInterface> ConstStackIterator;
-    // typedef const StackIterator ConstStackIterator;
+    /**
+     * this is the full type of the declared ParticleInterface: typedef typename
+     */
     typedef typename StackIterator::ParticleInterfaceType ParticleType;
+
     friend class StackIteratorInterface<StackData, ParticleInterface>;
     friend class ConstStackIteratorInterface<StackData, ParticleInterface>;
 
+  protected:
+    using StackData::Copy;
+    using StackData::Swap;
+
   public:
     using StackData::GetCapacity;
     using StackData::GetSize;
 
     using StackData::Clear;
-    using StackData::Copy;
 
     using StackData::DecrementSize;
     using StackData::IncrementSize;
@@ -88,9 +119,14 @@ namespace corsika::stack {
       IncrementSize();
       return StackIterator(*this, GetSize() - 1, parent, v...);
     }
-    void Copy(StackIterator& a, StackIterator& b) { Copy(a.GetIndex(), b.GetIndex()); }
+    void Swap(StackIterator a, StackIterator b) { Swap(a.GetIndex(), b.GetIndex()); }
+    void Swap(ConstStackIterator a, ConstStackIterator b) {
+      Swap(a.GetIndex(), b.GetIndex());
+    }
+    void Copy(StackIterator a, StackIterator b) { Copy(a.GetIndex(), b.GetIndex()); }
+    void Copy(ConstStackIterator a, StackIterator b) { Copy(a.GetIndex(), b.GetIndex()); }
     /// delete this particle
-    void Delete(StackIterator& p) {
+    void Delete(StackIterator p) {
       if (GetSize() == 0) { /*error*/
         throw std::runtime_error("Stack, cannot delete entry since size is zero");
       }
@@ -98,7 +134,7 @@ namespace corsika::stack {
       DeleteLast();
       // p.SetInvalid();
     }
-    void Delete(ParticleType& p) { Delete(p.GetIterator()); }
+    void Delete(ParticleType p) { Delete(p.GetIterator()); }
     /// delete last particle on stack by decrementing stack size
     void DeleteLast() { DecrementSize(); }
     /// check if there are no further particles on stack
diff --git a/Framework/StackInterface/StackIteratorInterface.h b/Framework/StackInterface/StackIteratorInterface.h
index ad00c6c0939a4e99c020164ef4cd85cbc81a8788..74f1a7802e1d55c0739686496d0ba568557785ac 100644
--- a/Framework/StackInterface/StackIteratorInterface.h
+++ b/Framework/StackInterface/StackIteratorInterface.h
@@ -51,11 +51,13 @@ namespace corsika::stack {
   class StackIteratorInterface
       : public ParticleInterface<StackIteratorInterface<StackData, ParticleInterface>> {
 
+  public:
     typedef Stack<StackData, ParticleInterface> StackType;
     /*typedef
         typename std::conditional<std::is_const<StackData>::value,
                                   const Stack<const StackData, ParticleInterface>&,
                                   Stack<StackData, ParticleInterface>&>::type StackType;*/
+
     typedef ParticleInterface<StackIteratorInterface<StackData, ParticleInterface>>
         ParticleInterfaceType;
 
@@ -63,7 +65,7 @@ namespace corsika::stack {
     friend class ParticleBase<StackIteratorInterface>; // for access to GetStackData
 
   private:
-    int fIndex = 0;
+    unsigned int fIndex = 0;
     StackType* fData = 0; // info: Particles and StackIterators become invalid when parent
                           // Stack is copied or deleted!
 
@@ -75,7 +77,7 @@ namespace corsika::stack {
         @param data reference to the stack [rw]
         @param index index on stack
      */
-    StackIteratorInterface(StackType& data, const int index)
+    StackIteratorInterface(StackType& data, const unsigned int index)
         : fIndex(index)
         , fData(&data) {}
 
@@ -87,7 +89,7 @@ namespace corsika::stack {
        ParticleInterfaceType::SetParticleData(...) function
      */
     template <typename... Args>
-    StackIteratorInterface(StackType& data, const int index, const Args... args)
+    StackIteratorInterface(StackType& data, const unsigned int index, const Args... args)
         : fIndex(index)
         , fData(&data) {
       (**this).SetParticleData(args...);
@@ -104,7 +106,7 @@ namespace corsika::stack {
        ParticleInterfaceType::SetParticleData(...) function
     */
     template <typename... Args>
-    StackIteratorInterface(StackType& data, const int index,
+    StackIteratorInterface(StackType& data, const unsigned int index,
                            StackIteratorInterface& parent, const Args... args)
         : fIndex(index)
         , fData(&data) {
@@ -124,6 +126,9 @@ namespace corsika::stack {
       ++fIndex;
       return tmp;
     }
+    StackIteratorInterface operator+(int delta) {
+      return StackIteratorInterface(*fData, fIndex + delta);
+    }
     bool operator==(const StackIteratorInterface& rhs) { return fIndex == rhs.fIndex; }
     bool operator!=(const StackIteratorInterface& rhs) { return fIndex != rhs.fIndex; }
     /// Convert to value type
@@ -141,7 +146,7 @@ namespace corsika::stack {
      */
     ///@{
     /// Get current particle index
-    inline int GetIndex() const { return fIndex; }
+    inline unsigned int GetIndex() const { return fIndex; }
     /// Get current particle Stack object
     inline StackType& GetStack() { return *fData; }
     /// Get current particle const Stack object
@@ -164,6 +169,7 @@ namespace corsika::stack {
       : public ParticleInterface<
             ConstStackIteratorInterface<StackData, ParticleInterface>> {
 
+  public:
     typedef Stack<StackData, ParticleInterface> StackType;
     typedef ParticleInterface<ConstStackIteratorInterface<StackData, ParticleInterface>>
         ParticleInterfaceType;
@@ -172,7 +178,7 @@ namespace corsika::stack {
     friend class ParticleBase<ConstStackIteratorInterface>; // for access to GetStackData
 
   private:
-    int fIndex = 0;
+    unsigned int fIndex = 0;
     const StackType* fData = 0; // info: Particles and StackIterators become invalid when
                                 // parent Stack is copied or deleted!
 
@@ -180,7 +186,7 @@ namespace corsika::stack {
     ConstStackIteratorInterface() = delete;
 
   public:
-    ConstStackIteratorInterface(const StackType& data, const int index)
+    ConstStackIteratorInterface(const StackType& data, const unsigned int index)
         : fIndex(index)
         , fData(&data) {}
 
@@ -210,6 +216,9 @@ namespace corsika::stack {
       ++fIndex;
       return tmp;
     }
+    ConstStackIteratorInterface operator+(int delta) {
+      return ConstStackIteratorInterface(*fData, fIndex + delta);
+    }
     bool operator==(const ConstStackIteratorInterface& rhs) {
       return fIndex == rhs.fIndex;
     }
@@ -227,7 +236,7 @@ namespace corsika::stack {
         Only the const versions for read-only access
      */
     ///@{
-    inline int GetIndex() const { return fIndex; }
+    inline unsigned int GetIndex() const { return fIndex; }
     inline const StackType& GetStack() const { return *fData; }
     inline const StackData& GetStackData() const { return fData->GetStackData(); }
     ///@}
diff --git a/Framework/StackInterface/testStackInterface.cc b/Framework/StackInterface/testStackInterface.cc
index cde5a76e2b72dc7563efbd554a3024c4d28fff9f..3a5e945c46e4faa552724ee743241de935f48a41 100644
--- a/Framework/StackInterface/testStackInterface.cc
+++ b/Framework/StackInterface/testStackInterface.cc
@@ -29,8 +29,8 @@ public:
   // these functions are needed for the Stack interface
   void Init() {}
   void Clear() { fData.clear(); }
-  int GetSize() const { return fData.size(); }
-  int GetCapacity() const { return fData.size(); }
+  unsigned int GetSize() const { return fData.size(); }
+  unsigned int GetCapacity() const { return fData.size(); }
   void Copy(const int i1, const int i2) { fData[i2] = fData[i1]; }
   void Swap(const int i1, const int i2) {
     double tmp0 = fData[i1];
@@ -97,8 +97,8 @@ TEST_CASE("Stack", "[Stack]") {
     s.Init();
     s.Clear();
     s.IncrementSize();
-    s.Copy(0, 0);
-    s.Swap(0, 0);
+    s.Copy(s.cbegin(), s.begin());
+    s.Swap(s.begin(), s.begin());
     s.GetCapacity();
     REQUIRE(s.GetSize() == 1);
     s.DecrementSize();
diff --git a/Processes/Sibyll/SibStack.h b/Processes/Sibyll/SibStack.h
index 22b446c57c4c8febbeef1fe9312ea7810bd99e3b..95f7011e9a922eace00bf7e0b7c07e95a5d8b181 100644
--- a/Processes/Sibyll/SibStack.h
+++ b/Processes/Sibyll/SibStack.h
@@ -29,34 +29,34 @@ namespace corsika::process::sibyll {
     void Init();
 
     void Clear() { s_plist_.np = 0; }
-    int GetSize() const { return s_plist_.np; }
-    int GetCapacity() const { return 8000; }
+    unsigned int GetSize() const { return s_plist_.np; }
+    unsigned int GetCapacity() const { return 8000; }
 
-    void SetId(const int i, const int v) { s_plist_.llist[i] = v; }
-    void SetEnergy(const int i, const corsika::units::si::HEPEnergyType v) {
+    void SetId(const unsigned int i, const int v) { s_plist_.llist[i] = v; }
+    void SetEnergy(const unsigned int i, const corsika::units::si::HEPEnergyType v) {
       using namespace corsika::units::si;
       s_plist_.p[3][i] = v / 1_GeV;
     }
-    void SetMass(const int i, const corsika::units::si::HEPMassType v) {
+    void SetMass(const unsigned int i, const corsika::units::si::HEPMassType v) {
       using namespace corsika::units::si;
       s_plist_.p[4][i] = v / 1_GeV;
     }
-    void SetMomentum(const int i, const MomentumVector& v) {
+    void SetMomentum(const unsigned int i, const MomentumVector& v) {
       using namespace corsika::units::si;
       auto tmp = v.GetComponents();
       for (int idx = 0; idx < 3; ++idx) s_plist_.p[idx][i] = tmp[idx] / 1_GeV;
     }
 
-    int GetId(const int i) const { return s_plist_.llist[i]; }
+    int GetId(const unsigned int i) const { return s_plist_.llist[i]; }
     corsika::units::si::HEPEnergyType GetEnergy(const int i) const {
       using namespace corsika::units::si;
       return s_plist_.p[3][i] * 1_GeV;
     }
-    corsika::units::si::HEPEnergyType GetMass(const int i) const {
+    corsika::units::si::HEPEnergyType GetMass(const unsigned int i) const {
       using namespace corsika::units::si;
       return s_plist_.p[4][i] * 1_GeV;
     }
-    MomentumVector GetMomentum(const int i) const {
+    MomentumVector GetMomentum(const unsigned int i) const {
       using corsika::geometry::CoordinateSystem;
       using corsika::geometry::QuantityVector;
       using corsika::geometry::RootCoordinateSystem;
@@ -69,9 +69,15 @@ namespace corsika::process::sibyll {
       return v1;
     }
 
-    void Copy(const int i1, const int i2) {
-      s_plist_.llist[i1] = s_plist_.llist[i2];
-      s_plist_.p[3][i1] = s_plist_.p[3][i2];
+    void Copy(const unsigned int i1, const unsigned int i2) {
+      s_plist_.llist[i2] = s_plist_.llist[i1];
+      for (unsigned int i = 0; i < 5; ++i) s_plist_.p[i][i2] = s_plist_.p[i][i1];
+    }
+
+    void Swap(const unsigned int i1, const unsigned int i2) {
+      std::swap(s_plist_.llist[i1], s_plist_.llist[i2]);
+      for (unsigned int i = 0; i < 5; ++i)
+        std::swap(s_plist_.p[i][i1], s_plist_.p[i][i2]);
     }
 
   protected:
diff --git a/Setup/CMakeLists.txt b/Setup/CMakeLists.txt
index 2f9a4ddf36a57d5f73d523fa5771ee17193a83bd..5b73fa434c99f42a5372c2f7488463a6d9c956de 100644
--- a/Setup/CMakeLists.txt
+++ b/Setup/CMakeLists.txt
@@ -20,6 +20,7 @@ target_link_libraries (
   INTERFACE
   CORSIKAgeometry
   SuperStupidStack
+  NuclearStackExtension
   )
 
 target_include_directories (
diff --git a/Setup/SetupStack.h b/Setup/SetupStack.h
index 0ab05fffbf177b72a21880e5b2b2043a4aba448d..bae7163e259538571b5192338158471865975fb5 100644
--- a/Setup/SetupStack.h
+++ b/Setup/SetupStack.h
@@ -12,11 +12,19 @@
 #ifndef _corsika_setup_setupstack_h_
 #define _corsika_setup_setupstack_h_
 
+#include <corsika/stack/nuclear_extension/NuclearStackExtension.h>
 #include <corsika/stack/super_stupid/SuperStupidStack.h>
 
+// this is an auxiliary help typedef, which I don't know how to put
+// into NuclearStackExtension.h where it belongs...
+template<typename StackIter> using ExtendedParticleInterfaceType =
+  corsika::stack::nuclear_extension::NuclearParticleInterface<corsika::stack::super_stupid::SuperStupidStack::PIType, StackIter>;
+
 namespace corsika::setup {
 
-  typedef corsika::stack::super_stupid::SuperStupidStack Stack;
+  using Stack = corsika::stack::nuclear_extension::NuclearStackExtension<corsika::stack::super_stupid::SuperStupidStack, ExtendedParticleInterfaceType>;
+
+  // typedef corsika::stack::super_stupid::SuperStupidStack Stack;
 }
 
 #endif
diff --git a/Stack/CMakeLists.txt b/Stack/CMakeLists.txt
index 4cf17a6163a1503f61ad50d4bb1658c4d35727c4..a9645130bd9f81cd4415fbfeccc88e5a736b7160 100644
--- a/Stack/CMakeLists.txt
+++ b/Stack/CMakeLists.txt
@@ -1,3 +1,4 @@
 
 add_subdirectory (DummyStack)
 add_subdirectory (SuperStupidStack)
+add_subdirectory (NuclearStackExtension)
diff --git a/Stack/DummyStack/SuperStupidStack.h b/Stack/DummyStack/SuperStupidStack.h
deleted file mode 100644
index 3abbb3139abb243f5782241db711edf852d98424..0000000000000000000000000000000000000000
--- a/Stack/DummyStack/SuperStupidStack.h
+++ /dev/null
@@ -1,108 +0,0 @@
-
-/*
- * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
- *
- * See file AUTHORS for a list of contributors.
- *
- * This software is distributed under the terms of the GNU General Public
- * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
- * the license.
- */
-
-#ifndef _include_superstupidstack_h_
-#define _include_superstupidstack_h_
-
-#include <corsika/particles/ParticleProperties.h>
-#include <corsika/stack/Stack.h>
-#include <corsika/units/PhysicalUnits.h>
-
-#include <string>
-#include <vector>
-
-namespace corsika::stack {
-
-  namespace super_stupid {
-
-    using corsika::particles::Code;
-    using corsika::units::si::EnergyType;
-    using corsika::units::si::operator""_GeV; // literals;
-
-    /**
-     * Example of a particle object on the stack.
-     */
-
-    template <typename _Stack>
-    class ParticleRead : public StackIteratorInfo<_Stack, ParticleRead<_Stack> > {
-
-      using StackIteratorInfo<_Stack, ParticleRead>::GetIndex;
-      using StackIteratorInfo<_Stack, ParticleRead>::GetStack;
-
-    public:
-      void SetId(const Code id) { GetStack().SetId(GetIndex(), id); }
-      void SetEnergy(const EnergyType& e) { GetStack().SetEnergy(GetIndex(), e); }
-
-      Code GetId() const { return GetStack().GetId(GetIndex()); }
-      const EnergyType& GetEnergy() const { return GetStack().GetEnergy(GetIndex()); }
-    };
-
-    /**
-     *
-     * Memory implementation of the most simple (stupid) particle stack object.
-     */
-
-    class SuperStupidStackImpl {
-
-    public:
-      void Init() {}
-
-      void Clear() {
-        fDataE.clear();
-        fDataId.clear();
-      }
-
-      int GetSize() const { return fDataId.size(); }
-      int GetCapacity() const { return fDataId.size(); }
-
-      void SetId(const int i, const Code id) { fDataId[i] = id; }
-      void SetEnergy(const int i, const EnergyType& e) { fDataE[i] = e; }
-
-      const Code GetId(const int i) const { return fDataId[i]; }
-      const EnergyType& GetEnergy(const int i) const { return fDataE[i]; }
-
-      /**
-       *   Function to copy particle at location i2 in stack to i1
-       */
-      void Copy(const int i1, const int i2) {
-        fDataE[i2] = fDataE[i1];
-        fDataId[i2] = fDataId[i1];
-      }
-
-    protected:
-      void IncrementSize() {
-        fDataE.push_back(0_GeV);
-        fDataId.push_back(Code::unknown);
-      }
-      void DecrementSize() {
-        if (fDataE.size() > 0) {
-          fDataE.pop_back();
-          fDataId.pop_back();
-        }
-      }
-
-    private:
-      /// the actual memory to store particle data
-
-      std::vector<Code> fDataId;
-      std::vector<EnergyType> fDataE;
-
-    }; // end class SuperStupidStackImpl
-
-    typedef StackIterator<SuperStupidStackImpl, ParticleRead<SuperStupidStackImpl> >
-        Particle;
-    typedef Stack<SuperStupidStackImpl, Particle> SuperStupidStack;
-
-  } // namespace super_stupid
-
-} // namespace corsika::stack
-
-#endif
diff --git a/Stack/NuclearStackExtension/CMakeLists.txt b/Stack/NuclearStackExtension/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..f2f34cfcfa4cf6108651282b464c476de9c125a5
--- /dev/null
+++ b/Stack/NuclearStackExtension/CMakeLists.txt
@@ -0,0 +1,49 @@
+
+set (NuclearStackExtension_HEADERS NuclearStackExtension.h)
+set (NuclearStackExtension_NAMESPACE corsika/stack/nuclear_extension)
+
+add_library (NuclearStackExtension INTERFACE)
+
+CORSIKA_COPY_HEADERS_TO_NAMESPACE (NuclearStackExtension ${NuclearStackExtension_NAMESPACE} ${NuclearStackExtension_HEADERS})
+
+target_link_libraries (
+  NuclearStackExtension
+  INTERFACE
+  CORSIKAstackinterface
+  CORSIKAunits
+  CORSIKAparticles
+  CORSIKAgeometry
+  )
+
+target_include_directories (
+  NuclearStackExtension
+  INTERFACE
+  $<BUILD_INTERFACE:${PROJECT_BINARY_DIR}/include>
+  $<INSTALL_INTERFACE:include>
+  )
+
+install (
+  FILES
+  ${NuclearStackExtension_HEADERS}
+  DESTINATION
+  include/${NuclearStackExtension_NAMESPACE}
+  )
+
+# ----------------
+# code unit testing
+add_executable (
+  testNuclearStackExtension
+  testNuclearStackExtension.cc
+  )
+
+target_link_libraries (
+  testNuclearStackExtension
+  SuperStupidStack
+  NuclearStackExtension
+  #  CORSIKAutls
+  ProcessStackInspector
+  CORSIKAgeometry
+  CORSIKAunits
+  CORSIKAthirdparty # for catch2
+  )
+CORSIKA_ADD_TEST(testNuclearStackExtension)
diff --git a/Stack/NuclearStackExtension/NuclearStackExtension.h b/Stack/NuclearStackExtension/NuclearStackExtension.h
new file mode 100644
index 0000000000000000000000000000000000000000..4cdd876925d748c4ddbb6688f1dae5c56512b8d2
--- /dev/null
+++ b/Stack/NuclearStackExtension/NuclearStackExtension.h
@@ -0,0 +1,293 @@
+
+/*
+ * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+
+#ifndef _include_stack_nuclearstackextension_h_
+#define _include_stack_nuclearstackextension_h_
+
+#include <corsika/stack/super_stupid/SuperStupidStack.h>
+
+#include <corsika/particles/ParticleProperties.h>
+#include <corsika/stack/Stack.h>
+#include <corsika/units/PhysicalUnits.h>
+
+#include <corsika/geometry/Point.h>
+#include <corsika/geometry/Vector.h>
+
+#include <algorithm>
+#include <vector>
+
+namespace corsika::stack {
+
+  namespace nuclear_extension {
+
+    typedef corsika::geometry::Vector<corsika::units::si::hepmomentum_d> MomentumVector;
+
+    /**
+     * Define ParticleInterface for NuclearStackExtension Stack derived from
+     * ParticleInterface of Inner stack class
+     */
+    template <template <typename> typename InnerParticleInterface,
+              typename StackIteratorInterface>
+    class NuclearParticleInterface
+        : public InnerParticleInterface<StackIteratorInterface> {
+
+    public:
+      using ExtendedParticleInterface =
+          NuclearParticleInterface<InnerParticleInterface, StackIteratorInterface>;
+
+    protected:
+      using InnerParticleInterface<StackIteratorInterface>::GetStackData;
+      using InnerParticleInterface<StackIteratorInterface>::GetIndex;
+
+    public:
+      void SetParticleData(const corsika::particles::Code vDataPID,
+                           const corsika::units::si::HEPEnergyType vDataE,
+                           const MomentumVector& vMomentum,
+                           const corsika::geometry::Point& vPosition,
+                           const corsika::units::si::TimeType vTime, const int vA = 0,
+                           const int vZ = 0) {
+        if (vDataPID == corsika::particles::Code::Nucleus) {
+          if (vA == 0 || vZ == 0) {
+            std::ostringstream err;
+            err << "NuclearStackExtension: no A and Z specified for new Nucleus!";
+            throw std::runtime_error(err.str());
+          }
+          SetNuclearRef(
+              GetStackData().GetNucleusNextRef()); // store this nucleus data ref
+          SetNuclearA(vA);
+          SetNuclearZ(vZ);
+        } else {
+          SetNuclearRef(-1); // this is not a nucleus
+        }
+        InnerParticleInterface<StackIteratorInterface>::
+            SetParticleData(vDataPID, vDataE, vMomentum, vPosition, vTime);
+      }
+
+      void SetParticleData(InnerParticleInterface<StackIteratorInterface>&,
+                           const corsika::particles::Code vDataPID,
+                           const corsika::units::si::HEPEnergyType vDataE,
+                           const MomentumVector& vMomentum,
+                           const corsika::geometry::Point& vPosition,
+                           const corsika::units::si::TimeType vTime, const int vA = 0,
+                           const int vZ = 0) {
+        SetParticleData(vDataPID, vDataE, vMomentum, vPosition, vTime, vA, vZ);
+      }
+
+      /**
+       * @name individual setters
+       * @{
+       */
+      void SetNuclearA(const int vA) { GetStackData().SetNuclearA(GetIndex(), vA); }
+      void SetNuclearZ(const int vZ) { GetStackData().SetNuclearZ(GetIndex(), vZ); }
+      /// @}
+
+      /**
+       * @name individual getters
+       * @{
+       */
+      int GetNuclearA() const { return GetStackData().GetNuclearA(GetIndex()); }
+      int GetNuclearZ() const { return GetStackData().GetNuclearZ(GetIndex()); }
+      /// @}
+
+    protected:
+      void SetNuclearRef(const int vR) { GetStackData().SetNuclearRef(GetIndex(), vR); }
+    };
+
+    /**
+     * Memory implementation of the most simple (stupid) particle stack object.
+     */
+    template <typename InnerStackImpl>
+    class NuclearStackExtensionImpl : public InnerStackImpl {
+
+    public:
+      void Init() { InnerStackImpl::Init(); }
+
+      void Clear() {
+        InnerStackImpl::Clear();
+        fNuclearRef.clear();
+        fNuclearA.clear();
+        fNuclearZ.clear();
+      }
+
+      unsigned int GetSize() const { return fNuclearRef.size(); }
+      unsigned int GetCapacity() const { return fNuclearRef.size(); }
+
+      void SetNuclearA(const unsigned int i, const int vA) {
+        fNuclearA[GetNucleusRef(i)] = vA;
+      }
+      void SetNuclearZ(const unsigned int i, const int vZ) {
+        fNuclearZ[GetNucleusRef(i)] = vZ;
+      }
+      void SetNuclearRef(const unsigned int i, const int v) { fNuclearRef[i] = v; }
+
+      int GetNuclearA(const unsigned int i) const { return fNuclearA[GetNucleusRef(i)]; }
+      int GetNuclearZ(const unsigned int i) const { return fNuclearZ[GetNucleusRef(i)]; }
+      // this function will create new storage for Nuclear Properties, and return the
+      // reference to it
+      int GetNucleusNextRef() {
+        fNuclearA.push_back(0);
+        fNuclearZ.push_back(0);
+        return fNuclearA.size() - 1;
+      }
+
+      int GetNucleusRef(const unsigned int i) const {
+        if (fNuclearRef[i] >= 0) return fNuclearRef[i];
+        std::ostringstream err;
+        err << "NuclearStackExtension: no nucleus at ref=" << i;
+        throw std::runtime_error(err.str());
+      }
+
+      /**
+       *   Function to copy particle at location i1 in stack to i2
+       */
+      void Copy(const unsigned int i1, const unsigned int i2) {
+        if (i1 >= GetSize() || i2 >= GetSize()) {
+          std::ostringstream err;
+          err << "NuclearStackExtension: trying to access data beyond size of stack!";
+          throw std::runtime_error(err.str());
+        }
+        InnerStackImpl::Copy(i1, i2);
+        const int ref1 = fNuclearRef[i1];
+        const int ref2 = fNuclearRef[i2];
+        if (ref2 < 0) {
+          if (ref1 >= 0) {
+            // i1 is nucleus, i2 is not
+            fNuclearRef[i2] = GetNucleusNextRef();
+            fNuclearA[fNuclearRef[i2]] = fNuclearA[ref1];
+            fNuclearZ[fNuclearRef[i2]] = fNuclearZ[ref1];
+          } else {
+            // neither i1 nor i2 are nuclei
+          }
+        } else {
+          if (ref1 >= 0) {
+            // both are nuclei, i2 is overwritten with nucleus i1
+            fNuclearA[ref2] = fNuclearA[ref1];
+            fNuclearZ[ref2] = fNuclearZ[ref1];
+          } else {
+            // i2 is overwritten with non-nucleus i1
+            fNuclearA.erase(fNuclearA.cbegin() + ref2);
+            fNuclearZ.erase(fNuclearZ.cbegin() + ref2);
+            const int n = fNuclearRef.size();
+            for (int i = 0; i < n; ++i) {
+              if (fNuclearRef[i] >= ref2) { fNuclearRef[i] -= 1; }
+            }
+          }
+        }
+      }
+
+      /**
+       *   Function to copy particle at location i2 in stack to i1
+       */
+      void Swap(const unsigned int i1, const unsigned int i2) {
+        if (i1 >= GetSize() || i2 >= GetSize()) {
+          std::ostringstream err;
+          err << "NuclearStackExtension: trying to access data beyond size of stack!";
+          throw std::runtime_error(err.str());
+        }
+        InnerStackImpl::Swap(i1, i2);
+        const int ref1 = fNuclearRef[i1];
+        const int ref2 = fNuclearRef[i2];
+        if (ref2 < 0) {
+          if (ref1 >= 0) {
+            // i1 is nucleus, i2 is not
+            std::swap(fNuclearRef[i2], fNuclearRef[i1]);
+          } else {
+            // neither i1 nor i2 are nuclei
+          }
+        } else {
+          if (ref1 >= 0) {
+            // both are nuclei, i2 is overwritten with nucleus i1
+            std::swap(fNuclearRef[i2], fNuclearRef[i1]);
+          } else {
+            // i2 is overwritten with non-nucleus i1
+            std::swap(fNuclearRef[i2], fNuclearRef[i1]);
+          }
+        }
+      }
+
+    protected:
+      void IncrementSize() {
+        InnerStackImpl::IncrementSize();
+        fNuclearRef.push_back(-1);
+      }
+
+      void DecrementSize() {
+        InnerStackImpl::DecrementSize();
+        if (fNuclearRef.size() > 0) {
+          const int ref = fNuclearRef.back();
+          fNuclearRef.pop_back();
+          if (ref >= 0) {
+            fNuclearA.erase(fNuclearA.begin() + ref);
+            fNuclearZ.erase(fNuclearZ.begin() + ref);
+            const int n = fNuclearRef.size();
+            for (int i = 0; i < n; ++i) {
+              if (fNuclearRef[i] >= ref) { fNuclearRef[i] -= 1; }
+            }
+          }
+        }
+      }
+
+    private:
+      /// the actual memory to store particle data
+
+      std::vector<int> fNuclearRef;
+      std::vector<int> fNuclearA;
+      std::vector<int> fNuclearZ;
+
+    }; // end class NuclearStackExtensionImpl
+
+    //    template<typename StackIteratorInterface>
+    // using NuclearParticleInterfaceType<StackIteratorInterface> =
+    // NuclearParticleInterface< ,StackIteratorInterface>
+
+    // works, but requires stupd _PI class
+    // template<typename SS> using TEST =
+    // NuclearParticleInterface<corsika::stack::super_stupid::SuperStupidStack::PIType,
+    // SS>;
+    template <typename InnerStack, template <typename> typename _PI>
+    using NuclearStackExtension =
+        Stack<NuclearStackExtensionImpl<typename InnerStack::StackImpl>, _PI>;
+
+    // ----
+
+    // I'm dont't manage to do this properly.......
+    /*
+    template<typename TT, typename SS> using TESTi = typename
+    NuclearParticleInterface<TT::template PIType, SS>::ExtendedParticleInterface;
+    template<typename TT, typename SS> using TEST1 = TESTi<TT, SS>;
+    template<typename SS> using TEST2 = TEST1<typename
+    corsika::stack::super_stupid::SuperStupidStack, SS>;
+
+    using NuclearStackExtension = Stack<NuclearStackExtensionImpl<typename
+    InnerStack::StackImpl>, TEST2>;
+    */
+    /*
+      // .... this should be fixed ....
+
+    template <typename InnerStack, typename SS=StackIteratorInterface>
+      //using NuclearStackExtension = Stack<NuclearStackExtensionImpl<typename
+    InnerStack::StackImpl>, NuclearParticleInterface<typename InnerStack::template PIType,
+    StackIteratorInterface>::ExtendedParticleInterface>; using NuclearStackExtension =
+    Stack<NuclearStackExtensionImpl<typename InnerStack::StackImpl>, TEST1<typename
+    corsika::stack::super_stupid::SuperStupidStack, SS> >;
+
+    //template <typename InnerStack>
+      //  using NuclearStackExtension = Stack<NuclearStackExtensionImpl<typename
+    InnerStack::StackImpl>, TEST<typename
+    corsika::stack::super_stupid::SuperStupidStack::PIType>>;
+    //using NuclearStackExtension = Stack<NuclearStackExtensionImpl<typename
+    InnerStack::StackImpl>, TEST>;
+    */
+
+  } // namespace nuclear_extension
+} // namespace corsika::stack
+
+#endif
diff --git a/Stack/NuclearStackExtension/testNuclearStackExtension.cc b/Stack/NuclearStackExtension/testNuclearStackExtension.cc
new file mode 100644
index 0000000000000000000000000000000000000000..05983794970285750c462439bd3f1825515c56f7
--- /dev/null
+++ b/Stack/NuclearStackExtension/testNuclearStackExtension.cc
@@ -0,0 +1,221 @@
+
+/*
+ * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+
+#include <corsika/geometry/RootCoordinateSystem.h>
+#include <corsika/stack/nuclear_extension/NuclearStackExtension.h>
+#include <corsika/stack/super_stupid/SuperStupidStack.h>
+#include <corsika/units/PhysicalUnits.h>
+
+#include <boost/type_index.hpp>
+using boost::typeindex::type_id_with_cvr;
+
+using namespace corsika::geometry;
+using namespace corsika::units::si;
+
+#define CATCH_CONFIG_MAIN // This tells Catch to provide a main() - only do this in one
+                          // cpp file
+#include <catch2/catch.hpp>
+
+using namespace corsika;
+using namespace corsika::stack::nuclear_extension;
+
+// this is an auxiliary help typedef, which I don't know how to put
+// into NuclearStackExtension.h where it belongs...
+template <typename StackIter>
+using ExtendedParticleInterfaceType =
+    corsika::stack::nuclear_extension::NuclearParticleInterface<
+        corsika::stack::super_stupid::SuperStupidStack::PIType, StackIter>;
+
+using ExtStack = NuclearStackExtension<corsika::stack::super_stupid::SuperStupidStack,
+                                       ExtendedParticleInterfaceType>;
+
+#include <iostream>
+using namespace std;
+
+TEST_CASE("NuclearStackExtension", "[stack]") {
+
+  geometry::CoordinateSystem& dummyCS =
+      geometry::RootCoordinateSystem::GetInstance().GetRootCoordinateSystem();
+
+  // cout << "ParticleType=" << type_id_with_cvr<ParticleType>().pretty_name() << endl;
+
+  SECTION("write non nucleus") {
+    NuclearStackExtension<corsika::stack::super_stupid::SuperStupidStack,
+                          ExtendedParticleInterfaceType>
+        s;
+    s.AddParticle(particles::Code::Electron, 1.5_GeV,
+                  MomentumVector(dummyCS, {1_GeV, 1_GeV, 1_GeV}),
+                  Point(dummyCS, {1 * meter, 1 * meter, 1 * meter}), 100_s);
+    REQUIRE(s.GetSize() == 1);
+  }
+
+  SECTION("write nucleus") {
+    NuclearStackExtension<corsika::stack::super_stupid::SuperStupidStack,
+                          ExtendedParticleInterfaceType>
+        s;
+    s.AddParticle(particles::Code::Nucleus, 1.5_GeV,
+                  MomentumVector(dummyCS, {1_GeV, 1_GeV, 1_GeV}),
+                  Point(dummyCS, {1 * meter, 1 * meter, 1 * meter}), 100_s, 10, 10);
+    REQUIRE(s.GetSize() == 1);
+  }
+
+  SECTION("write invalid nucleus") {
+    ExtStack s;
+    REQUIRE_THROWS(s.AddParticle(
+        particles::Code::Nucleus, 1.5_GeV, MomentumVector(dummyCS, {1_GeV, 1_GeV, 1_GeV}),
+        Point(dummyCS, {1 * meter, 1 * meter, 1 * meter}), 100_s, 0, 0));
+  }
+
+  SECTION("read non nucleus") {
+    ExtStack s;
+    s.AddParticle(particles::Code::Electron, 1.5_GeV,
+                  MomentumVector(dummyCS, {1_GeV, 1_GeV, 1_GeV}),
+                  Point(dummyCS, {1 * meter, 1 * meter, 1 * meter}), 100_s);
+    const auto pout = s.GetNextParticle();
+    REQUIRE(pout.GetPID() == particles::Code::Electron);
+    REQUIRE(pout.GetEnergy() == 1.5_GeV);
+    REQUIRE(pout.GetTime() == 100_s);
+  }
+
+  SECTION("read nucleus") {
+    ExtStack s;
+    s.AddParticle(particles::Code::Nucleus, 1.5_GeV,
+                  MomentumVector(dummyCS, {1_GeV, 1_GeV, 1_GeV}),
+                  Point(dummyCS, {1 * meter, 1 * meter, 1 * meter}), 100_s, 10, 9);
+    const auto pout = s.GetNextParticle();
+    REQUIRE(pout.GetPID() == particles::Code::Nucleus);
+    REQUIRE(pout.GetEnergy() == 1.5_GeV);
+    REQUIRE(pout.GetTime() == 100_s);
+    REQUIRE(pout.GetNuclearA() == 10);
+    REQUIRE(pout.GetNuclearZ() == 9);
+  }
+
+  SECTION("read invalid nucleus") {
+    ExtStack s;
+    s.AddParticle(particles::Code::Electron, 1.5_GeV,
+                  MomentumVector(dummyCS, {1_GeV, 1_GeV, 1_GeV}),
+                  Point(dummyCS, {1 * meter, 1 * meter, 1 * meter}), 100_s);
+    const auto pout = s.GetNextParticle();
+    REQUIRE_THROWS(pout.GetNuclearA());
+    REQUIRE_THROWS(pout.GetNuclearZ());
+  }
+
+  SECTION("stack fill and cleanup") {
+
+    ExtStack s;
+    // add 99 particles, each 10th particle is a nucleus with A=i and Z=A/2!
+    for (int i = 0; i < 99; ++i) {
+      if ((i + 1) % 10 == 0) {
+        s.AddParticle(particles::Code::Nucleus, 1.5_GeV,
+                      MomentumVector(dummyCS, {1_GeV, 1_GeV, 1_GeV}),
+                      Point(dummyCS, {1 * meter, 1 * meter, 1 * meter}), 100_s, i, i / 2);
+      } else {
+        s.AddParticle(particles::Code::Electron, 1.5_GeV,
+                      MomentumVector(dummyCS, {1_GeV, 1_GeV, 1_GeV}),
+                      Point(dummyCS, {1 * meter, 1 * meter, 1 * meter}), 100_s);
+      }
+    }
+
+    REQUIRE(s.GetSize() == 99);
+    for (int i = 0; i < 99; ++i) s.GetNextParticle().Delete();
+    REQUIRE(s.GetSize() == 0);
+  }
+
+  SECTION("stack operations") {
+
+    ExtStack s;
+    // add 99 particles, each 10th particle is a nucleus with A=i and Z=A/2!
+    for (int i = 0; i < 99; ++i) {
+      if ((i + 1) % 10 == 0) {
+        s.AddParticle(particles::Code::Nucleus, i * 15_GeV,
+                      MomentumVector(dummyCS, {1_GeV, 1_GeV, 1_GeV}),
+                      Point(dummyCS, {1 * meter, 1 * meter, 1 * meter}), 100_s, i, i / 2);
+      } else {
+        s.AddParticle(particles::Code::Electron, i * 1.5_GeV,
+                      MomentumVector(dummyCS, {1_GeV, 1_GeV, 1_GeV}),
+                      Point(dummyCS, {1 * meter, 1 * meter, 1 * meter}), 100_s);
+      }
+    }
+
+    // copy
+    {
+      s.Copy(s.begin() + 9, s.begin() + 10); // nuclei to non-nuclei
+      const auto& p9 = s.cbegin() + 9;
+      const auto& p10 = s.cbegin() + 10;
+
+      REQUIRE(p9.GetPID() == particles::Code::Nucleus);
+      REQUIRE(p9.GetEnergy() == 9 * 15_GeV);
+      REQUIRE(p9.GetTime() == 100_s);
+      REQUIRE(p9.GetNuclearA() == 9);
+      REQUIRE(p9.GetNuclearZ() == 9 / 2);
+
+      REQUIRE(p10.GetPID() == particles::Code::Nucleus);
+      REQUIRE(p10.GetEnergy() == 9 * 15_GeV);
+      REQUIRE(p10.GetTime() == 100_s);
+      REQUIRE(p10.GetNuclearA() == 9);
+      REQUIRE(p10.GetNuclearZ() == 9 / 2);
+    }
+
+    // copy
+    {
+      s.Copy(s.begin() + 93, s.begin() + 9); // non-nuclei to nuclei
+      const auto& p93 = s.cbegin() + 93;
+      const auto& p9 = s.cbegin() + 9;
+
+      REQUIRE(p9.GetPID() == particles::Code::Electron);
+      REQUIRE(p9.GetEnergy() == 93 * 1.5_GeV);
+      REQUIRE(p9.GetTime() == 100_s);
+
+      REQUIRE(p93.GetPID() == particles::Code::Electron);
+      REQUIRE(p93.GetEnergy() == 93 * 1.5_GeV);
+      REQUIRE(p93.GetTime() == 100_s);
+    }
+
+    // swap
+    {
+      s.Swap(s.begin() + 11, s.begin() + 10);
+      const auto& p11 = s.cbegin() + 11; // now: nucleus
+      const auto& p10 = s.cbegin() + 10; // now: electron
+
+      REQUIRE(p11.GetPID() == particles::Code::Nucleus);
+      REQUIRE(p11.GetEnergy() == 9 * 15_GeV);
+      REQUIRE(p11.GetTime() == 100_s);
+      REQUIRE(p11.GetNuclearA() == 9);
+      REQUIRE(p11.GetNuclearZ() == 9 / 2);
+
+      REQUIRE(p10.GetPID() == particles::Code::Electron);
+      REQUIRE(p10.GetEnergy() == 11 * 1.5_GeV);
+      REQUIRE(p10.GetTime() == 100_s);
+    }
+
+    // swap two nuclei
+    {
+      s.Swap(s.begin() + 29, s.begin() + 59);
+      const auto& p29 = s.cbegin() + 29;
+      const auto& p59 = s.cbegin() + 59;
+
+      REQUIRE(p29.GetPID() == particles::Code::Nucleus);
+      REQUIRE(p29.GetEnergy() == 59 * 15_GeV);
+      REQUIRE(p29.GetTime() == 100_s);
+      REQUIRE(p29.GetNuclearA() == 59);
+      REQUIRE(p29.GetNuclearZ() == 59 / 2);
+
+      REQUIRE(p59.GetPID() == particles::Code::Nucleus);
+      REQUIRE(p59.GetEnergy() == 29 * 15_GeV);
+      REQUIRE(p59.GetTime() == 100_s);
+      REQUIRE(p59.GetNuclearA() == 29);
+      REQUIRE(p59.GetNuclearZ() == 29 / 2);
+    }
+
+    for (int i = 0; i < 99; ++i) s.DeleteLast();
+    REQUIRE(s.GetSize() == 0);
+  }
+}
diff --git a/Stack/SuperStupidStack/SuperStupidStack.h b/Stack/SuperStupidStack/SuperStupidStack.h
index c3d4e633685195682038dcca39af5d2475bea32c..3daa825abcabef8d81bbb058cf62d2338933d1df 100644
--- a/Stack/SuperStupidStack/SuperStupidStack.h
+++ b/Stack/SuperStupidStack/SuperStupidStack.h
@@ -23,8 +23,6 @@
 #include <algorithm>
 #include <vector>
 
-using namespace corsika;
-
 namespace corsika::stack {
 
   namespace super_stupid {
@@ -38,20 +36,12 @@ namespace corsika::stack {
     template <typename StackIteratorInterface>
     class ParticleInterface : public ParticleBase<StackIteratorInterface> {
 
+    protected:
       using corsika::stack::ParticleBase<StackIteratorInterface>::GetStack;
       using corsika::stack::ParticleBase<StackIteratorInterface>::GetStackData;
       using corsika::stack::ParticleBase<StackIteratorInterface>::GetIndex;
 
     public:
-      /// the factory function, this is how to create a new particle:
-      /*void AddSecondary(const corsika::particles::Code vDataPID,
-                        const corsika::units::si::HEPEnergyType vDataE,
-                        const MomentumVector& vMomentum,
-                        const corsika::geometry::Point& vPosition,
-                        const corsika::units::si::TimeType vTime) {
-        GetStack().AddParticle(vDataPID, vDataE, vMomentum, vPosition, vTime);
-        }*/
-      //
       void SetParticleData(const corsika::particles::Code vDataPID,
                            const corsika::units::si::HEPEnergyType vDataE,
                            const MomentumVector& vMomentum,
@@ -117,7 +107,6 @@ namespace corsika::stack {
     };
 
     /**
-     *
      * Memory implementation of the most simple (stupid) particle stack object.
      */
 
@@ -134,29 +123,41 @@ namespace corsika::stack {
         fTime.clear();
       }
 
-      int GetSize() const { return fDataPID.size(); }
-      int GetCapacity() const { return fDataPID.size(); }
+      unsigned int GetSize() const { return fDataPID.size(); }
+      unsigned int GetCapacity() const { return fDataPID.size(); }
 
-      void SetPID(const int i, const corsika::particles::Code id) { fDataPID[i] = id; }
-      void SetEnergy(const int i, const corsika::units::si::HEPEnergyType e) {
+      void SetPID(const unsigned int i, const corsika::particles::Code id) {
+        fDataPID[i] = id;
+      }
+      void SetEnergy(const unsigned int i, const corsika::units::si::HEPEnergyType e) {
         fDataE[i] = e;
       }
-      void SetMomentum(const int i, const MomentumVector& v) { fMomentum[i] = v; }
-      void SetPosition(const int i, const corsika::geometry::Point& v) {
+      void SetMomentum(const unsigned int i, const MomentumVector& v) {
+        fMomentum[i] = v;
+      }
+      void SetPosition(const unsigned int i, const corsika::geometry::Point& v) {
         fPosition[i] = v;
       }
-      void SetTime(const int i, const corsika::units::si::TimeType& v) { fTime[i] = v; }
+      void SetTime(const unsigned int i, const corsika::units::si::TimeType& v) {
+        fTime[i] = v;
+      }
 
-      corsika::particles::Code GetPID(const int i) const { return fDataPID[i]; }
-      corsika::units::si::HEPEnergyType GetEnergy(const int i) const { return fDataE[i]; }
-      MomentumVector GetMomentum(const int i) const { return fMomentum[i]; }
-      corsika::geometry::Point GetPosition(const int i) const { return fPosition[i]; }
-      corsika::units::si::TimeType GetTime(const int i) const { return fTime[i]; }
+      corsika::particles::Code GetPID(const unsigned int i) const { return fDataPID[i]; }
+      corsika::units::si::HEPEnergyType GetEnergy(const unsigned int i) const {
+        return fDataE[i];
+      }
+      MomentumVector GetMomentum(const unsigned int i) const { return fMomentum[i]; }
+      corsika::geometry::Point GetPosition(const unsigned int i) const {
+        return fPosition[i];
+      }
+      corsika::units::si::TimeType GetTime(const unsigned int i) const {
+        return fTime[i];
+      }
 
       /**
        *   Function to copy particle at location i2 in stack to i1
        */
-      void Copy(const int i1, const int i2) {
+      void Copy(const unsigned int i1, const unsigned int i2) {
         fDataPID[i2] = fDataPID[i1];
         fDataE[i2] = fDataE[i1];
         fMomentum[i2] = fMomentum[i1];
@@ -167,7 +168,7 @@ namespace corsika::stack {
       /**
        *   Function to copy particle at location i2 in stack to i1
        */
-      void Swap(const int i1, const int i2) {
+      void Swap(const unsigned int i1, const unsigned int i2) {
         std::swap(fDataPID[i2], fDataPID[i1]);
         std::swap(fDataE[i2], fDataE[i1]);
         std::swap(fMomentum[i2], fMomentum[i1]);
@@ -181,7 +182,6 @@ namespace corsika::stack {
         using corsika::particles::Code;
         fDataPID.push_back(Code::Unknown);
         fDataE.push_back(0 * corsika::units::si::electronvolt);
-        //#TODO this here makes no sense: see issue #48
         geometry::CoordinateSystem& dummyCS =
             geometry::RootCoordinateSystem::GetInstance().GetRootCoordinateSystem();
         fMomentum.push_back(MomentumVector(
diff --git a/do-copyright.py b/do-copyright.py
index 10f60e389b1cfe27d476bb3f1936fc1c11cd57e1..da81b79521919a69805f63df66bffc1766c0af4f 100755
--- a/do-copyright.py
+++ b/do-copyright.py
@@ -21,43 +21,48 @@ extensions = [".cc", ".h", ".test"]
 
 def checkNote(filename):
 
-    startNote = -1
-    endNote = -1
-    isCopyright = False
-    lines = []
+    startNote = []
+    endNote = []
 
+    # read input file into lines
+    lines = []
     with open(filename, "r") as file:
         for line in file.readlines():
             lines.append(line)            
         file.close()
 
-            
+    searchStatus = 0 # 0:before comment block, #1 in comment block, #2 found copyright
+    blockStart = 0
     for iLine in range(len(lines)):
         line = lines[iLine]
-        if "/**" in line and startNote == -1:
-            startNote = iLine
-        if "copyright" in line.lower() and startNote>=0 and endNote==-1:
-            isCopyright = True
-        if "*/" in line and startNote>=0 and endNote==-1:
-            endNote = iLine
+        if "/*" in line and searchStatus==0:
+            searchStatus = 1
+            blockStart = iLine
+        if "copyright" in line.lower() and searchStatus>0:
+            searchStatus = 2
+        if "*/" in line:
+            if searchStatus>=2:
+                startNote.append(blockStart)
+                endNote.append(iLine)
+            searchStatus = 0
         iLine += 1
-
-    # now check if copyright notice is already there and identical...
+        
+    # now check if first copyright notices is already identical...
     isSame = False
-    if startNote>=0 and endNote>=0 and isCopyright:
+    if len(startNote)>0: 
         isSame = True
         noteLines = text.split('\n')        
         for iLine in range(len(noteLines)-2):
-            if startNote+iLine >= len(lines):
+            if startNote[0]+iLine >= len(lines):
                 isSame = False
-                break            
-            if noteLines[iLine+1].strip(" \n") != lines[startNote+iLine].strip(" \n"):
+                break
+            if noteLines[iLine+1].strip(" \n") != lines[startNote[0]+iLine].strip(" \n"):
                 isSame = False
-                print "not same: " + filename + " new=\'" + noteLines[iLine+1] + "\' vs old=\'" + lines[startNote+iLine].rstrip('\n') + "\'"
+                print "need update: " + filename + " new=\'" + noteLines[iLine+1] + "\' vs old=\'" + lines[startNote+iLine].rstrip('\n') + "\'"
                 break
-
-    # check if notice is the same 
-    if isSame:
+    
+    # check if notice is the same, or we need to remove multiple notices...
+    if isSame and len(startNote)<=1:
         return                
     
     # add (new) copyright notice here:
@@ -68,12 +73,24 @@ def checkNote(filename):
 
         file.write(text)
 
-        firstLine = 0
-        if startNote>=0 and endNote>=0 and isCopyright:
-            firstLine = endNote + 2
+        skip = False
+        for iLine in range(len(lines)):
+
+            inBlock = False
+            for iBlock in range(len(startNote)):
+                if iLine>=startNote[iBlock] and iLine<=endNote[iBlock]:
+                    print "   [remove " + str(iBlock) + "] " + (lines[iLine]).strip()
+                    inBlock = True
+                    skip = True
 
-        for iLine in range(firstLine, len(lines)):
-            file.write(lines[iLine])
+            if inBlock:
+                continue
+            
+            if lines[iLine].strip() != "": # if line after comment is empty, also remove it
+                skip = False
+                    
+            if not skip:
+                file.write(lines[iLine])
         
         file.close()
 
@@ -84,6 +101,8 @@ def next_file(x, dir_name, files):
             return
     for check in files :
         filename, file_extension = os.path.splitext(check)
+        if '#' in check or '~' in check:
+            return
         for check2 in excludeFiles :
             if check2 in check:
                 return
diff --git a/test.cc b/test.cc
new file mode 100644
index 0000000000000000000000000000000000000000..6839def48c2f6fddfafd0beb5fd4ca46058ebd09
--- /dev/null
+++ b/test.cc
@@ -0,0 +1,10 @@
+
+/*
+ * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */