diff --git a/CMakeLists.txt b/CMakeLists.txt
index 919c73ff..8af78472 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -6,7 +6,7 @@ cmake_minimum_required( VERSION 3.1 )
 
 # Define the project
 cmake_policy( SET CMP0048 NEW ) # version in project()
-project( locust_mc VERSION 1.16.0)
+project( locust_mc VERSION 1.16.1)
 
 list( APPEND CMAKE_MODULE_PATH ${PROJECT_SOURCE_DIR}/Scarab/cmake )
 include( PackageBuilder )
@@ -41,10 +41,19 @@ set (LOCUST_MC_INCLUDE_DIRECTORIES
     ${CMAKE_CURRENT_SOURCE_DIR}/Source/RxComponents/
     ${CMAKE_CURRENT_SOURCE_DIR}/Source/Transforms/
     ${CMAKE_CURRENT_SOURCE_DIR}/Source/Core/
+    ${CMAKE_CURRENT_SOURCE_DIR}/Source/IO/
     ${CMAKE_CURRENT_SOURCE_DIR}/Source/Distributions/
     ${CMAKE_CURRENT_SOURCE_DIR}/Source/Kassiopeia/
 )
 
+#if (locust_mc_BUILD_WITH_ROOT)
+#    set( LOCUST_MC_INCLUDE_DIRECTORIES ${LOCUST_MC_INCUDE_DIRECTORIES}
+#        ${CMAKE_CURRENT_SOURCE_DIR}/Source/IO/
+#        )
+#endif (locust_mc_BUILD_WITH_ROOT)
+
+        
+
 include_directories( BEFORE ${LOCUST_MC_INCLUDE_DIRECTORIES} )
 
 
diff --git a/Dockerfile b/Dockerfile
index 72007014..4124f727 100644
--- a/Dockerfile
+++ b/Dockerfile
@@ -3,7 +3,7 @@ FROM project8/p8compute_dependencies:v0.9.0 as locust_common
 ARG build_type=Release
 ENV LOCUST_BUILD_TYPE=$build_type
 
-ENV LOCUST_TAG=v1.16.0
+ENV LOCUST_TAG=v1.16.1
 ENV LOCUST_BUILD_PREFIX=/usr/local/p8/locust/$LOCUST_TAG
 
 RUN mkdir -p $LOCUST_BUILD_PREFIX &&\
diff --git a/Source/CMakeLists.txt b/Source/CMakeLists.txt
index 5f183122..8aae9d2c 100644
--- a/Source/CMakeLists.txt
+++ b/Source/CMakeLists.txt
@@ -79,7 +79,7 @@ set( LOCUST_MC_HEADER_FILES
     RxComponents/LMCPatchAntenna.hh
     RxComponents/LMCSlotAntenna.hh
     RxComponents/LMCChannel.hh
-    
+
 )
 
 
@@ -139,7 +139,7 @@ set( LOCUST_MC_SOURCE_FILES
     RxComponents/LMCReceiver.cc
     RxComponents/LMCPatchAntenna.cc
     RxComponents/LMCSlotAntenna.cc
-    RxComponents/LMCChannel.cc
+    RxComponents/LMCChannel.cc    
 
 )
 
@@ -197,14 +197,24 @@ if (locust_mc_BUILD_WITH_ROOT)
     set( LOCUST_MC_HEADER_FILES ${LOCUST_MC_HEADER_FILES}
         Core/LMCEvent.hh
         Core/LMCTrack.hh
-        Generators/LMCFakeTrackSignalGenerator.hh
+        Core/LMCRunParameters.hh
+        Generators/LMCFakeTrackSignalGenerator.hh    
+	    IO/LMCFileWriter.hh
+    	IO/LMCRootTreeWriter.hh
+    	IO/LMCRootHistoWriter.hh
+    	IO/LMCRootGraphWriter.hh
        )
 
     set( LOCUST_MC_SOURCE_FILES ${LOCUST_MC_SOURCE_FILES}
         Core/LMCEvent.cc
         Core/LMCTrack.cc
-        Generators/LMCFakeTrackSignalGenerator.cc
-       )
+        Core/LMCRunParameters.cc
+        Generators/LMCFakeTrackSignalGenerator.cc   
+	    IO/LMCFileWriter.cc
+    	IO/LMCRootTreeWriter.cc
+    	IO/LMCRootHistoWriter.cc
+    	IO/LMCRootGraphWriter.cc
+        )
 
 
     if (CMAKE_COMPILER_IS_GNUCXX)
