diff --git a/apps/showcases/ChargedParticles/ChargedParticles.cpp b/apps/showcases/ChargedParticles/ChargedParticles.cpp
index a09c4bed91d8b30e588a4c613e51a863cfd488d1..9641c2b268078f6e1fda710ead534ccd26a5fdcc 100644
--- a/apps/showcases/ChargedParticles/ChargedParticles.cpp
+++ b/apps/showcases/ChargedParticles/ChargedParticles.cpp
@@ -278,32 +278,35 @@ void createPlaneSetup(const shared_ptr< mesa_pd::data::ParticleStorage >& ps,
 
 struct ParticleInfo
 {
-   real_t averageVelocity = 0_r;
-   real_t maximumVelocity = 0_r;
-   uint_t numParticles    = 0;
-   real_t maximumHeight   = 0_r;
-   real_t particleVolume  = 0_r;
-   real_t heightOfMass    = 0_r;
+   real_t averageVelocity           = 0_r;
+   real_t ensembledAverageVelocityZ = 0_r;
+   real_t maximumVelocity           = 0_r;
+   uint_t numParticles              = 0;
+   real_t maximumHeight             = 0_r;
+   real_t particleVolume            = 0_r;
+   real_t heightOfMass              = 0_r;
 
    void allReduce()
    {
       walberla::mpi::allReduceInplace(numParticles, walberla::mpi::SUM);
       walberla::mpi::allReduceInplace(averageVelocity, walberla::mpi::SUM);
+      walberla::mpi::allReduceInplace(ensembledAverageVelocityZ, walberla::mpi::SUM);
       walberla::mpi::allReduceInplace(maximumVelocity, walberla::mpi::MAX);
       walberla::mpi::allReduceInplace(maximumHeight, walberla::mpi::MAX);
       walberla::mpi::allReduceInplace(particleVolume, walberla::mpi::SUM);
       walberla::mpi::allReduceInplace(heightOfMass, walberla::mpi::SUM);
 
       averageVelocity /= real_c(numParticles);
+      ensembledAverageVelocityZ /= real_c(numParticles);
       heightOfMass /= particleVolume;
    }
 };
 
 std::ostream& operator<<(std::ostream& os, ParticleInfo const& m)
 {
-   return os << "Particle Info: uAvg = " << m.averageVelocity << ", uMax = " << m.maximumVelocity
-             << ", numParticles = " << m.numParticles << ", zMax = " << m.maximumHeight << ", Vp = " << m.particleVolume
-             << ", zMass = " << m.heightOfMass;
+   return os << "Ensembled Avg Z = " << m.ensembledAverageVelocityZ << ", Particle Info: uAvg = " << m.averageVelocity
+             << ", uMax = " << m.maximumVelocity << ", numParticles = " << m.numParticles
+             << ", zMax = " << m.maximumHeight << ", Vp = " << m.particleVolume << ", zMass = " << m.heightOfMass;
 }
 
 template< typename Accessor_T >
@@ -319,9 +322,11 @@ ParticleInfo evaluateParticleInfo(const Accessor_T& ac)
 
       ++info.numParticles;
       real_t velMagnitude   = ac.getLinearVelocity(i).length();
+      real_t velocityZ      = ac.getLinearVelocity(i)[2];
       real_t particleVolume = ac.getShape(i)->getVolume();
       real_t height         = ac.getPosition(i)[2];
       info.averageVelocity += velMagnitude;
+      info.ensembledAverageVelocityZ += velocityZ;
       info.maximumVelocity = std::max(info.maximumVelocity, velMagnitude);
       info.maximumHeight   = std::max(info.maximumHeight, height);
       info.particleVolume += particleVolume;
@@ -1241,6 +1246,8 @@ int main(int argc, char** argv)
       auto particleInfo = evaluateParticleInfo(*accessor);
       auto fluidInfo = evaluateFluidInfo< BoundaryHandling_T >(blocks, pdfFieldID, boundaryHandlingID);
 
+      if (timeStep % 1000 == 0) { WriteEnsembledVelocityToFile(timeStep, particleInfo); }
+
       if (timeStep % infoSpacing == 0)
       {
          timeloopTiming["Evaluate infos"].start();
diff --git a/apps/showcases/ChargedParticles/postProcessingUtilities.h b/apps/showcases/ChargedParticles/postProcessingUtilities.h
index bb4b841da2c8ddb4d21efe7f4456f889193dc362..4eda4967bfad3c6c7da676aa3973b92b99610c19 100644
--- a/apps/showcases/ChargedParticles/postProcessingUtilities.h
+++ b/apps/showcases/ChargedParticles/postProcessingUtilities.h
@@ -211,6 +211,32 @@ void writeStressesToFile(std::vector<real_t>& hydroStress, std::vector<real_t>&
    }
 }
 
+template< typename ParticleInfo >
+void WriteEnsembledVelocityToFile(const uint_t timestep, ParticleInfo info)
+{
+   WALBERLA_ROOT_SECTION()
+   {
+      std::ofstream file("ensembledVelocity.txt", std::ios_base::app);
+      if (file.is_open())
+      {
+         if (file.tellp() == 0)
+         {
+            file << "timestep"
+                 << ", "
+                 << "ensembledVelocity "
+                 << "\n";
+         }
+         file << timestep << "," << info.ensembledAverageVelocityZ << "\n";
+
+         file.close();
+      }
+      else { std::cerr << "Error opening file for writing.\n"; }
+      file.close();
+   }
+}
+
+
+
 }