@@ -235,6 +245,7 @@ if (locust_mc_BUILD_WITH_ROOT)
     endforeach()
     add_custom_command(OUTPUT ${CMAKE_CURRENT_BINARY_DIR}/LMCDict.cxx 
         COMMAND rootcling -f LMCDict.cxx -c -inlineInputHeader ${includedirs}
+            ${PROJECT_SOURCE_DIR}/Source/Core/LMCRunParameters.hh 
             ${PROJECT_SOURCE_DIR}/Source/Core/LMCEvent.hh 
             ${PROJECT_SOURCE_DIR}/Source/Core/LMCTrack.hh 
             ${PROJECT_SOURCE_DIR}/Source/Core/LinkDef/LMCLinkDef.hh 
diff --git a/Source/Core/LMCRunParameters.cc b/Source/Core/LMCRunParameters.cc
new file mode 100644
index 00000000..6a389aa4
--- /dev/null
+++ b/Source/Core/LMCRunParameters.cc
@@ -0,0 +1,23 @@
+/*
+ * LMCRunParameters.cc
+ *
+ * This class has access to both locust and ROOT libraries.  The syntax is
+ * consistent with KTROOTData.cc and the instructions in
+ * https://root.cern.ch/root/Using.html .   It is also mentioned in LMCEventLinkDef.hh .
+ *
+ *  Created on: Jul 7, 2020
+ *      Author: pslocum
+ */
+
+
+#include "LMCRunParameters.hh"
+#include <iostream>
+
+ClassImp(locust::RunParameters);
+
+namespace locust
+{
+    RunParameters::RunParameters() {}
+    RunParameters::~RunParameters() {}
+
+}
diff --git a/Source/Core/LMCRunParameters.hh b/Source/Core/LMCRunParameters.hh
new file mode 100644
index 00000000..5bde8157
--- /dev/null
+++ b/Source/Core/LMCRunParameters.hh
@@ -0,0 +1,37 @@
+/*
+ * LMCRunParameters.hh
+ *
+ * This class has access to both locust and ROOT libraries.  The syntax is
+ * consistent with KTROOTData.hh and the instructions in
+ * https://root.cern.ch/root/Using.html .  It is also mentioned in LMCEventLinkDef.hh .
+ *  Created on: Jul 7, 2020
+ *      Author: pslocum
+ */
+
+
+
+#ifndef LMCRUNPARAMETERS_HH_
+#define LMCRUNPARAMETERS_HH_
+
+#include "TObject.h"
+#include "LMCRunParameters.hh"
+
+namespace locust
+{
+
+    class RunParameters : public TObject
+    {
+
+        public:
+            RunParameters();
+            virtual ~RunParameters();
+
+            double fNoise;
+            double fLOfrequency;
+
+            ClassDef(RunParameters,1)  // Root syntax.
+
+    };
+
+}
+#endif /* LMCRUNPARAMETERS_HH_ */
diff --git a/Source/Core/LMCSingleton.hh b/Source/Core/LMCSingleton.hh
new file mode 100644
index 00000000..6da50fbe
--- /dev/null
+++ b/Source/Core/LMCSingleton.hh
@@ -0,0 +1,147 @@
+/*
+ * singleton.hh
+ *
+ *  Created on: Nov 7, 2011
+ *      Author: nsoblath
+ */
+
+#ifndef SCARAB_SINGLETON_HH_
+#define SCARAB_SINGLETON_HH_
+
+#include "destroyer.hh"
+#include "error.hh"
+
+#include <cstddef>
+#include <mutex>
+
+namespace scarab
+{
+    /*!
+     \def allow_singleton_access( class_name )
+     Gives friend access to your class from scarab::singleton< class_name > and scarab::destroyer< class_name > 
+     so that those classes can control the life of your singleton class.
+    */
+#define allow_singleton_access( class_name ) \
+    friend class scarab::singleton< class_name >; \
+    friend class scarab::destroyer< class_name >;
+
+    /*!
+      @class singleton
+      @author N.S. Oblath
+      @brief Base class that turns a class into a singleton
+
+      @details
+      To use:
+        1. Inherit your class from singleton< your_class >
+        2. Make your constructor and destructor protected (or private)
+        3. Add the `allow_singleton_access( your_class )` macro to your class definition to allow the base classes to access your class
+
+      The mutex f_mutex provides thread safety for creation and destruction of an instance of the singleton.  It's available to the 
+      derived class for thread-safe access to the object in general.
+    */
+    template< class x_type >
+    class singleton
+    {
+        public:
+            static x_type* get_instance();
+            static void kill_instance();
+
+            template< class... x_args >
+            static x_type* create_instance( x_args... args );
+
+        private:
+            static void construct_instance();
+            static void delete_instance();
+
+        private:
+            static x_type* f_instance;
+            static destroyer< x_type > f_destroyer;
+
+        protected:
+            static std::mutex f_mutex;
+
+        protected:
+            singleton();
+
+            friend class destroyer< x_type >;
+            ~singleton();
+    };
+
+    template< class x_type >
+    x_type* singleton< x_type >::f_instance = nullptr;
+
+    template< class x_type >
+    destroyer< x_type > singleton< x_type >::f_destroyer;
+
+    template< class x_type >
+    std::mutex singleton< x_type >::f_mutex;
+
+    template< class x_type >
+    x_type* singleton< x_type >::get_instance()
+    {
+        if( f_instance == nullptr )
+        {
+            std::unique_lock< std::mutex > t_lock( f_mutex );
+            construct_instance();
+        }
+        return f_instance;
+    }
+
+    template< class x_type >
+    void singleton< x_type >::kill_instance()
+    {
+        if( f_instance != nullptr )
+        {
+            std::unique_lock< std::mutex > t_lock( f_mutex );
+            delete_instance();
+        }
+        return;
+    }
+
+    template< class x_type >
+    template< class... x_args >
+    x_type* singleton< x_type >::create_instance( x_args... args )
+    {
+        if( f_instance != nullptr )
+        {
+            throw error() << "Instance already exists; create_instance can only be called before the instance exists";
+        }
+        std::unique_lock< std::mutex > t_lock( f_mutex );
+        f_instance = new x_type( args... );
+        f_destroyer.set_doomed( f_instance );
+        return f_instance;
+    }
+
+    template< class x_type >
+    void singleton< x_type >::construct_instance()
+    {
+        if( f_instance == nullptr )
+        {
+            f_instance = new x_type();
+            f_destroyer.set_doomed( f_instance );
+        }
+    }
+
+    template< class x_type >
+    void singleton< x_type >::delete_instance()
+    {
+        if( f_instance != nullptr )
+        {
+            delete f_instance;
+            f_instance = nullptr;
+            f_destroyer.set_doomed( nullptr );
+        }
+    }
+
+    template< class x_type >
+    singleton< x_type >::singleton()
+    {
+    }
+    template< class x_type >
+    singleton< x_type >::~singleton()
+    {
+    }
+
+} /* namespace scarab */
+
+#endif /* SCARAB_SINGLETON_HH_ */
diff --git a/Source/Core/LinkDef/LMCLinkDef.hh b/Source/Core/LinkDef/LMCLinkDef.hh
index 714332e6..34641280 100644
--- a/Source/Core/LinkDef/LMCLinkDef.hh
+++ b/Source/Core/LinkDef/LMCLinkDef.hh
@@ -9,6 +9,7 @@
 
 #pragma link C++ namespace locust;
 
+#pragma link C++ class locust::RunParameters+;
 #pragma link C++ class locust::Event+;
 #pragma link C++ class locust::Track+;
 
diff --git a/Source/Distributions/LMCDistributionInterface.cc b/Source/Distributions/LMCDistributionInterface.cc
index 1a0b1788..ed89f8dc 100644
--- a/Source/Distributions/LMCDistributionInterface.cc
+++ b/Source/Distributions/LMCDistributionInterface.cc
@@ -43,7 +43,7 @@ namespace locust
         if(dist_name == "dirac" || dist_name == "fixed" )
             fDistributionList.push_back( std::make_shared< DiracDistribution >(aParam) );
 
-        if(dist_name == "exponential")
+        else if(dist_name == "exponential")
             fDistributionList.push_back( std::make_shared< ExponentialDistribution >(aParam) );
 
         else if(dist_name == "gaussian")
diff --git a/Source/Distributions/LMCExponentialDistribution.cc b/Source/Distributions/LMCExponentialDistribution.cc
index 0259e31c..d7dc1e16 100644
--- a/Source/Distributions/LMCExponentialDistribution.cc
+++ b/Source/Distributions/LMCExponentialDistribution.cc
@@ -12,24 +12,28 @@ namespace locust
     LOGGER( lmclog, "LMCExponentialDistribution" );
 
     ExponentialDistribution::ExponentialDistribution(const scarab::param_node &aParam) :
-        fLambda( 1. )
+        fTau( 1. ),
+        fShift( 0. )
     {
-        if(aParam.has("lambda"))
-            fLambda = aParam.get_value< double >( "lambda", fLambda );
+        if(aParam.has("tau"))
+            fTau = aParam.get_value< double >( "tau", fTau );
+        if(aParam.has("shift"))
+            fShift = aParam.get_value< double >( "shift", fShift );
 
-        fDistribution = std::exponential_distribution<double>(fLambda);
-        LDEBUG( lmclog, "Created exponential distribution. lambda: " <<fLambda);
+        fDistribution = std::exponential_distribution<double>(1. / fTau);
+        LDEBUG( lmclog, "Created exponential distribution. tau: " <<fTau<<" shift: "<<fShift);
     }
 
-    ExponentialDistribution::ExponentialDistribution(const double &aLambda) :
-        fLambda( aLambda ),
-        fDistribution( aLambda )
+    ExponentialDistribution::ExponentialDistribution(const double &aTau, const double &aShift) :
+        fTau( aTau ),
+        fShift( aShift ),
+        fDistribution( aTau )
     {
     }
 
     double ExponentialDistribution::Generate()
     {
-        return fDistribution(*fRNEngine);
+        return fDistribution(*fRNEngine) + fShift;
     }
 
 } /* namespace locust */
diff --git a/Source/Distributions/LMCExponentialDistribution.hh b/Source/Distributions/LMCExponentialDistribution.hh
index 869f9c46..fdf74b6f 100644
--- a/Source/Distributions/LMCExponentialDistribution.hh
+++ b/Source/Distributions/LMCExponentialDistribution.hh
@@ -26,13 +26,14 @@ namespace locust
 
         public:
             ExponentialDistribution(const scarab::param_node &aParam);
-            ExponentialDistribution(const double &aLambda);
+            ExponentialDistribution(const double &aTau, const double &aShift);
 
             double Generate();
 
         private:
             std::exponential_distribution<double> fDistribution;
-            double fLambda;
+            double fTau;
+            double fShift; //shift of exponential distribution from 0
 };
 
 
diff --git a/Source/Distributions/LMCKrComplexLineDistribution.cc b/Source/Distributions/LMCKrComplexLineDistribution.cc
index a12ad1d7..734a2ec5 100644
--- a/Source/Distributions/LMCKrComplexLineDistribution.cc
+++ b/Source/Distributions/LMCKrComplexLineDistribution.cc
@@ -22,7 +22,7 @@ namespace locust
         fFWHM( 5. ),
         fLinePosition( 17826. ),
         fAmplitude{1,0},
-        fScatterProbability{0.,0.},
+        fScatterProbability(0.),
         fNPointsSELA(10000),
         fGases{"H2","Kr"},
         fEmittedPeak("shake")
@@ -33,6 +33,8 @@ namespace locust
             fLinePosition = aParam.get_value< double >( "line-position", fLinePosition );
         if(aParam.has("emitted-peak"))
             fEmittedPeak = aParam.get_value< std::string >( "emitted-peak", fEmittedPeak );
+        if(aParam.has("scatter-probability"))
+            fScatterProbability = aParam.get_value< double >( "scatter-probability", fScatterProbability );
 
         for(unsigned i=0; i < fGases.size(); ++i)
             fGasIndex.insert( std::pair<std::string, unsigned>(fGases[i], i) );
@@ -45,14 +47,6 @@ namespace locust
                     fAmplitude[fGasIndex[fGases[i]]] = aAmplitudeNode.get_value<double >( fGases[i], fAmplitude[i]) ;
         }
 
-        if(aParam.has("scatter-probability"))
-        {
-            scarab::param_node aScatterNode = aParam["scatter-probability"].as_node();
-            for(unsigned i=0;i<fGases.size();++i)
-                if(aScatterNode.has(fGases[i]))
-                    fScatterProbability[fGasIndex[fGases[i]]] = aScatterNode.get_value<double >( fGases[i], fScatterProbability[i]) ;
-        }
-
         //create random number generator instances
         const double kr_line_width = 2.83; // eV
         fXArray = linspace(-1000,1000,fNPointsSELA);
@@ -61,10 +55,7 @@ namespace locust
         fLorentzian = std::cauchy_distribution<double>(0, kr_line_width / 2. ); //pos, hwhm
         fUniform = std::uniform_real_distribution<double>(0,1);
 
-
-        for(unsigned i=0; i < fGases.size(); ++i)
-            fGeometric.push_back(std::geometric_distribution<int>(1. - fScatterProbability[i]));
-
+        fGeometric = std::geometric_distribution<int>(1. - fScatterProbability);
 
         fShakeAccelerator = gsl_interp_accel_alloc();
 
@@ -89,9 +80,9 @@ namespace locust
         for(unsigned i=0;i<fGases.size();++i)
         {
             amp_log_string += std::to_string(fAmplitude[i]) + "  ";
-            scatter_log_string += std::to_string(fScatterProbability[i]) + "  ";
         }
-        LDEBUG( lmclog, "Created kr-complex-line distribution. fwhm: " <<fFWHM<<" line-position: "<<fLinePosition<<" emitted-peak: "<<fEmittedPeak);
+
+        LDEBUG( lmclog, "Created kr-complex-line distribution. fwhm: " <<fFWHM<<" line-position: "<<fLinePosition<<" emitted-peak: "<<fEmittedPeak<<" scatter-prob: "<<fScatterProbability);
         LDEBUG( lmclog, "amplitudes: {" <<amp_log_string<<"} scatter-probabilities: {"<<scatter_log_string<<"}");
     }
 
@@ -357,9 +348,9 @@ namespace locust
         return generate_from_cdf(u, fEnergyLossInterpolator[gas_index], fEnergyLossAccelerator[gas_index]);
     }
 
-    int KrComplexLineDistribution::generate_nscatters(std::string &gas_species)
+    int KrComplexLineDistribution::generate_nscatters()
     {
-        return fGeometric[fGasIndex[gas_species]](*fRNEngine);
+        return fGeometric(*fRNEngine);
     }
 
     double KrComplexLineDistribution::Generate()
@@ -369,18 +360,22 @@ namespace locust
 
         if(fEmittedPeak == "lorentzian")
             generated_energy += fLorentzian(*fRNEngine);
-        if(fEmittedPeak == "shake")
+        else if(fEmittedPeak == "shake")
             generated_energy += generate_shake();
-        
+        else if(fEmittedPeak == "dirac" || fEmittedPeak == "fixed") //just to be explicit
+            generated_energy += 0;
+
         //choose gas species
         std::string gas_species = generate_gas_species();
         
         //generate number of scatters
-        int nScatters = generate_nscatters(gas_species);
+        int nScatters = generate_nscatters();
 
-        //calculate n missed tracks on energy loss
         for(int i=0; i < nScatters; ++i)
+        {
+            gas_species = generate_gas_species();
             generated_energy -= generate_energy_loss(gas_species);
+        }
 
         //include gaussian smearing from finite detector resolution
         generated_energy += fNormal(*fRNEngine);
diff --git a/Source/Distributions/LMCKrComplexLineDistribution.hh b/Source/Distributions/LMCKrComplexLineDistribution.hh
index 1640ba3d..6a6d7d6e 100644
--- a/Source/Distributions/LMCKrComplexLineDistribution.hh
+++ b/Source/Distributions/LMCKrComplexLineDistribution.hh
@@ -41,7 +41,7 @@ namespace locust
             double fFWHM;
             double fLinePosition;
             std::vector<double> fAmplitude;
-            std::vector<double> fScatterProbability;
+            double fScatterProbability;
             std::string fEmittedPeak;
 
             unsigned fNPointsSELA;
@@ -70,7 +70,7 @@ namespace locust
             std::normal_distribution<double> fNormal;
             std::cauchy_distribution<double> fLorentzian;
 
-            std::vector<std::geometric_distribution<int>> fGeometric;
+            std::geometric_distribution<int> fGeometric;
 
             void read_shake_data();
             std::valarray<double> nprime(const double &E_b, const std::valarray<double> &W);
@@ -93,7 +93,7 @@ namespace locust
             double generate_shake();
             std::string generate_gas_species();
             double generate_energy_loss(std::string gas_species);
-            int generate_nscatters(std::string &gas_species);
+            int generate_nscatters();
 
             std::vector<std::vector<double>> transpose_vector(const std::vector<std::vector<double>> aVector);
             std::valarray<double> linspace(double a, double b, unsigned N);
diff --git a/Source/Generators/LMCFakeTrackSignalGenerator.cc b/Source/Generators/LMCFakeTrackSignalGenerator.cc
index cbb16f7e..1cf0751e 100644
--- a/Source/Generators/LMCFakeTrackSignalGenerator.cc
+++ b/Source/Generators/LMCFakeTrackSignalGenerator.cc
@@ -36,7 +36,6 @@ namespace locust
         fStartPitchMax( 90. ),
         fPitchMin( 0. ),
         fLO_frequency( 0. ),
-        fTrackLengthMean( 0. ),
         fNTracksMean( 0. ),
         fBField(1.0),
         fRandomSeed(0),
@@ -45,7 +44,7 @@ namespace locust
         fRandomEngine(0),
         fHydrogenFraction(1),
         fTrapLength(0.1784),  //Phase II harmonic trap L0 (A. Ashtari Esfahani et al.- Phys. Rev. C 99, 055501 )
-        fRoot_filename("LocustEvent.root"),
+        fRootFilename("LocustEvent.root"),
         fSlope( 0. ),
         fPitch( 0. ),
         fTrackLength( 0. ),
@@ -106,6 +105,31 @@ namespace locust
             fSlopeDistribution = fDistributionInterface.get_dist(default_setting);
         }
 
+        if( aParam.has( "track-length" ) )
+        {
+            fTrackLengthDistribution = fDistributionInterface.get_dist(aParam["track-length"].as_node());
+        }
+        else
+        {
+            LWARN( lmclog, "Using default distribution: Track Length = 1e-4 ");
+            scarab::param_node default_setting;
+            default_setting.add("name","dirac");
+            default_setting.add("value","1e-4");
+            fTrackLengthDistribution = fDistributionInterface.get_dist(default_setting);
+        }
+
+        if( aParam.has( "z0" ) )
+        {
+            fz0Distribution = fDistributionInterface.get_dist(aParam["z0"].as_node());
+        }
+        else
+        {
+            LWARN( lmclog, "Using default distribution: z0 = 0 ");
+            scarab::param_node default_setting;
+            default_setting.add("name","dirac");
+            fz0Distribution = fDistributionInterface.get_dist(default_setting);
+        }
+
         if( aParam.has( "signal-power" ) )
             SetSignalPower( aParam.get_value< double >( "signal-power", fSignalPower ) );
 
@@ -130,8 +154,6 @@ namespace locust
         if( aParam.has( "lo-frequency" ) )
             SetFrequency( aParam.get_value< double >( "lo-frequency", fLO_frequency ) );
 
-        if( aParam.has( "track-length-mean" ) )
-            SetTrackLengthMean( aParam.get_value< double >( "track-length-mean", fTrackLengthMean ) );
 
         if (aParam.has( "ntracks-mean") )
             SetNTracksMean( aParam.get_value< double >( "ntracks-mean",fNTracksMean) );
@@ -148,6 +170,9 @@ namespace locust
         if (aParam.has( "pitch-correction") )
             SetPitchCorrection(  aParam.get_value< bool >( "pitch-correction", fPitchCorrection) );
 
+        if (aParam.has( "trap-length") )
+            SetTrapLength(  aParam.get_value< double >( "trap-length", fTrapLength) );
+
         if (aParam.has( "hydrogen-fraction") )
         {
             SetHydrogenFraction(  aParam.get_value< int >( "hydrogen-fraction",fHydrogenFraction) );
@@ -158,7 +183,7 @@ namespace locust
 
         if( aParam.has( "root-filename" ) )
         {
-            fRoot_filename = aParam["root-filename"]().as_string();
+            fRootFilename = aParam["root-filename"]().as_string();
         }
 
         if(fRandomSeed)
@@ -282,17 +307,6 @@ namespace locust
         return;
     }
 
-    double FakeTrackSignalGenerator::GetTrackLengthMean() const
-    {
-        return fTrackLengthMean;
-    }
-
-    void FakeTrackSignalGenerator::SetTrackLengthMean( double aTrackLengthMean )
-    {
-        fTrackLengthMean = aTrackLengthMean;
-        return;
-    }
-
     double FakeTrackSignalGenerator::GetStartTimeMin() const
     {
         return fStartTimeMin;
@@ -394,6 +408,17 @@ namespace locust
         return;
     }
 
+    double FakeTrackSignalGenerator::GetTrapLength() const
+    {
+        return fTrapLength;
+    }
+
+    void FakeTrackSignalGenerator::SetTrapLength( double aTrapLength )
+    {
+        fTrapLength = aTrapLength;
+        return;
+    }
+
 
     Signal::State FakeTrackSignalGenerator::GetDomain() const
     {
@@ -598,7 +623,6 @@ namespace locust
         double theta_scatter;
         const double deg_to_rad = LMCConst::Pi() / 180.;
 
-        std::exponential_distribution<double> tracklength_distribution(1./fTrackLengthMean);
         std::uniform_real_distribution<double> starttime_distribution(fStartTimeMin,fStartTimeMax);
         std::uniform_real_distribution<double> startpitch_distribution(cos(fStartPitchMin * deg_to_rad),cos(fStartPitchMax * deg_to_rad));
         std::uniform_real_distribution<double> dist(0,1);
@@ -623,10 +647,18 @@ namespace locust
                 fStartFrequency = rel_cyc(fStartEnergyDistribution->Generate(), fBField);
             }
 
-            fPitch = acos(startpitch_distribution(fRandomEngine));
+            do
+            {
+                fPitch = acos(startpitch_distribution(fRandomEngine));
+                double z0 = fz0Distribution->Generate();
+                fPitch = GetPitchAngleZ(fPitch, GetBField(z0), fBField);
+
+            } while(fPitch < fStartPitchMin * deg_to_rad );
+
             aTrack.StartTime = fStartTime;
             aTrack.StartFrequency = fStartFrequency;
         }
+
        else
        {
             fStartTime = fEndTime + 0.;  // old track endtime + margin=0.
@@ -653,7 +685,7 @@ namespace locust
         }
 
         fSlope = fSlopeDistribution->Generate();
-        fTrackLength = tracklength_distribution(fRandomEngine);
+        fTrackLength = fTrackLengthDistribution->Generate();
         fEndTime = fStartTime + fTrackLength;  // reset endtime.
         aTrack.Slope = fSlope;
         aTrack.TrackLength = fTrackLength;
@@ -686,32 +718,12 @@ namespace locust
     }
 
 
-    void WriteRootFile(Event* anEvent, TFile* hfile)
-    {
-    	char buffer[100];
-        int n=sprintf(buffer, "Event_%d", anEvent->fEventID);
-    	char* treename = buffer;
-
-        TTree *aTree = new TTree(treename,"Locust Tree");
-        aTree->Branch("EventID", &anEvent->fEventID, "EventID/I");
-        aTree->Branch("ntracks", &anEvent->fNTracks, "ntracks/I");
-        aTree->Branch("StartFrequencies", "std::vector<double>", &anEvent->fStartFrequencies);
-        aTree->Branch("StartTimes", "std::vector<double>", &anEvent->fStartTimes);
-        aTree->Branch("EndTimes", "std::vector<double>", &anEvent->fEndTimes);
-        aTree->Branch("TrackLengths", "std::vector<double>", &anEvent->fTrackLengths);
-        aTree->Branch("Slopes", "std::vector<double>", &anEvent->fSlopes);
-        aTree->Branch("LOFrequency", &anEvent->fLOFrequency, "LOFrequency/D");
-        aTree->Branch("RandomSeed", &anEvent->fRandomSeed, "RandomSeed/I");
-        aTree->Branch("TrackPower", "std::vector<double>", &anEvent->fTrackPowers);
-        aTree->Branch("PitchAngles", "std::vector<double>", &anEvent->fPitchAngles);
-        aTree->Fill();
-        aTree->Write();
-        delete aTree;
-    }
-
     bool FakeTrackSignalGenerator::DoGenerateTime( Signal* aSignal )
     {
-        TFile* hfile = new TFile(fRoot_filename.c_str(),"RECREATE");
+
+    	FileWriter* aRootTreeWriter = RootTreeWriter::get_instance();
+    	aRootTreeWriter->SetFilename(fRootFilename);
+        aRootTreeWriter->OpenFile("RECREATE");
 
         const unsigned nChannels = fNChannels;
         const double tLocustStep = 1./aSignal->DecimationFactor()/(fAcquisitionRate*1.e6);
@@ -766,11 +778,12 @@ namespace locust
                 }
 
             } //track loop
-            WriteRootFile(anEvent, hfile);
+
+            aRootTreeWriter->WriteEvent(anEvent);
             delete anEvent;
         } //event loop
 
-        hfile->Close();
+        aRootTreeWriter->CloseFile();
         return true;
     }
 
diff --git a/Source/Generators/LMCFakeTrackSignalGenerator.hh b/Source/Generators/LMCFakeTrackSignalGenerator.hh
index 83780f32..71903701 100644
--- a/Source/Generators/LMCFakeTrackSignalGenerator.hh
+++ b/Source/Generators/LMCFakeTrackSignalGenerator.hh
@@ -8,12 +8,10 @@
 #ifndef LMCFAKETRACKSIGNALGENERATOR_HH_
 #define LMCFAKETRACKSIGNALGENERATOR_HH_
 
-#include "TFile.h"  // order of includes matters.
-#include "TTree.h"  // include these first.
-
 #include "LMCGenerator.hh"
-#include "LMCRunLengthCalculator.hh"
 #include "LMCEvent.hh"
+#include "LMCRunParameters.hh"
+#include "LMCRootTreeWriter.hh"
 #include "LMCDistributionInterface.hh"
 
 #include <gsl/gsl_errno.h>
@@ -30,7 +28,6 @@ namespace scarab
 
 namespace locust
 {
-    class Digitizer;
 
     /*!
       @class FakeTrackSignalGenerator
@@ -54,8 +51,10 @@ namespace locust
       - "magnetic-field": double -- Magnetic field used to convert from frequency to energy (for jumpsize) (T).
       - "n-events": int -- Number of events per simulation, spaced by 0.5 ms (hardcoded).
       - "random-seed": integer -- integer seed for random number generator for above pdfs, if set to 0 random_device will be used.
-      - "root-filename": str -- Name of root file containing event information to be written at output. 
       - "pitch-correction": bool -- Flag to switch pitch angle corrections on [default] or off.
+      - "root-filename": string -- Name of output Root file.  This can have the same name as other
+          	 	 	 	  generators' output Root files, in which case all of the Root objects will
+          	 	 	 	  be written to the same output file.
 
 
 */
@@ -84,9 +83,6 @@ namespace locust
             double GetPitchMin() const;
             void SetPitchMin( double aPitchMin );
 
-            double GetTrackLengthMean() const;
-            void SetTrackLengthMean( double aTrackLengthMean );
-
             double GetStartVPhase() const;
             void SetStartVPhase( double aPhase );
 
@@ -114,6 +110,9 @@ namespace locust
             bool GetPitchCorrection() const;
             void SetPitchCorrection(  bool aPitchCorrection );
 
+            double GetTrapLength() const;
+            void SetTrapLength(  double aTrapLength );
+
 
             Signal::State GetDomain() const;
             void SetDomain( Signal::State aDomain );
@@ -158,24 +157,25 @@ namespace locust
             std::shared_ptr< BaseDistribution> fStartEnergyDistribution;
             std::shared_ptr< BaseDistribution> fStartFrequencyDistribution;
             std::shared_ptr< BaseDistribution> fSlopeDistribution;
+            std::shared_ptr< BaseDistribution> fTrackLengthDistribution;
+            std::shared_ptr< BaseDistribution> fz0Distribution;
             double fStartTimeMax;
             double fStartTimeMin;
             double fStartPitchMin;
             double fStartPitchMax;
             double fPitchMin;
             double fLO_frequency;
-            double fTrackLengthMean;
             double fNTracksMean;
             double fBField;
             int fRandomSeed;
             int fNEvents;
             bool fPitchCorrection;
             double fHydrogenFraction;
-            std::string fRoot_filename;
+            std::string fRootFilename;
             std::default_random_engine fRandomEngine;
             std::vector<gsl_spline*> fInterpolators;
             std::vector<gsl_interp_accel*> fAccelerators;
-            const double fTrapLength;
+            double fTrapLength;
             bool fUseEnergyDistribution;
             bool fUseFrequencyDistribution;
 
diff --git a/Source/Generators/LMCGaussianNoiseGenerator.cc b/Source/Generators/LMCGaussianNoiseGenerator.cc
index 02bf475a..66d442b4 100644
--- a/Source/Generators/LMCGaussianNoiseGenerator.cc
+++ b/Source/Generators/LMCGaussianNoiseGenerator.cc
@@ -25,8 +25,10 @@ namespace locust
         fDoGenerateFunc( &GaussianNoiseGenerator::DoGenerateFreq ),
         fMean( 0. ),
         fSigma( 1. ),
-        fRandomSeed(0),
-        fNormDist( fMean, fSigma )
+        fRandomSeed( 0 ),
+        fNormDist( fMean, fSigma ),
+		fWriteRootTree( false ),
+        fRootFilename( "LocustNoise.root")
     {
         fRequiredSignalState = Signal::kFreq;
     }
@@ -40,7 +42,6 @@ namespace locust
 
         if( aParam.has( "noise-floor-psd" ) )
         {
-            //tSigma = sqrt( aParam->get_value< double >( "noise-floor" ) * fAcquisitionRate * 1.e6);  // sampling rate
             fSigma = sqrt( aParam["noise-floor-psd"]().as_double() );
             if(aParam.has("noise-temperature"))
             {
@@ -59,6 +60,19 @@ namespace locust
             return false;
         }
 
+        if (aParam.has( "write-root-tree" ))
+        {
+            if (aParam["write-root-tree"]().as_bool())
+            {
+            	fWriteRootTree = true;
+            }
+        }
+
+    	if( aParam.has( "root-filename" ) )
+        {
+            fRootFilename = aParam["root-filename"]().as_string();
+        }
+
 
         if (aParam.has( "random-seed") )
         {
@@ -123,6 +137,7 @@ namespace locust
         fSigma = aSigma;
         return;
     }
+
     int GaussianNoiseGenerator::GetRandomSeed() const
     {
         return fRandomSeed;
@@ -159,6 +174,22 @@ namespace locust
     }
 
 
+    bool GaussianNoiseGenerator::WriteRootTree()
+    {
+		#ifdef LMCFILEWRITER_HH_
+    	FileWriter* aRootTreeWriter = RootTreeWriter::get_instance();
+    	aRootTreeWriter->SetFilename(fRootFilename);
+    	aRootTreeWriter->OpenFile("UPDATE");
+        RunParameters* aRunParameter = new RunParameters();
+        aRunParameter->fNoise = fSigma*fSigma;
+        aRootTreeWriter->WriteRunParameters(aRunParameter, "Noise");
+        aRootTreeWriter->CloseFile();
+        delete aRunParameter;
+		#endif
+        return true;
+    }
+
+
     bool GaussianNoiseGenerator::DoGenerate( Signal* aSignal )
     {
         return (this->*fDoGenerateFunc)( aSignal );
@@ -180,10 +211,10 @@ namespace locust
         std::default_random_engine generator(random_seed_val);
 
     	SetMeanAndSigma( fMean, fSigma, fSigma * sqrt(fAcquisitionRate * 1.e6) );
+    	if (fWriteRootTree) WriteRootTree();
 
         double gain=1.;
         const unsigned nchannels = fNChannels;
-        //double phi = 0.;  // voltage phase
         double mag_r = 0.;  // voltage mag
         double mag_i = 0.;
         const double tResistance = 50.;
@@ -192,7 +223,6 @@ namespace locust
         {
             for( unsigned index = 0; index < aSignal->TimeSize(); ++index )
             {
-                //phi = fUniDist( fRNG );
                 mag_r = fNormDist( generator ) * sqrt(0.5);
                 mag_i = fNormDist( generator ) * sqrt(0.5);
                 aSignal->SignalTimeComplex()[ch*aSignal->TimeSize() + index][0] += gain*sqrt(tResistance)* mag_r;
diff --git a/Source/Generators/LMCGaussianNoiseGenerator.hh b/Source/Generators/LMCGaussianNoiseGenerator.hh
index 08351fb6..5feacc56 100644
--- a/Source/Generators/LMCGaussianNoiseGenerator.hh
+++ b/Source/Generators/LMCGaussianNoiseGenerator.hh
@@ -12,6 +12,10 @@
 #include "LMCRunLengthCalculator.hh"
 #include "LMCConst.hh"
 
+#ifdef LMCFILEWRITER_HH_
+    #include "LMCRunParameters.hh"
+    #include "LMCRootTreeWriter.hh"
+#endif
 
 #include <random>
 
@@ -30,12 +34,17 @@ namespace locust
      Configuration name: "gaussian-noise"
 
      Available configuration options:
-     - "noise-floor" -- Measured noise floor of system in W/Hz.  Overrides "sigma"
-     - "acquisition-rate" -- Digitization acquisition rate in MHz; Used in the calculation of sigma from the noise floor. Default is 100 MHz.
-     - "mean": double -- Mean of the Gaussian noise (usually remains at 0)
-     - "sigma": double -- Standard deviation of the Gaussian noise.  Ignored if "noise-floor" is present.
+     - "noise-floor-psd": double -- Noise power in W/Hz.
+     - "noise-temperature": double -- Noise temperature in K.
      - "domain": string -- Determines whether the noise is generated in the time or frequency domain
                            Available options: "time" and "freq" [default]
+     - "write-root-tree": boolean -- Flag to control whether or not value of noise power is written
+     	 	 			  to output Root file as a LMCRunParameter.
+     - "random-seed": int -- Random seed used to generate random noise.  If this is omitted then the
+     	 	 	 	 	  noise spectrum is reproducible.
+     - "root-filename": string -- Name of output Root file.  This can have the same name as other
+          	 	 	 	  generators' output Root files, in which case all of the Root objects will
+          	 	 	 	  be written to the same output file.
 
     */
     class GaussianNoiseGenerator : public Generator
@@ -67,6 +76,7 @@ namespace locust
             void SetDomain( Signal::State aDomain );
 
         private:
+            bool WriteRootTree();
             bool DoGenerate( Signal* aSignal );
 
             bool DoGenerateTime( Signal* aSignal );
@@ -77,6 +87,8 @@ namespace locust
             double fMean;
             double fSigma;
             int fRandomSeed;
+            bool fWriteRootTree;
+            std::string fRootFilename;
 
             mutable std::normal_distribution< double > fNormDist;
     };
diff --git a/Source/Generators/LMCTestSignalGenerator.cc b/Source/Generators/LMCTestSignalGenerator.cc
index b12f32fb..1e202647 100644
--- a/Source/Generators/LMCTestSignalGenerator.cc
+++ b/Source/Generators/LMCTestSignalGenerator.cc
@@ -25,7 +25,11 @@ namespace locust
         fDoGenerateFunc( &TestSignalGenerator::DoGenerateTime ),
         fLO_frequency( 20.05e9 ),
         fRF_frequency( 20.1e9 ),
-        fAmplitude( 5.e-8 )
+        fAmplitude( 5.e-8 ),
+		fMixingProduct( false ),
+    	fWriteRootHisto( false ),
+		fWriteRootGraph( false ),
+		fRootFilename( "LocustTestSignal.root")
     {
         fRequiredSignalState = Signal::kTime;
     }
@@ -39,18 +43,37 @@ namespace locust
     {
         if( aParam.has( "rf-frequency" ) )
         {
-        SetRFFrequency( aParam.get_value< double >( "rf-frequency", fRF_frequency ) );
+            SetRFFrequency( aParam.get_value< double >( "rf-frequency", fRF_frequency ) );
         }
 
         if( aParam.has( "lo-frequency" ) )
         {
-        SetLOFrequency( aParam.get_value< double >( "lo-frequency", fLO_frequency ) );
+            SetLOFrequency( aParam.get_value< double >( "lo-frequency", fLO_frequency ) );
         }
 
-
         if( aParam.has( "amplitude" ) )
         {
-        SetAmplitude( aParam.get_value< double >( "amplitude", fAmplitude ) );
+            SetAmplitude( aParam.get_value< double >( "amplitude", fAmplitude ) );
+        }
+
+        if( aParam.has( "mixing-product" ) )
+        {
+        	SetMixingProduct( aParam.get_value< bool >( "mixing-product", fMixingProduct ));
+        }
+
+        if( aParam.has( "write-root-histo" ) )
+        {
+        	fWriteRootHisto = aParam.get_value< bool >( "write-root-histo", fWriteRootHisto );
+        }
+
+        if( aParam.has( "write-root-graph" ) )
+        {
+        	fWriteRootGraph = aParam.get_value< bool >( "write-root-graph", fWriteRootGraph );
+        }
+
+    	if( aParam.has( "root-filename" ) )
+        {
+            fRootFilename = aParam["root-filename"]().as_string();
         }
 
 
@@ -120,6 +143,16 @@ namespace locust
         return;
     }
 
+    bool TestSignalGenerator::GetMixingProduct() const
+    {
+        return fMixingProduct;
+    }
+
+    void TestSignalGenerator::SetMixingProduct( bool aMixingProduct )
+    {
+        fMixingProduct = aMixingProduct;
+        return;
+    }
 
 
     Signal::State TestSignalGenerator::GetDomain() const
@@ -146,6 +179,39 @@ namespace locust
         return;
     }
 
+    bool TestSignalGenerator::WriteRootHisto()
+    {
+		#ifdef LMCFILEWRITER_HH_
+    	FileWriter* aRootHistoWriter = RootHistoWriter::get_instance();
+    	aRootHistoWriter->SetFilename(fRootFilename);
+    	aRootHistoWriter->OpenFile("UPDATE");
+    	TH1D* aHisto = new TH1D("testhisto", "histotitle", 200, 0., 400.);
+    	for (unsigned i=0; i<200; i++)
+    	{
+    		aHisto->SetBinContent(i+1, (double)i);
+    	}
+    	aRootHistoWriter->Write1DHisto(aHisto);
+        aRootHistoWriter->CloseFile();
+		#endif
+        return true;
+    }
+
+    bool TestSignalGenerator::WriteRootGraph()
+    {
+		#ifdef LMCFILEWRITER_HH_
+    	FileWriter* aRootGraphWriter = RootGraphWriter::get_instance();
+    	aRootGraphWriter->SetFilename(fRootFilename);
+    	aRootGraphWriter->OpenFile("UPDATE");
+    	std::vector<double> xVector{1,2,3,4,5};
+    	std::vector<double> yVector{1,2,3,4,5};
+    	aRootGraphWriter->WriteVectorGraph(xVector, yVector);
+        aRootGraphWriter->CloseFile();
+		#endif
+        return true;
+    }
+
+
+
 
     bool TestSignalGenerator::DoGenerate( Signal* aSignal )
     {
@@ -155,6 +221,9 @@ namespace locust
     bool TestSignalGenerator::DoGenerateTime( Signal* aSignal )
     {
 
+    	if (fWriteRootHisto) WriteRootHisto();
+    	if (fWriteRootGraph) WriteRootGraph();
+
         const unsigned nchannels = fNChannels;
 
         double LO_phase = 0.;
@@ -168,8 +237,16 @@ namespace locust
                 LO_phase = 2.*LMCConst::Pi()*fLO_frequency*(double)index/aSignal->DecimationFactor()/(fAcquisitionRate*1.e6);
                 voltage_phase = 2.*LMCConst::Pi()*fRF_frequency*(double)index/aSignal->DecimationFactor()/(fAcquisitionRate*1.e6);
 
-                aSignal->LongSignalTimeComplex()[ch*aSignal->TimeSize()*aSignal->DecimationFactor() + index][0] += sqrt(50.)*fAmplitude*cos(voltage_phase-LO_phase);
-                aSignal->LongSignalTimeComplex()[ch*aSignal->TimeSize()*aSignal->DecimationFactor() + index][1] += sqrt(50.)*fAmplitude*cos(-LMCConst::Pi()/2. + voltage_phase-LO_phase);
+                if (fMixingProduct)  // keep upper and lower sidebands RF+LO and RF-LO
+                {
+        		    aSignal->LongSignalTimeComplex()[ch*aSignal->TimeSize()*aSignal->DecimationFactor() + index][0] +=  2. * sqrt(50.)*fAmplitude*cos(voltage_phase) * sin(LO_phase);
+        		    aSignal->LongSignalTimeComplex()[ch*aSignal->TimeSize()*aSignal->DecimationFactor() + index][1] +=  2. * sqrt(50.)*fAmplitude*cos(voltage_phase) * cos(LO_phase);
+                }
+                else // keep only lower sideband RF-LO
+                {
+                	aSignal->LongSignalTimeComplex()[ch*aSignal->TimeSize()*aSignal->DecimationFactor() + index][0] += sqrt(50.)*fAmplitude*cos(voltage_phase-LO_phase);
+                	aSignal->LongSignalTimeComplex()[ch*aSignal->TimeSize()*aSignal->DecimationFactor() + index][1] += sqrt(50.)*fAmplitude*cos(-LMCConst::Pi()/2. + voltage_phase-LO_phase);
+                }
 
 //printf("signal %d is with acqrate %g, lo %g and rf %g is %g\n", index, fAcquisitionRate, fLO_frequency, fRF_frequency, aSignal->LongSignalTimeComplex()[ch*aSignal->TimeSize()*aSignal->DecimationFactor() + index][0]); getchar();
 
diff --git a/Source/Generators/LMCTestSignalGenerator.hh b/Source/Generators/LMCTestSignalGenerator.hh
index dac41623..c91b4efe 100644
--- a/Source/Generators/LMCTestSignalGenerator.hh
+++ b/Source/Generators/LMCTestSignalGenerator.hh
@@ -11,6 +11,10 @@
 #include "LMCGenerator.hh"
 #include "LMCRunLengthCalculator.hh"
 
+#ifdef LMCFILEWRITER_HH_
+    #include "LMCRootGraphWriter.hh"
+    #include "LMCRootHistoWriter.hh"
+#endif
 
 namespace scarab
 {
@@ -35,8 +39,17 @@ namespace locust
      Available configuration options:
      - "frequency": double -- Frequency of the sine wave.
      - "amplitude": double -- Amplitude of the sine wave.
+     - "mixing-product": bool -- If true, use product of RF and LO for downmixing.  If false,
+     	 calculate lower sideband at RF-LO and replace mixing product with it.  The latter is best used
+     	 for sinusoidal RF signals.  Toggling this parameter is useful to see whether the upper sideband
+     	 RF+LO is downmixing into the measurement band, or not.
+     - "write-root-histo": bool -- Flag to determine whether an example Root histogram is generated
+         and written to file.
+     - "write-root-graph": bool -- Flag to determine whether an example Root TGraph is generated
+         and written to file.
+     - "root-filename": string -- name of output Root file.
      - "domain": string -- Determines whether the sinusoidal test signal is generated in the time 
-            or frequency domain
+     	 or frequency domain
     
      Available options: "time" and "freq" [default]
 
@@ -58,6 +71,8 @@ namespace locust
             double GetLOFrequency() const;
             void SetLOFrequency( double aFrequency );
 
+            bool GetMixingProduct() const;
+            void SetMixingProduct( bool aMixingProduct );
 
             double GetAmplitude() const;
             void SetAmplitude( double aAmplitude );
@@ -65,6 +80,10 @@ namespace locust
             Signal::State GetDomain() const;
             void SetDomain( Signal::State aDomain );
 
+            bool WriteRootHisto();
+            bool WriteRootGraph();
+
+
 
         private:
             bool DoGenerate( Signal* aSignal );
@@ -77,6 +96,10 @@ namespace locust
             double fRF_frequency;
             double fLO_frequency;
             double fAmplitude;
+            bool fMixingProduct;
+            bool fWriteRootHisto;
+            bool fWriteRootGraph;
+            std::string fRootFilename;
 
             
     };
diff --git a/Source/IO/LMCFileWriter.cc b/Source/IO/LMCFileWriter.cc
new file mode 100644
index 00000000..9dbb3b4b
--- /dev/null
+++ b/Source/IO/LMCFileWriter.cc
@@ -0,0 +1,46 @@
+/*
+ * LMCFileWriter.cc
+ *
+ *  Created on: Mar 21, 2020
+ *      Author: pslocum
+ */
+
+#include "LMCFileWriter.hh"
+
+namespace locust
+{
+    FileWriter::FileWriter()
+	{
+	}
+ 
+    FileWriter::~FileWriter() {}
+
+    bool FileWriter::Configure( const scarab::param_node& aParam )
+    {
+        return true;
+    }
+
+    void FileWriter::SetFilename(std::string aFilename)
+    {
+    	fRoot_filename = aFilename;
+    }
+
+    std::string FileWriter::GetFilename()
+    {
+    	return fRoot_filename;
+    }
+
+    void FileWriter::OpenFile(std::string aFileFlag)
+    {
+        fFile = new TFile(fRoot_filename.c_str(), aFileFlag.c_str());
+    }
+
+    void FileWriter::CloseFile()
+    {
+    	fFile->Close();
+    }
+
+
+
+} /* namespace locust */
+
diff --git a/Source/IO/LMCFileWriter.hh b/Source/IO/LMCFileWriter.hh
new file mode 100644
index 00000000..8f134cfb
--- /dev/null
+++ b/Source/IO/LMCFileWriter.hh
@@ -0,0 +1,70 @@
+/*
+ * LMCFileWriter.hh
+ *
+ *  Created on: Mar 21, 2020
+ *      Author: pslocum
+ */
+
+#ifndef LMCFILEWRITER_HH_
+#define LMCFILEWRITER_HH_
+
+#include "TFile.h"  // order of includes matters.
+#include "TTree.h"  // include these first.
+#include "TH1.h"
+#include "TH2.h"
+#include "TGraph.h"
+#include "LMCEvent.hh"
+#include "LMCRunParameters.hh"
+
+#include "param.hh"
+#include "logger.hh"
+#include "singleton.hh"
+
+namespace locust
+{
+ /*!
+ @class FileWriter
+ @author P. Slocum
+ @brief Base class to contain LMCFileWriter subclasses, such as LMCRootTreeFileWriter.
+ @details
+ Available configuration options:
+ No input parameters
+ */
+
+//    class FileWriter : public scarab::singleton< locust::FileWriter >
+    class FileWriter
+    {
+    	protected:
+    	FileWriter();
+    	virtual ~FileWriter();
+        virtual bool Configure( const scarab::param_node& aNode );
+
+    	public:
+        virtual double GetTestVar() {};
+        virtual void SetTestVar(double aValue) {};
+        virtual void WriteEvent(Event* anEvent) {};
+        virtual void WriteRunParameters( RunParameters* aRunParameter, const char* aParameterName ) {};
+        virtual void Write1DHisto(TH1D* aHisto) {};
+        virtual void Write2DHisto(TH2D* aHisto) {};
+        virtual void WriteVector1DHisto(std::vector<double> aVector, double xmin, double xmax) {};
+        virtual void Write2DGraph(TGraph* aGraph) {};
+        virtual void WriteVectorGraph(std::vector<double> xVector, std::vector<double> yVector) {};
+
+        virtual void OpenFile(std::string aFileFlag);
+        virtual void CloseFile();
+
+        virtual void SetFilename(std::string aFilename);
+        virtual std::string GetFilename();
+
+
+    	private:
+        TFile* fFile;
+        std::string fRoot_filename;
+
+
+    };
+
+
+} /* namespace locust */
+
+#endif /* LMCFILEWRITER_HH_ */
diff --git a/Source/IO/LMCRootGraphWriter.cc b/Source/IO/LMCRootGraphWriter.cc
new file mode 100644
index 00000000..9f5aa600
--- /dev/null
+++ b/Source/IO/LMCRootGraphWriter.cc
@@ -0,0 +1,69 @@
+/*
+ * LMCRootGraphWriter.cc
+ *
+ *  Created on: Mar. 21, 2020
+ *      Author: pslocum
+ */
+
+#include "LMCRootGraphWriter.hh"
+#include "logger.hh"
+
+
+
+namespace locust
+{
+    LOGGER( lmclog, "RootGraphWriter" );
+
+
+    RootGraphWriter::RootGraphWriter()
+    {
+    }
+
+    RootGraphWriter::~RootGraphWriter()
+    {
+    }
+
+    bool RootGraphWriter::Configure( const scarab::param_node& aParam )
+    {
+
+    	if( !FileWriter::Configure(aParam))
+    	{
+    		LERROR(lmclog,"Error configuring FileWriter class from RootGraphWriter child class");
+    	}
+
+    	return true;
+    }
+
+    void RootGraphWriter::Write2DGraph(TGraph* aGraph)
+    {
+    	aGraph->Write();
+    }
+
+    void RootGraphWriter::WriteVectorGraph(std::vector<double> xVector, std::vector<double> yVector)
+    {
+    	int n = xVector.size();
+    	double* xArray = new double[n];
+    	double* yArray = new double [n];
+
+    	unsigned count = 0;
+    	for (auto it = xVector.begin(); it!=xVector.end(); ++it)
+    	{
+    	    xArray[count] = *it;
+    	    count ++;
+    	}
+    	count = 0;
+    	for (auto it = yVector.begin(); it!=yVector.end(); ++it)
+    	{
+    	    yArray[count] = *it;
+    	    count ++;
+    	}
+
+    	TGraph* aGraph = new TGraph(n, xArray, yArray);
+    	aGraph->Write();
+
+
+    }
+
+
+
+} /* namespace locust */
diff --git a/Source/IO/LMCRootGraphWriter.hh b/Source/IO/LMCRootGraphWriter.hh
new file mode 100644
index 00000000..d76f759e
--- /dev/null
+++ b/Source/IO/LMCRootGraphWriter.hh
@@ -0,0 +1,48 @@
+/*
+ * LMCRootGraphWriter.hh
+ *
+ *  Created on: Jul. 10, 2020
+ *      Author: pslocum
+ */
+
+#ifndef LMCROOTGRAPHWRITER_HH_
+#define LMCROOTGRAPHWRITER_HH_
+
+#include "LMCFileWriter.hh"
+#include <string>
+#include <vector>
+
+
+namespace locust
+{
+ /*!
+ @class RootGraphWriter
+ @author P. Slocum
+ @brief Derived class to configure the Root TGraph that will be written to file.
+ @details
+ Available configuration options:
+ No input parameters
+ */
+    class RootGraphWriter : public FileWriter, public scarab::singleton< locust::RootGraphWriter >
+    {
+
+	public:
+
+    	RootGraphWriter();
+        virtual ~RootGraphWriter();
+
+        virtual bool Configure( const scarab::param_node& aNode );
+        allow_singleton_access( RootGraphWriter );
+
+    private:
+
+        void Write2DGraph(TGraph* aGraph);
+        void WriteVectorGraph(std::vector<double> xVector, std::vector<double> yVector);
+
+
+    };
+
+
+} /* namespace locust */
+
+#endif /* LMCROOTGRAPHWRITER_HH_ */
diff --git a/Source/IO/LMCRootHistoWriter.cc b/Source/IO/LMCRootHistoWriter.cc
new file mode 100644
index 00000000..cf66c170
--- /dev/null
+++ b/Source/IO/LMCRootHistoWriter.cc
@@ -0,0 +1,60 @@
+/*
+ * LMCRootHistoWriter.cc
+ *
+ *  Created on: Mar. 21, 2020
+ *      Author: pslocum
+ */
+
+#include "LMCRootHistoWriter.hh"
+#include "logger.hh"
+
+
+
+namespace locust
+{
+    LOGGER( lmclog, "RootHistoWriter" );
+
+
+    RootHistoWriter::RootHistoWriter()
+    {
+    }
+
+    RootHistoWriter::~RootHistoWriter()
+    {
+    }
+
+    bool RootHistoWriter::Configure( const scarab::param_node& aParam )
+    {
+
+
+    	if( !FileWriter::Configure(aParam))
+    	{
+    		LERROR(lmclog,"Error configuring FileWriter class from RootHistoWriter child class");
+    	}
+
+    	return true;
+    }
+
+    void RootHistoWriter::Write1DHisto(TH1D* aHisto)
+    {
+    	aHisto->Write();
+    }
+
+    void RootHistoWriter::Write2DHisto(TH2D* aHisto)
+    {
+    	aHisto->Write();
+    }
+
+    void RootHistoWriter::WriteVector1DHisto(std::vector<double> aVector, double xmin, double xmax)
+    {
+    	TH1D* aHisto = new TH1D("histo", "title", aVector.size(), xmin, xmax);
+    	unsigned count = 1;
+    	for (auto it = aVector.begin(); it!=aVector.end(); ++it)
+    	{
+    	    aHisto->SetBinContent(count,*it);
+    	    count ++;
+    	}
+    	aHisto->Write();
+    }
+
+} /* namespace locust */
diff --git a/Source/IO/LMCRootHistoWriter.hh b/Source/IO/LMCRootHistoWriter.hh
new file mode 100644
index 00000000..99a878f9
--- /dev/null
+++ b/Source/IO/LMCRootHistoWriter.hh
@@ -0,0 +1,49 @@
+/*
+ * LMCRootHistoWriter.hh
+ *
+ *  Created on: Jul. 8, 2020
+ *      Author: pslocum
+ */
+
+#ifndef LMCROOTHISTOWRITER_HH_
+#define LMCROOTHISTOWRITER_HH_
+
+#include "LMCFileWriter.hh"
+#include <string>
+#include <vector>
+
+
+namespace locust
+{
+ /*!
+ @class RootHistoWriter
+ @author P. Slocum
+ @brief Derived class to configure the Root histogram that will be written to file.
+ @details
+ Available configuration options:
+ No input parameters
+ */
+    class RootHistoWriter : public FileWriter, public scarab::singleton< locust::RootHistoWriter >
+    {
+
+	public:
+
+    	RootHistoWriter();
+        virtual ~RootHistoWriter();
+
+        virtual bool Configure( const scarab::param_node& aNode );
+        allow_singleton_access( RootHistoWriter );
+
+        private:
+
+        virtual void Write1DHisto(TH1D* aHisto);
+        virtual void Write2DHisto(TH2D* aHisto);
+        virtual void WriteVector1DHisto(std::vector<double> aVector, double xmin, double xmax);
+
+
+    };
+
+
+} /* namespace locust */
+
+#endif /* LMCROOTHISTOWRITER_HH_ */
diff --git a/Source/IO/LMCRootTreeWriter.cc b/Source/IO/LMCRootTreeWriter.cc
new file mode 100644
index 00000000..7e9304f7
--- /dev/null
+++ b/Source/IO/LMCRootTreeWriter.cc
@@ -0,0 +1,81 @@
+/*
+ * LMCRootTreeWriter.cc
+ *
+ *  Created on: Mar. 21, 2020
+ *      Author: pslocum
+ */
+
+#include "LMCRootTreeWriter.hh"
+#include "logger.hh"
+using std::string;
+
+
+namespace locust
+{
+    LOGGER( lmclog, "RootTreeWriter" );
+
+
+    RootTreeWriter::RootTreeWriter()
+    {
+    }
+
+    RootTreeWriter::~RootTreeWriter()
+    {
+    }
+
+    bool RootTreeWriter::Configure( const scarab::param_node& aParam )
+    {
+
+    	if( !FileWriter::Configure(aParam))
+    	{
+    		LERROR(lmclog,"Error configuring FileWriter class from RootTreeWriter child class");
+    	}
+
+
+    	return true;
+    }
+
+    void RootTreeWriter::WriteRunParameters( RunParameters* aRunParameter, const char* aParameterName )
+    {
+        TTree *aTree = new TTree("Run Parameters","Locust Tree");
+
+        if (aParameterName=="Noise")
+        	{
+        	aTree->Branch("Noise", &aRunParameter->fNoise, "Noise/D");
+        	}
+        if (aParameterName=="LOfrequency")
+        	{
+        	aTree->Branch("LO frequency", &aRunParameter->fLOfrequency, "LOfrequency/D");
+        	}
+
+        aTree->Fill();
+        aTree->Write();
+        delete aTree;
+    }
+
+    void RootTreeWriter::WriteEvent(Event* anEvent)
+    {
+    	char buffer[100];
+        int n=sprintf(buffer, "Event_%d", anEvent->fEventID);
+    	char* treename = buffer;
+
+        TTree *aTree = new TTree(treename,"Locust Tree");
+        aTree->Branch("EventID", &anEvent->fEventID, "EventID/I");
+        aTree->Branch("ntracks", &anEvent->fNTracks, "ntracks/I");
+        aTree->Branch("StartFrequencies", "std::vector<double>", &anEvent->fStartFrequencies);
+        aTree->Branch("StartTimes", "std::vector<double>", &anEvent->fStartTimes);
+        aTree->Branch("EndTimes", "std::vector<double>", &anEvent->fEndTimes);
+        aTree->Branch("TrackLengths", "std::vector<double>", &anEvent->fTrackLengths);
+        aTree->Branch("Slopes", "std::vector<double>", &anEvent->fSlopes);
+        aTree->Branch("LOFrequency", &anEvent->fLOFrequency, "LOFrequency/D");
+        aTree->Branch("RandomSeed", &anEvent->fRandomSeed, "RandomSeed/I");
+        aTree->Branch("TrackPower", "std::vector<double>", &anEvent->fTrackPowers);
+        aTree->Branch("PitchAngles", "std::vector<double>", &anEvent->fPitchAngles);
+        aTree->Fill();
+        aTree->Write();
+        delete aTree;
+    }
+
+
+
+} /* namespace locust */
diff --git a/Source/IO/LMCRootTreeWriter.hh b/Source/IO/LMCRootTreeWriter.hh
new file mode 100644
index 00000000..ee642069
--- /dev/null
+++ b/Source/IO/LMCRootTreeWriter.hh
@@ -0,0 +1,47 @@
+/*
+ * LMCRootTreeWriter.hh
+ *
+ *  Created on: Mar. 21, 2020
+ *      Author: pslocum
+ */
+
+#ifndef LMCROOTTREEWRITER_HH_
+#define LMCROOTTREEWRITER_HH_
+
+#include "LMCFileWriter.hh"
+#include <string>
+
+
+namespace locust
+{
+ /*!
+ @class RootTreeWriter
+ @author P. Slocum
+ @brief Derived class to configure the Root TTree that will be written to file.
+ @details
+ Available configuration options:
+ No input parameters
+ */
+    class RootTreeWriter : public FileWriter, public scarab::singleton< locust::RootTreeWriter >
+    {
+
+	public:
+
+    	RootTreeWriter();
+        virtual ~RootTreeWriter();
+
+        virtual bool Configure( const scarab::param_node& aNode );
+        allow_singleton_access( RootTreeWriter );
+
+        virtual void WriteEvent(Event* anEvent);
+        virtual void WriteRunParameters( RunParameters* aRunParameter, const char* aParameterName );
+
+        private:
+
+
+    };
+
+
+} /* namespace locust */
+
+#endif /* LMCROOTTREEWRITER_HH_ */
diff --git a/Source/RxComponents/LMCPowerCombiner.cc b/Source/RxComponents/LMCPowerCombiner.cc
index caea1f2d..815de629 100644
--- a/Source/RxComponents/LMCPowerCombiner.cc
+++ b/Source/RxComponents/LMCPowerCombiner.cc
@@ -22,7 +22,8 @@ namespace locust
             fpatchLoss( 0.6 ),
             famplifierLoss( 0.66 ),
             fendPatchLoss( 1.0 ),
-            fjunctionResistance( 0.3 )
+            fjunctionResistance( 0.3 ),
+			fvoltageCheck( false )
 
     {}
     PowerCombiner::~PowerCombiner() {}
@@ -30,6 +31,12 @@ namespace locust
 
     bool PowerCombiner::Configure( const scarab::param_node& aParam )
     {
+
+    	if ( aParam.has( "voltage-check" ) )
+    	{
+    		fvoltageCheck = aParam["voltage-check"]().as_bool();
+    	}
+
         if( aParam.has( "nelements-per-strip" ) )
         {
             fnElementsPerStrip = aParam["nelements-per-strip"]().as_int();
@@ -63,7 +70,6 @@ namespace locust
     }
 
 
-
 	bool PowerCombiner::AddOneVoltageToStripSum(Signal* aSignal, double VoltageFIRSample, double phi_LO, unsigned z_index, unsigned sampleIndex)
 	{
 
@@ -71,6 +77,9 @@ namespace locust
 		aSignal->LongSignalTimeComplex()[sampleIndex][0] += 2.*VoltageFIRSample * sin(phi_LO);
 		aSignal->LongSignalTimeComplex()[sampleIndex][1] += 2.*VoltageFIRSample * cos(phi_LO);
 
+		if ( (sampleIndex%100 < 1) && (fvoltageCheck==true) )
+			LWARN( lmclog, "Voltage " << sampleIndex << " is <" << aSignal->LongSignalTimeComplex()[sampleIndex][1] << ">" );
+
 		return true;
 	}
 
diff --git a/Source/RxComponents/LMCPowerCombiner.hh b/Source/RxComponents/LMCPowerCombiner.hh
index 18ef59ab..51d84ec6 100644
--- a/Source/RxComponents/LMCPowerCombiner.hh
+++ b/Source/RxComponents/LMCPowerCombiner.hh
@@ -67,6 +67,7 @@ namespace locust
             double famplifierLoss;
             double fendPatchLoss;
             double fjunctionResistance;
+            bool fvoltageCheck;
 
 };
 
diff --git a/Source/Transmitters/LMCPlaneWaveTransmitter.cc b/Source/Transmitters/LMCPlaneWaveTransmitter.cc
index ed0b2332..f152c4a5 100644
--- a/Source/Transmitters/LMCPlaneWaveTransmitter.cc
+++ b/Source/Transmitters/LMCPlaneWaveTransmitter.cc
@@ -50,11 +50,19 @@ namespace locust
     void PlaneWaveTransmitter::AddPropagationPhaseDelay(LMCThreeVector pointOfInterest)
     {
 	//Assuming the element strip is always along Z
-	double meanZ = GetMeanofFieldPoints(2);
+
+      // Some previous attempts to get phase delay:
+      
+	/*double meanZ = GetMeanofFieldPoints(2);
+	  double distanceFromCenter = meanZ - pointOfInterest.GetZ();*/
 	/*int z_index = fieldPointIndex%nElementsPerStrip;
     	double stripLength = (nElementsPerStrip-1)*elementSpacing;
     	double distanceFromCenter = stripLength/2. - z_index*elementSpacing;*/
-    	double distanceFromCenter = meanZ - pointOfInterest.GetZ();
+
+      // Take the end of the array as the first place that the phase front of the planewave touches.
+      double endOfStrip = GetFieldPoint(0).GetZ();
+      double distanceFromCenter = pointOfInterest.GetZ() - endOfStrip;
+    
     	double phaseDelay = 2*LMCConst::Pi()*distanceFromCenter*sin(fAOI)*fRF_Frequency/LMCConst::C();
 	Transmitter::AddPropagationPhaseDelay(phaseDelay);
     }