diff --git a/docs/lib_guide/lib_guide.tex b/docs/lib_guide/lib_guide.tex
index 9e94f004..20bd583b 100644
--- a/docs/lib_guide/lib_guide.tex
+++ b/docs/lib_guide/lib_guide.tex
@@ -183,7 +183,7 @@ \subsection{Particle Transforms}
 apr.particle_intensities.map(apr,absolute_value_of_intensities,[](const uint16_t &a) { return abs(a); });
 \end{lstlisting}
 \section{Pixel Images and Reconstruction}
-Pixel images are handled using the MeshData class and tiff library is used for input-output.
+Pixel images are handled using the PixelData class and tiff library is used for input-output.
 \subsection{Reconstruction}
 Any particle property, stored using ExtraParticleData class, can be interpolated to reconstruct an image, allowing both a piecewise constant and smooth reconstruction. For examples see Example\_reconstruct\_image.
 \subsection{Sampling particles from an image}
diff --git a/examples/Example_compress_apr.cpp b/examples/Example_compress_apr.cpp
index 028105ea..7aa0ac46 100644
--- a/examples/Example_compress_apr.cpp
+++ b/examples/Example_compress_apr.cpp
@@ -70,7 +70,7 @@ int main(int argc, char **argv) {
     apr.read_apr(options.directory + name + "_compress_apr.h5");
     timer.stop_timer();
 
-    MeshData<uint16_t> img;
+    PixelData<uint16_t> img;
     apr.interp_img(img,apr.particles_intensities);
     std::string output = options.directory + name + "_compress.tif";
     TiffUtils::saveMeshAsTiff(output, img);
diff --git a/examples/Example_compute_gradient.cpp b/examples/Example_compute_gradient.cpp
index c58ca7b9..03524a4b 100644
--- a/examples/Example_compute_gradient.cpp
+++ b/examples/Example_compute_gradient.cpp
@@ -79,7 +79,7 @@ int main(int argc, char **argv) {
     gradient.map(apr,gradient_magnitude,[](const std::vector<float> &a) { return 20.0f*sqrt(pow(a[0], 2.0f) + pow(a[1], 2.0f) + pow(a[2], 2.0f)); });
 
     // write result to image
-    MeshData<float> gradient_magnitude_image;
+    PixelData<float> gradient_magnitude_image;
     apr.interp_img(gradient_magnitude_image,gradient_magnitude);
     //apr.interp_img(gradient_magnitude_image,apr.particles_intensities);
 
@@ -96,9 +96,9 @@ int main(int argc, char **argv) {
     if(options.original_image.size() > 0) {
 
         TiffUtils::TiffInfo inputTiff(options.directory + options.original_image);
-        MeshData<uint16_t> original_image = TiffUtils::getMesh<uint16_t>(inputTiff);
+        PixelData<uint16_t> original_image = TiffUtils::getMesh<uint16_t>(inputTiff);
 
-        std::vector<MeshData<float>> gradient_mesh;
+        std::vector<PixelData<float>> gradient_mesh;
 
         MeshNumerics::compute_gradient(original_image, gradient_mesh);
 
diff --git a/examples/Example_get_apr.cpp b/examples/Example_get_apr.cpp
index 59e6f67b..4b2f0db5 100644
--- a/examples/Example_get_apr.cpp
+++ b/examples/Example_get_apr.cpp
@@ -82,7 +82,7 @@ int main(int argc, char **argv) {
 
         timer.verbose_flag = true;
 
-        MeshData<uint16_t> level;
+        PixelData<uint16_t> level;
 
         apr.interp_depth_ds(level);
 
@@ -125,14 +125,14 @@ int main(int argc, char **argv) {
             unsigned int blosc_comp_level = options.compress_level;
             unsigned int blosc_shuffle = 1;
 
-            MeshData<uint16_t> recon_image;
+            PixelData<uint16_t> recon_image;
 
             apr.interp_img(recon_image, apr.particles_intensities);
 
             TiffUtils::TiffInfo inputTiff(options.directory + options.input);
-            MeshData<uint16_t> inputImage = TiffUtils::getMesh<uint16_t>(inputTiff);
+            PixelData<uint16_t> inputImage = TiffUtils::getMesh<uint16_t>(inputTiff);
 
-            MeshData<int16_t> diff_image(inputImage.y_num,inputImage.x_num,inputImage.z_num,0);
+            PixelData<int16_t> diff_image(inputImage.y_num,inputImage.x_num,inputImage.z_num,0);
 
 #ifdef HAVE_OPENMP
 #pragma omp parallel for schedule(static)
diff --git a/examples/Example_get_apr.h b/examples/Example_get_apr.h
index d12e95cc..baa153c6 100644
--- a/examples/Example_get_apr.h
+++ b/examples/Example_get_apr.h
@@ -5,7 +5,7 @@
 #include <string>
 
 #include "algorithm/APRParameters.hpp"
-#include "data_structures/Mesh/MeshData.hpp"
+#include "data_structures/Mesh/PixelData.hpp"
 #include "algorithm/APRConverter.hpp"
 #include "data_structures/APR/APR.hpp"
 
diff --git a/examples/Example_ray_cast.cpp b/examples/Example_ray_cast.cpp
index ba06e2cb..e9081dde 100644
--- a/examples/Example_ray_cast.cpp
+++ b/examples/Example_ray_cast.cpp
@@ -68,7 +68,7 @@ int main(int argc, char **argv) {
 
     apr_raycaster.name = apr.name;
 
-    MeshData<uint16_t> views;
+    PixelData<uint16_t> views;
 
     /////////////
     ///
@@ -92,9 +92,9 @@ int main(int argc, char **argv) {
     if(options.original_image.size() > 0){
 
         TiffUtils::TiffInfo inputTiff(options.directory + options.original_image);
-        MeshData<uint16_t> original_image = TiffUtils::getMesh<uint16_t>(inputTiff);
+        PixelData<uint16_t> original_image = TiffUtils::getMesh<uint16_t>(inputTiff);
 
-        MeshData<uint16_t> mesh_views;
+        PixelData<uint16_t> mesh_views;
 
         apr_raycaster.perpsective_mesh_raycast(original_image,mesh_views);
 
diff --git a/examples/Example_reconstruct_image.cpp b/examples/Example_reconstruct_image.cpp
index 3d845b5e..3f83bc25 100644
--- a/examples/Example_reconstruct_image.cpp
+++ b/examples/Example_reconstruct_image.cpp
@@ -127,7 +127,7 @@ int main(int argc, char **argv) {
 
         if(options.output_pc_recon) {
             //create mesh data structure for reconstruction
-            MeshData<uint16_t> recon_pc;
+            PixelData<uint16_t> recon_pc;
 
             timer.start_timer("pc interp");
             //perform piece-wise constant interpolation
@@ -179,7 +179,7 @@ int main(int argc, char **argv) {
 
         // Intentionaly block-scoped since local type_recon will be destructed when block ends and release memory.
         {
-            MeshData<uint16_t> type_recon;
+            PixelData<uint16_t> type_recon;
 
             apr.interp_img(type_recon, type);
             TiffUtils::saveMeshAsTiff(options.directory + apr.name + "_type.tif", type_recon);
@@ -205,7 +205,7 @@ int main(int argc, char **argv) {
     if(options.output_smooth_recon) {
 
         //smooth reconstruction - requires float
-        MeshData<float> recon_smooth;
+        PixelData<float> recon_smooth;
         std::vector<float> scale_d = {2, 2, 2};
 
         timer.start_timer("smooth reconstrution");
diff --git a/examples/Example_reconstruct_patch.cpp b/examples/Example_reconstruct_patch.cpp
index 270bdcfe..7d292cc4 100644
--- a/examples/Example_reconstruct_patch.cpp
+++ b/examples/Example_reconstruct_patch.cpp
@@ -191,7 +191,7 @@ int main(int argc, char **argv) {
 
         if(options.output_pc_recon) {
             //create mesh data structure for reconstruction
-            MeshData<uint16_t> recon_pc;
+            PixelData<uint16_t> recon_pc;
 
 
 
@@ -245,7 +245,7 @@ int main(int argc, char **argv) {
 
         // Intentionaly block-scoped since local type_recon will be destructed when block ends and release memory.
         {
-            MeshData<uint16_t> type_recon;
+            PixelData<uint16_t> type_recon;
 
             aprReconstruction.interp_image_patch(apr,type_recon, type,reconPatch);
             TiffUtils::saveMeshAsTiff(options.directory + apr.name + "_type.tif", type_recon);
@@ -291,7 +291,7 @@ int main(int argc, char **argv) {
 
         // Intentionaly block-scoped since local type_recon will be destructed when block ends and release memory.
         {
-            MeshData<uint8_t> type_recon;
+            PixelData<uint8_t> type_recon;
 
             //level
             aprReconstruction.interp_image_patch(apr,type_recon, level,reconPatch);
diff --git a/libapr.i b/libapr.i
index d4688154..f0a6e374 100644
--- a/libapr.i
+++ b/libapr.i
@@ -68,7 +68,7 @@ namespace std {
 #include "src/numerics/APRNumerics.hpp"
 %}
 
-%include "src/data_structures/Mesh/MeshData.hpp"
+%include "src/data_structures/Mesh/PixelData.hpp"
 %include "src/data_structures/APR/APR.hpp"
 %include "src/data_structures/APR/APRIterator.hpp"
 %include "src/numerics/APRNumerics.hpp"
diff --git a/src/algorithm/APRConverter.hpp b/src/algorithm/APRConverter.hpp
index df8834e7..2805f02f 100644
--- a/src/algorithm/APRConverter.hpp
+++ b/src/algorithm/APRConverter.hpp
@@ -9,7 +9,7 @@
 #ifndef PARTPLAY_APR_CONVERTER_HPP
 #define PARTPLAY_APR_CONVERTER_HPP
 
-#include "../data_structures/Mesh/MeshData.hpp"
+#include "data_structures/Mesh/PixelData.hpp"
 #include "../io/TiffUtils.hpp"
 #include "../data_structures/APR/APR.hpp"
 
@@ -51,30 +51,30 @@ class APRConverter: public LocalIntensityScale, public ComputeGradient, public L
 private:
     //get apr without setting parameters, and with an already loaded image.
     template<typename T>
-    bool get_apr_method(APR<ImageType> &aAPR, MeshData<T> &input_image);
+    bool get_apr_method(APR<ImageType> &aAPR, PixelData<T> &input_image);
 
     //pointer to the APR structure so member functions can have access if they need
     const APR<ImageType> *apr;
 
     template<typename T>
-    void init_apr(APR<ImageType>& aAPR, MeshData<T>& input_image);
+    void init_apr(APR<ImageType>& aAPR, PixelData<T>& input_image);
 
     template<typename T>
-    void auto_parameters(const MeshData<T> &input_img);
+    void auto_parameters(const PixelData<T> &input_img);
 
     template<typename T>
     bool get_apr_method_from_file(APR<ImageType> &aAPR, const TiffUtils::TiffInfo &aTiffFile);
 
-    void get_gradient(MeshData<ImageType> &image_temp, MeshData<ImageType> &grad_temp, MeshData<float> &local_scale_temp, MeshData<float> &local_scale_temp2, float bspline_offset);
-    void get_local_intensity_scale(MeshData<float> &local_scale_temp, MeshData<float> &local_scale_temp2);
-    void get_local_particle_cell_set(MeshData<ImageType> &grad_temp, MeshData<float> &local_scale_temp, MeshData<float> &local_scale_temp2);
+    void get_gradient(PixelData<ImageType> &image_temp, PixelData<ImageType> &grad_temp, PixelData<float> &local_scale_temp, PixelData<float> &local_scale_temp2, float bspline_offset);
+    void get_local_intensity_scale(PixelData<float> &local_scale_temp, PixelData<float> &local_scale_temp2);
+    void get_local_particle_cell_set(PixelData<ImageType> &grad_temp, PixelData<float> &local_scale_temp, PixelData<float> &local_scale_temp2);
 };
 
 template <typename T>
 struct MinMax{T min; T max; };
 
 template <typename T>
-MinMax<T> getMinMax(const MeshData<T>& input_image) {
+MinMax<T> getMinMax(const PixelData<T>& input_image) {
     T minVal = std::numeric_limits<T>::max();
     T maxVal = std::numeric_limits<T>::min();
 
@@ -96,7 +96,7 @@ MinMax<T> getMinMax(const MeshData<T>& input_image) {
 template<typename ImageType> template<typename T>
 bool APRConverter<ImageType>::get_apr_method_from_file(APR<ImageType> &aAPR, const TiffUtils::TiffInfo &aTiffFile) {
     allocation_timer.start_timer("read tif input image");
-    MeshData<T> inputImage = TiffUtils::getMesh<T>(aTiffFile);
+    PixelData<T> inputImage = TiffUtils::getMesh<T>(aTiffFile);
     allocation_timer.stop_timer();
 
     method_timer.start_timer("calculate automatic parameters");
@@ -141,7 +141,7 @@ bool APRConverter<ImageType>::get_apr_method_from_file(APR<ImageType> &aAPR, con
  * Main method for constructing the APR from an input image
  */
 template<typename ImageType> template<typename T>
-bool APRConverter<ImageType>::get_apr_method(APR<ImageType> &aAPR, MeshData<T>& input_image) {
+bool APRConverter<ImageType>::get_apr_method(APR<ImageType> &aAPR, PixelData<T>& input_image) {
     apr = &aAPR; // in case it was called directly
 
     total_timer.start_timer("Total_pipeline_excluding_IO");
@@ -155,12 +155,12 @@ bool APRConverter<ImageType>::get_apr_method(APR<ImageType> &aAPR, MeshData<T>&
     //assuming uint16, the total memory cost shoudl be approximately (1 + 1 + 1/8 + 2/8 + 2/8) = 2 5/8 original image size in u16bit
     //storage of the particle cell tree for computing the pulling scheme
     allocation_timer.start_timer("init and copy image");
-    MeshData<ImageType> image_temp(input_image, false /* don't copy */); // global image variable useful for passing between methods, or re-using memory (should be the only full sized copy of the image)
-    MeshData<ImageType> grad_temp; // should be a down-sampled image
+    PixelData<ImageType> image_temp(input_image, false /* don't copy */); // global image variable useful for passing between methods, or re-using memory (should be the only full sized copy of the image)
+    PixelData<ImageType> grad_temp; // should be a down-sampled image
     grad_temp.initDownsampled(input_image.y_num, input_image.x_num, input_image.z_num, 0);
-    MeshData<float> local_scale_temp; // Used as down-sampled images for some averaging steps where it is useful to not lose precision, or get over-flow errors
+    PixelData<float> local_scale_temp; // Used as down-sampled images for some averaging steps where it is useful to not lose precision, or get over-flow errors
     local_scale_temp.initDownsampled(input_image.y_num, input_image.x_num, input_image.z_num);
-    MeshData<float> local_scale_temp2;
+    PixelData<float> local_scale_temp2;
     local_scale_temp2.initDownsampled(input_image.y_num, input_image.x_num, input_image.z_num);
     allocation_timer.stop_timer();
 
@@ -206,7 +206,7 @@ bool APRConverter<ImageType>::get_apr_method(APR<ImageType> &aAPR, MeshData<T>&
     method_timer.stop_timer();
 
     method_timer.start_timer("downsample_pyramid");
-    std::vector<MeshData<T>> downsampled_img;
+    std::vector<PixelData<T>> downsampled_img;
     //Down-sample the image for particle intensity estimation
     downsamplePyrmaid(input_image, downsampled_img, aAPR.level_max(), aAPR.level_min());
     method_timer.stop_timer();
@@ -229,7 +229,7 @@ bool APRConverter<ImageType>::get_apr_method(APR<ImageType> &aAPR, MeshData<T>&
 }
 
 template<typename ImageType>
-void APRConverter<ImageType>::get_local_particle_cell_set(MeshData<ImageType> &grad_temp, MeshData<float> &local_scale_temp, MeshData<float> &local_scale_temp2) {
+void APRConverter<ImageType>::get_local_particle_cell_set(PixelData<ImageType> &grad_temp, PixelData<float> &local_scale_temp, PixelData<float> &local_scale_temp2) {
     //
     //  Computes the Local Particle Cell Set from a down-sampled local intensity scale (\sigma) and gradient magnitude
     //
@@ -274,7 +274,7 @@ void APRConverter<ImageType>::get_local_particle_cell_set(MeshData<ImageType> &g
 }
 
 template<typename ImageType>
-void APRConverter<ImageType>::get_gradient(MeshData<ImageType> &image_temp, MeshData<ImageType> &grad_temp, MeshData<float> &local_scale_temp, MeshData<float> &local_scale_temp2, float bspline_offset) {
+void APRConverter<ImageType>::get_gradient(PixelData<ImageType> &image_temp, PixelData<ImageType> &grad_temp, PixelData<float> &local_scale_temp, PixelData<float> &local_scale_temp2, float bspline_offset) {
     //  Bevan Cheeseman 2018
     //  Calculate the gradient from the input image. (You could replace this method with your own)
     //  Input: full sized image.
@@ -335,7 +335,7 @@ void APRConverter<ImageType>::get_gradient(MeshData<ImageType> &image_temp, Mesh
 }
 
 template<typename ImageType>
-void APRConverter<ImageType>::get_local_intensity_scale(MeshData<float> &local_scale_temp, MeshData<float> &local_scale_temp2) {
+void APRConverter<ImageType>::get_local_intensity_scale(PixelData<float> &local_scale_temp, PixelData<float> &local_scale_temp2) {
     //
     //  Calculate the Local Intensity Scale (You could replace this method with your own)
     //
@@ -385,7 +385,7 @@ void APRConverter<ImageType>::get_local_intensity_scale(MeshData<float> &local_s
 
 
 template<typename ImageType> template<typename T>
-void APRConverter<ImageType>::init_apr(APR<ImageType>& aAPR,MeshData<T>& input_image){
+void APRConverter<ImageType>::init_apr(APR<ImageType>& aAPR,PixelData<T>& input_image){
     //
     //  Initializing the size of the APR, min and maximum level (in the data structures it is called depth)
     //
@@ -406,7 +406,7 @@ void APRConverter<ImageType>::init_apr(APR<ImageType>& aAPR,MeshData<T>& input_i
 }
 
 template<typename ImageType> template<typename T>
-void APRConverter<ImageType>::auto_parameters(const MeshData<T>& input_img){
+void APRConverter<ImageType>::auto_parameters(const PixelData<T>& input_img){
     //
     //  Simple automatic parameter selection for 3D APR Flouresence Images
     //
@@ -484,7 +484,7 @@ void APRConverter<ImageType>::auto_parameters(const MeshData<T>& input_img){
     float proportion_flat = freq[0]/(counter*1.0f);
     float proportion_next = freq[1]/(counter*1.0f);
 
-    MeshData<T> histogram;
+    PixelData<T> histogram;
     histogram.init(num_bins, 1, 1);
     std::copy(freq.begin(),freq.end(),histogram.mesh.begin());
 
diff --git a/src/algorithm/ComputeGradient.hpp b/src/algorithm/ComputeGradient.hpp
index 3c79c2d1..f14587d4 100644
--- a/src/algorithm/ComputeGradient.hpp
+++ b/src/algorithm/ComputeGradient.hpp
@@ -5,7 +5,7 @@
 #ifndef PARTPLAY_GRADIENT_HPP
 #define PARTPLAY_GRADIENT_HPP
 
-#include "../data_structures/Mesh/MeshData.hpp"
+#include "data_structures/Mesh/PixelData.hpp"
 #include "../io/TiffUtils.hpp"
 
 #ifdef HAVE_OPENMP
@@ -19,18 +19,18 @@ class ComputeGradient {
 public:
 
     template<typename T>
-    void get_smooth_bspline_3D(MeshData<T> &input, float lambda);
+    void get_smooth_bspline_3D(PixelData<T> &input, float lambda);
 
 // Calculate inverse B-Spline Transform
 
     template<typename T>
-    void calc_inv_bspline_y(MeshData<T> &input);
+    void calc_inv_bspline_y(PixelData<T> &input);
 
     template<typename T>
-    void calc_inv_bspline_x(MeshData<T> &input);
+    void calc_inv_bspline_x(PixelData<T> &input);
 
     template<typename T>
-    void calc_inv_bspline_z(MeshData<T> &input);
+    void calc_inv_bspline_z(PixelData<T> &input);
 
     struct three_temps {
         float temp_1, temp_2, temp_3;
@@ -40,22 +40,22 @@ class ComputeGradient {
 
     template<typename S>
     void
-    calc_bspline_fd_ds_mag(const MeshData<S> &input, MeshData<S> &grad, const float hx, const float hy, const float hz);
+    calc_bspline_fd_ds_mag(const PixelData<S> &input, PixelData<S> &grad, const float hx, const float hy, const float hz);
 
     template<typename T,typename S>
-    void mask_gradient(MeshData<T>& grad_ds,MeshData<S>& temp_ds,MeshData<T>& temp_full,APRParameters& par);
+    void mask_gradient(PixelData<T>& grad_ds,PixelData<S>& temp_ds,PixelData<T>& temp_full,APRParameters& par);
 
     template<typename T,typename S>
-    void threshold_gradient(MeshData<T> &grad, const MeshData<S> &img, const float Ip_th);
+    void threshold_gradient(PixelData<T> &grad, const PixelData<S> &img, const float Ip_th);
 
     template<typename T>
-    void bspline_filt_rec_y(MeshData<T> &image, float lambda, float tol);
+    void bspline_filt_rec_y(PixelData<T> &image, float lambda, float tol);
 
     template<typename T>
-    void bspline_filt_rec_x(MeshData<T> &image, float lambda, float tol);
+    void bspline_filt_rec_x(PixelData<T> &image, float lambda, float tol);
 
     template<typename T>
-    void bspline_filt_rec_z(MeshData<T> &image, float lambda, float tol);
+    void bspline_filt_rec_z(PixelData<T> &image, float lambda, float tol);
 
     inline float impulse_resp(float k, float rho, float omg);
 
@@ -65,7 +65,7 @@ class ComputeGradient {
 
 
 template<typename T,typename S>
-void ComputeGradient::mask_gradient(MeshData<T>& grad_ds,MeshData<S>& temp_ds,MeshData<T>& temp_full,APRParameters& par){
+void ComputeGradient::mask_gradient(PixelData<T>& grad_ds,PixelData<S>& temp_ds,PixelData<T>& temp_full,APRParameters& par){
     //
     //  Bevan Cheeseman 2018
     //
@@ -95,7 +95,7 @@ void ComputeGradient::mask_gradient(MeshData<T>& grad_ds,MeshData<S>& temp_ds,Me
 }
 
 template<typename T,typename S>
-void ComputeGradient::threshold_gradient(MeshData<T> &grad, const MeshData<S> &img, const float Ip_th){
+void ComputeGradient::threshold_gradient(PixelData<T> &grad, const PixelData<S> &img, const float Ip_th){
     //
     //  Bevan Cheeseman 2016
     //
@@ -109,7 +109,7 @@ void ComputeGradient::threshold_gradient(MeshData<T> &grad, const MeshData<S> &i
 }
 
 template<typename T>
-void ComputeGradient::get_smooth_bspline_3D(MeshData<T>& input, float lambda) {
+void ComputeGradient::get_smooth_bspline_3D(PixelData<T>& input, float lambda) {
     //
     //  Gets smoothing bspline co-efficients for 3D
     //
@@ -151,7 +151,7 @@ inline float ComputeGradient::impulse_resp_back(float k,float rho,float omg,floa
 }
 
 template<typename T>
-void ComputeGradient::bspline_filt_rec_y(MeshData<T>& image,float lambda,float tol){
+void ComputeGradient::bspline_filt_rec_y(PixelData<T>& image,float lambda,float tol){
     //
     //  Bevan Cheeseman 2016
     //
@@ -284,7 +284,7 @@ void ComputeGradient::bspline_filt_rec_y(MeshData<T>& image,float lambda,float t
 }
 
 template<typename T>
-void ComputeGradient::bspline_filt_rec_z(MeshData<T>& image,float lambda,float tol){
+void ComputeGradient::bspline_filt_rec_z(PixelData<T>& image,float lambda,float tol){
     //
     //  Bevan Cheeseman 2016
     //
@@ -439,7 +439,7 @@ void ComputeGradient::bspline_filt_rec_z(MeshData<T>& image,float lambda,float t
 }
 
 template<typename T>
-void ComputeGradient::bspline_filt_rec_x(MeshData<T>& image,float lambda,float tol){
+void ComputeGradient::bspline_filt_rec_x(PixelData<T>& image,float lambda,float tol){
     //
     //  Bevan Cheeseman 2016
     //
@@ -594,7 +594,7 @@ void ComputeGradient::bspline_filt_rec_x(MeshData<T>& image,float lambda,float t
  * Caclulation of signal value from B-Spline co-efficients
  */
 template<typename T>
-void ComputeGradient::calc_inv_bspline_y(MeshData<T>& input){
+void ComputeGradient::calc_inv_bspline_y(PixelData<T>& input){
     //  Bevan Cheeseman 2016
     //
     //  Inverse cubic bspline inverse filter in y direciton (Memory direction)
@@ -640,7 +640,7 @@ void ComputeGradient::calc_inv_bspline_y(MeshData<T>& input){
 }
 
 template<typename T>
-void ComputeGradient::calc_inv_bspline_z(MeshData<T>& input){
+void ComputeGradient::calc_inv_bspline_z(PixelData<T>& input){
     //  Bevan Cheeseman 2016
     //
     //  Inverse cubic bspline inverse filter in x direciton (Off-stride direction)
@@ -707,7 +707,7 @@ void ComputeGradient::calc_inv_bspline_z(MeshData<T>& input){
 
 
 template<typename T>
-void ComputeGradient::calc_inv_bspline_x(MeshData<T>& input) {
+void ComputeGradient::calc_inv_bspline_x(PixelData<T>& input) {
     //  Bevan Cheeseman 2016
     //
     //  Inverse cubic bspline inverse filter in x direciton (Off-stride direction)
@@ -779,7 +779,7 @@ void ComputeGradient::calc_inv_bspline_x(MeshData<T>& input) {
  * @param hz - step in z dir
  */
 template<typename S>
-void ComputeGradient::calc_bspline_fd_ds_mag(const MeshData<S> &input, MeshData<S> &grad, const float hx, const float hy,const float hz) {
+void ComputeGradient::calc_bspline_fd_ds_mag(const PixelData<S> &input, PixelData<S> &grad, const float hx, const float hy,const float hz) {
     const size_t z_num = input.z_num;
     const size_t x_num = input.x_num;
     const size_t y_num = input.y_num;
diff --git a/src/algorithm/LocalIntensityScale.hpp b/src/algorithm/LocalIntensityScale.hpp
index e6ae4e6d..fccf34ca 100644
--- a/src/algorithm/LocalIntensityScale.hpp
+++ b/src/algorithm/LocalIntensityScale.hpp
@@ -11,24 +11,24 @@ class LocalIntensityScale {
 
 public:
     template<typename T>
-    void calc_abs_diff(const MeshData<T> &input_image, MeshData<T> &var);
+    void calc_abs_diff(const PixelData<T> &input_image, PixelData<T> &var);
 
     template<typename T>
-    void calc_sat_mean_z(MeshData<T> &input, const size_t offset);
+    void calc_sat_mean_z(PixelData<T> &input, const size_t offset);
 
     template<typename T>
-    void calc_sat_mean_x(MeshData<T> &input, const size_t offset);
+    void calc_sat_mean_x(PixelData<T> &input, const size_t offset);
 
     template<typename T>
-    void calc_sat_mean_y(MeshData<T> &input, const size_t offset);
+    void calc_sat_mean_y(PixelData<T> &input, const size_t offset);
 
     void get_window(float &var_rescale, std::vector<int> &var_win, const APRParameters &par);
     template<typename T>
-    void rescale_var_and_threshold(MeshData<T>& var,const float var_rescale, const APRParameters& par);
+    void rescale_var_and_threshold(PixelData<T>& var,const float var_rescale, const APRParameters& par);
 };
 
 template<typename T>
-void LocalIntensityScale::rescale_var_and_threshold(MeshData<T> &var, const float var_rescale, const APRParameters &par) {
+void LocalIntensityScale::rescale_var_and_threshold(PixelData<T> &var, const float var_rescale, const APRParameters &par) {
     const float max_th = 60000.0;
 
     #ifdef HAVE_OPENMP
@@ -44,7 +44,7 @@ void LocalIntensityScale::rescale_var_and_threshold(MeshData<T> &var, const floa
 }
 
 template<typename T>
-void LocalIntensityScale::calc_abs_diff(const MeshData<T> &input_image, MeshData<T> &var) {
+void LocalIntensityScale::calc_abs_diff(const PixelData<T> &input_image, PixelData<T> &var) {
     #ifdef HAVE_OPENMP
 	#pragma omp parallel for default(shared)
     #endif
@@ -103,7 +103,7 @@ void LocalIntensityScale::get_window(float& var_rescale, std::vector<int>& var_w
  * @param offset
  */
 template<typename T>
-void LocalIntensityScale::calc_sat_mean_y(MeshData<T>& input, const size_t offset){
+void LocalIntensityScale::calc_sat_mean_y(PixelData<T>& input, const size_t offset){
     const size_t z_num = input.z_num;
     const size_t x_num = input.x_num;
     const size_t y_num = input.y_num;
@@ -162,7 +162,7 @@ void LocalIntensityScale::calc_sat_mean_y(MeshData<T>& input, const size_t offse
 }
 
 template<typename T>
-void LocalIntensityScale::calc_sat_mean_x(MeshData<T>& input, const size_t offset) {
+void LocalIntensityScale::calc_sat_mean_x(PixelData<T>& input, const size_t offset) {
     const size_t z_num = input.z_num;
     const size_t x_num = input.x_num;
     const size_t y_num = input.y_num;
@@ -223,7 +223,7 @@ void LocalIntensityScale::calc_sat_mean_x(MeshData<T>& input, const size_t offse
 }
 
 template<typename T>
-void LocalIntensityScale::calc_sat_mean_z(MeshData<T>& input,const size_t offset) {
+void LocalIntensityScale::calc_sat_mean_z(PixelData<T>& input,const size_t offset) {
     const size_t z_num = input.z_num;
     const size_t x_num = input.x_num;
     const size_t y_num = input.y_num;
diff --git a/src/algorithm/LocalParticleCellSet.hpp b/src/algorithm/LocalParticleCellSet.hpp
index c87b3e97..39dc9083 100644
--- a/src/algorithm/LocalParticleCellSet.hpp
+++ b/src/algorithm/LocalParticleCellSet.hpp
@@ -14,7 +14,7 @@ class LocalParticleCellSet {
     }
 
     template< typename T>
-    void compute_level_for_array(MeshData<T>& input, float k_factor, float rel_error) {
+    void compute_level_for_array(PixelData<T>& input, float k_factor, float rel_error) {
         //
         //  Takes the sqrt of the grad vector to caluclate the magnitude
         //
diff --git a/src/algorithm/PullingScheme.hpp b/src/algorithm/PullingScheme.hpp
index 679de785..c364afd4 100644
--- a/src/algorithm/PullingScheme.hpp
+++ b/src/algorithm/PullingScheme.hpp
@@ -13,7 +13,7 @@
 
 #include <cassert>
 #include "../data_structures/APR/APRIterator.hpp"
-#include "../data_structures/Mesh/MeshData.hpp"
+#include "data_structures/Mesh/PixelData.hpp"
 #include "../data_structures/APR/APR.hpp"
 
 #ifdef HAVE_OPENMP
@@ -57,12 +57,12 @@ class PullingScheme {
 
 public:
 
-    std::vector<MeshData<uint8_t>> particle_cell_tree;
+    std::vector<PixelData<uint8_t>> particle_cell_tree;
     unsigned int l_min;
     unsigned int l_max;
 
     template<typename T>
-    void fill(float k, const MeshData<T> &input);
+    void fill(float k, const PixelData<T> &input);
     void pulling_scheme_main();
     template<typename T>
     void initialize_particle_cell_tree(APR<T>& apr);
@@ -117,7 +117,7 @@ void PullingScheme::pulling_scheme_main() {
 }
 
 template<typename T>
-void PullingScheme::fill(const float k, const MeshData<T> &input) {
+void PullingScheme::fill(const float k, const PixelData<T> &input) {
     //  Bevan Cheeseman 2016
     //
     //  Updates the hash table from the down sampled images
diff --git a/src/data_structures/APR/APR.hpp b/src/data_structures/APR/APR.hpp
index 51697af3..3890354c 100644
--- a/src/data_structures/APR/APR.hpp
+++ b/src/data_structures/APR/APR.hpp
@@ -94,7 +94,7 @@ class APR {
     //////////////////////////
 
     template<typename U,typename V>
-    void interp_img(MeshData<U>& img,ExtraParticleData<V>& parts){
+    void interp_img(PixelData<U>& img,ExtraParticleData<V>& parts){
         //
         //  Bevan Cheeseman 2016
         //
@@ -105,7 +105,7 @@ class APR {
     }
 
     template<typename U>
-    void interp_depth_ds(MeshData<U>& img){
+    void interp_depth_ds(PixelData<U>& img){
         //
         //  Returns an image of the depth, this is down-sampled by one, as the Particle Cell solution reflects this
         //
@@ -114,7 +114,7 @@ class APR {
     }
 
     template<typename U>
-    void interp_depth(MeshData<U>& img){
+    void interp_depth(PixelData<U>& img){
         //
         //  Returns an image of the depth, this is down-sampled by one, as the Particle Cell solution reflects this
         //
@@ -123,7 +123,7 @@ class APR {
     }
 
     template<typename U>
-    void interp_type(MeshData<U>& img){
+    void interp_type(PixelData<U>& img){
         //
         //  Interpolates the APR
         //
@@ -132,7 +132,7 @@ class APR {
     }
 
     template<typename U,typename V>
-    void interp_parts_smooth(MeshData<U>& out_image,ExtraParticleData<V>& interp_data,std::vector<float> scale_d = {2,2,2}){
+    void interp_parts_smooth(PixelData<U>& out_image,ExtraParticleData<V>& interp_data,std::vector<float> scale_d = {2,2,2}){
         //
         //  Performs a smooth interpolation, based on the depth (level l) in each direction.
         //
@@ -141,7 +141,7 @@ class APR {
     }
 
     template<typename U,typename V>
-    void get_parts_from_img(std::vector<MeshData<U>>& img_by_level,ExtraParticleData<V>& parts){
+    void get_parts_from_img(std::vector<PixelData<U>>& img_by_level,ExtraParticleData<V>& parts){
         //
         //  Bevan Cheeseman 2016
         //
diff --git a/src/data_structures/APR/APRAccess.hpp b/src/data_structures/APR/APRAccess.hpp
index dff7a660..efc05783 100644
--- a/src/data_structures/APR/APRAccess.hpp
+++ b/src/data_structures/APR/APRAccess.hpp
@@ -9,7 +9,7 @@
 
 #include <map>
 #include <utility>
-#include "../../data_structures/Mesh/MeshData.hpp"
+#include "data_structures/Mesh/PixelData.hpp"
 
 //TODO: IT SHOULD NOT BE DEFINDED HERE SINCE IT DUPLICATES FROM PullingScheme
 #define SEED_TYPE 1
@@ -325,7 +325,7 @@ class APRAccess {
     }
 
     template<typename T>
-    void initialize_structure_from_particle_cell_tree(APR<T>& apr,std::vector<MeshData<uint8_t>>& layers){
+    void initialize_structure_from_particle_cell_tree(APR<T>& apr,std::vector<PixelData<uint8_t>>& layers){
        x_num.resize(level_max+1);
        y_num.resize(level_max+1);
        z_num.resize(level_max+1);
diff --git a/src/data_structures/APR/APRTree.hpp b/src/data_structures/APR/APRTree.hpp
index 1ad08e7b..7198f4e1 100644
--- a/src/data_structures/APR/APRTree.hpp
+++ b/src/data_structures/APR/APRTree.hpp
@@ -41,7 +41,7 @@ class APRTree {
         /// Set the iterator using random access by particle cell spatial index (x,y,z) and particle cell level
         ///
         ////////////////////////
-        std::vector<MeshData<uint8_t>> particle_cell_parent_tree;
+        std::vector<PixelData<uint8_t>> particle_cell_parent_tree;
         uint64_t l_max = apr.level_max()-1;
         uint64_t l_min = apr.level_min()-1; // extend one extra level
 
diff --git a/src/data_structures/Mesh/MeshData.hpp b/src/data_structures/Mesh/PixelData.hpp
similarity index 90%
rename from src/data_structures/Mesh/MeshData.hpp
rename to src/data_structures/Mesh/PixelData.hpp
index 4d4e3d51..b414812d 100644
--- a/src/data_structures/Mesh/MeshData.hpp
+++ b/src/data_structures/Mesh/PixelData.hpp
@@ -3,15 +3,13 @@
 //
 //  ImageGen 2016 Bevan Cheeseman
 //                Krzysztof Gonciarz
-//  Meshdata class for storing the image/mesh data
-//
-//
+//  Meshdata class for storing the image/pixel/mesh data
 //
 //
 ///////////////////////////////////////////////////////////////
 
-#ifndef PARTPLAY_MESHCLASS_H
-#define PARTPLAY_MESHCLASS_H
+#ifndef PIXEL_DATA_HPP
+#define PIXEL_DATA_HPP
 
 #include <vector>
 #include <cmath>
@@ -73,7 +71,7 @@ class ArrayWrapper
  * @tparam T type of mesh elements
  */
 template <typename T>
-class MeshData {
+class PixelData {
 public :
     // size of mesh and container for data
     size_t y_num;
@@ -85,7 +83,7 @@ public :
     /**
      * Constructor - initialize mesh with size of 0,0,0
      */
-    MeshData() { init(0, 0, 0); }
+    PixelData() { init(0, 0, 0); }
 
     /**
      * Constructor - initialize initial size of mesh to provided values
@@ -93,7 +91,7 @@ public :
      * @param aSizeOfX
      * @param aSizeOfZ
      */
-    MeshData(int aSizeOfY, int aSizeOfX, int aSizeOfZ) { init(aSizeOfY, aSizeOfX, aSizeOfZ); }
+    PixelData(int aSizeOfY, int aSizeOfX, int aSizeOfZ) { init(aSizeOfY, aSizeOfX, aSizeOfZ); }
 
     /**
      * Constructor - creates mesh with provided dimentions initialized to aInitVal
@@ -102,13 +100,13 @@ public :
      * @param aSizeOfZ
      * @param aInitVal - initial value of all elements
      */
-    MeshData(int aSizeOfY, int aSizeOfX, int aSizeOfZ, T aInitVal) { init(aSizeOfY, aSizeOfX, aSizeOfZ, aInitVal); }
+    PixelData(int aSizeOfY, int aSizeOfX, int aSizeOfZ, T aInitVal) { init(aSizeOfY, aSizeOfX, aSizeOfZ, aInitVal); }
 
     /**
      * Move constructor
      * @param aObj mesh to be moved
      */
-    MeshData(MeshData &&aObj) {
+    PixelData(PixelData &&aObj) {
         x_num = aObj.x_num;
         y_num = aObj.y_num;
         z_num = aObj.z_num;
@@ -120,7 +118,7 @@ public :
      * Move assignment operator
     * @param aObj
     */
-    MeshData& operator=(MeshData &&aObj) {
+    PixelData& operator=(PixelData &&aObj) {
         x_num = aObj.x_num;
         y_num = aObj.y_num;
         z_num = aObj.z_num;
@@ -134,7 +132,7 @@ public :
      * @param aMesh input mesh
      */
     template<typename U>
-    MeshData(const MeshData<U> &aMesh, bool aShouldCopyData) {
+    PixelData(const PixelData<U> &aMesh, bool aShouldCopyData) {
         init(aMesh.y_num, aMesh.x_num, aMesh.z_num);
         if (aShouldCopyData) std::copy(aMesh.mesh.begin(), aMesh.mesh.end(), mesh.begin());
     }
@@ -145,8 +143,8 @@ public :
      * @return created object by value
      */
     template <typename U>
-    MeshData<U> toType() const {
-        MeshData<U> new_value(y_num, x_num, z_num);
+    PixelData<U> toType() const {
+        PixelData<U> new_value(y_num, x_num, z_num);
         std::copy(mesh.begin(), mesh.end(), new_value.mesh.begin());
         return new_value;
     }
@@ -198,7 +196,7 @@ public :
      * @param aNumberOfBlocks in how many chunks copy will be done
      */
     template<typename U>
-    void copyFromMesh(const MeshData<U> &aInputMesh, unsigned int aNumberOfBlocks = 8) {
+    void copyFromMesh(const PixelData<U> &aInputMesh, unsigned int aNumberOfBlocks = 8) {
         aNumberOfBlocks = std::min((unsigned int)z_num, aNumberOfBlocks);
         unsigned int numOfElementsPerBlock = z_num/aNumberOfBlocks;
 
@@ -259,7 +257,7 @@ public :
      * @param aInputMesh
      */
     template<typename S>
-    void init(const MeshData<S> &aInputMesh) {
+    void init(const PixelData<S> &aInputMesh) {
         init(aInputMesh.y_num, aInputMesh.x_num, aInputMesh.z_num);
     }
 
@@ -313,7 +311,7 @@ public :
      * @param aMesh - mesh used to get dimensions
      */
     template <typename U>
-    void initDownsampled(const MeshData<U> &aMesh) {
+    void initDownsampled(const PixelData<U> &aMesh) {
         const int z_num_ds = ceil(1.0*aMesh.z_num/2.0);
         const int x_num_ds = ceil(1.0*aMesh.x_num/2.0);
         const int y_num_ds = ceil(1.0*aMesh.y_num/2.0);
@@ -327,7 +325,7 @@ public :
      * @param aInitVal
      */
     template <typename U>
-    void initDownsampled(const MeshData<U> &aMesh, T aInitVal) {
+    void initDownsampled(const PixelData<U> &aMesh, T aInitVal) {
         const int z_num_ds = ceil(1.0*aMesh.z_num/2.0);
         const int x_num_ds = ceil(1.0*aMesh.x_num/2.0);
         const int y_num_ds = ceil(1.0*aMesh.y_num/2.0);
@@ -339,7 +337,7 @@ public :
      * Swaps data of meshes this <-> aObj
      * @param aObj
      */
-    void swap(MeshData &aObj) {
+    void swap(PixelData &aObj) {
         std::swap(x_num, aObj.x_num);
         std::swap(y_num, aObj.y_num);
         std::swap(z_num, aObj.z_num);
@@ -356,7 +354,7 @@ public :
      * @param aNumberOfBlocks - in how many chunks copy will be done
      */
     template<typename U, typename R>
-    void copyFromMeshWithUnaryOp(const MeshData<U> &aInputMesh, R aOp, size_t aNumberOfBlocks = 10) {
+    void copyFromMeshWithUnaryOp(const PixelData<U> &aInputMesh, R aOp, size_t aNumberOfBlocks = 10) {
         aNumberOfBlocks = std::min(aInputMesh.z_num, (size_t)aNumberOfBlocks);
         size_t numOfElementsPerBlock = aInputMesh.z_num/aNumberOfBlocks;
 
@@ -427,21 +425,21 @@ public :
         }
     }
 
-    friend std::ostream & operator<<(std::ostream &os, const MeshData<T> &obj) {
-        os << "MeshData: size(Y/X/Z)=" << obj.y_num << "/" << obj.x_num << "/" << obj.z_num << " vSize:" << obj.mesh.size() << " vCapacity:" << obj.mesh.capacity() << " elementSize:" << sizeof(T);
+    friend std::ostream & operator<<(std::ostream &os, const PixelData<T> &obj) {
+        os << "PixelData: size(Y/X/Z)=" << obj.y_num << "/" << obj.x_num << "/" << obj.z_num << " vSize:" << obj.mesh.size() << " vCapacity:" << obj.mesh.capacity() << " elementSize:" << sizeof(T);
         return os;
     }
 
 private:
 
-    MeshData(const MeshData&) = delete; // make it noncopyable
-    MeshData& operator=(const MeshData&) = delete; // make it not assignable
+    PixelData(const PixelData&) = delete; // make it noncopyable
+    PixelData& operator=(const PixelData&) = delete; // make it not assignable
 
 };
 
 
 template<typename T, typename S, typename R, typename C>
-void downsample(const MeshData<T> &aInput, MeshData<S> &aOutput, R reduce, C constant_operator, bool aInitializeOutput = false) {
+void downsample(const PixelData<T> &aInput, PixelData<S> &aOutput, R reduce, C constant_operator, bool aInitializeOutput = false) {
     const size_t z_num = aInput.z_num;
     const size_t x_num = aInput.x_num;
     const size_t y_num = aInput.y_num;
@@ -495,7 +493,7 @@ void downsample(const MeshData<T> &aInput, MeshData<S> &aOutput, R reduce, C con
 }
 
 template<typename T>
-void downsamplePyrmaid(MeshData<T> &original_image, std::vector<MeshData<T>> &downsampled, size_t l_max, size_t l_min) {
+void downsamplePyrmaid(PixelData<T> &original_image, std::vector<PixelData<T>> &downsampled, size_t l_max, size_t l_min) {
     downsampled.resize(l_max + 1); // each level is kept at same index
     downsampled.back().swap(original_image); // put original image at l_max index
 
@@ -507,4 +505,4 @@ void downsamplePyrmaid(MeshData<T> &original_image, std::vector<MeshData<T>> &do
     }
 }
 
-#endif //PARTPLAY_MESHCLASS_H
+#endif //PIXEL_DATA_HPP
diff --git a/src/io/APRWriter.hpp b/src/io/APRWriter.hpp
index f8dfe53c..d8848ff4 100644
--- a/src/io/APRWriter.hpp
+++ b/src/io/APRWriter.hpp
@@ -368,7 +368,7 @@ class APRWriter {
     }
 
     template<typename ImageType>
-    float write_mesh_to_hdf5(MeshData<ImageType>& input_mesh,const std::string &save_loc, const std::string &file_name,unsigned int blosc_comp_type = BLOSC_ZSTD, unsigned int blosc_comp_level = 2, unsigned int blosc_shuffle=1){
+    float write_mesh_to_hdf5(PixelData<ImageType>& input_mesh,const std::string &save_loc, const std::string &file_name,unsigned int blosc_comp_type = BLOSC_ZSTD, unsigned int blosc_comp_level = 2, unsigned int blosc_shuffle=1){
         std::string hdf5_file_name = save_loc + file_name + "_pixels.h5";
 
         AprFile f{hdf5_file_name, AprFile::Operation::WRITE};
diff --git a/src/io/TiffUtils.hpp b/src/io/TiffUtils.hpp
index b4d18ea5..474afb17 100644
--- a/src/io/TiffUtils.hpp
+++ b/src/io/TiffUtils.hpp
@@ -10,7 +10,7 @@
 #include <string>
 #include <tiffio.h>
 #include <sstream>
-#include "../data_structures/Mesh/MeshData.hpp"
+#include "data_structures/Mesh/PixelData.hpp"
 
 
 namespace TiffUtils {
@@ -156,7 +156,7 @@ namespace TiffUtils {
      * @return mesh with tiff or empty mesh if reading file failed
      */
     template<typename T>
-    MeshData<T> getMesh(const std::string &aFileName) {
+    PixelData<T> getMesh(const std::string &aFileName) {
         TiffInfo tiffInfo(aFileName);
         return getMesh<T>(tiffInfo);
     }
@@ -168,9 +168,9 @@ namespace TiffUtils {
      * @return mesh with tiff or empty mesh if reading file failed
      */
     template<typename T>
-    MeshData<T> getMesh(const TiffInfo &aTiff) {
-        if (!aTiff.isFileOpened()) return MeshData<T>();
-        MeshData<T> mesh(aTiff.iImgHeight, aTiff.iImgWidth, aTiff.iNumberOfDirectories);
+    PixelData<T> getMesh(const TiffInfo &aTiff) {
+        if (!aTiff.isFileOpened()) return PixelData<T>();
+        PixelData<T> mesh(aTiff.iImgHeight, aTiff.iImgWidth, aTiff.iNumberOfDirectories);
         getMesh(aTiff, mesh);
         return mesh;
     }
@@ -183,7 +183,7 @@ namespace TiffUtils {
     * @return mesh with tiff or empty mesh if reading file failed
     */
     template<typename T>
-    void getMesh(const TiffInfo &aTiff, MeshData<T> &aInputMesh) {
+    void getMesh(const TiffInfo &aTiff, PixelData<T> &aInputMesh) {
         if (!aTiff.isFileOpened()) return;
 
         std::cout << "getMesh: " << aInputMesh << std::endl;
@@ -192,7 +192,7 @@ namespace TiffUtils {
         const long stripSize = TIFFStripSize(aTiff.iFile);
         std::cout << __func__ << ": ScanlineSize=" << TIFFScanlineSize(aTiff.iFile) << " StripSize=" << stripSize << " NumberOfStrips=" << TIFFNumberOfStrips(aTiff.iFile) << std::endl;
 
-        // Read TIF to MeshData
+        // Read TIF to PixelData
         size_t currentOffset = 0;
         for (uint32_t i = 0; i < aTiff.iNumberOfDirectories; ++i) {
             TIFFSetDirectory(aTiff.iFile, i);
@@ -217,7 +217,7 @@ namespace TiffUtils {
      * @param aData mesh with data
      */
     template<typename T>
-    void saveMeshAsTiff(const std::string &aFileName, const MeshData<T> &aData) {
+    void saveMeshAsTiff(const std::string &aFileName, const PixelData<T> &aData) {
         std::cout << __func__ << ": " << "FileName: [" << aFileName << "] " << aData << std::endl;
 
         // Set proper dimensions (x and y are exchanged giving transpose)
@@ -282,9 +282,9 @@ namespace TiffUtils {
      * @param aData mesh with data
      */
     template<typename T>
-    void saveMeshAsTiffUint16(const std::string &filename, const MeshData<T> &aData) {
+    void saveMeshAsTiffUint16(const std::string &filename, const PixelData<T> &aData) {
         //  Converts the data to uint16t then writes it (requires creation of a complete copy of the data)
-        MeshData<uint16_t> mesh16{aData, true /*copy data*/};
+        PixelData<uint16_t> mesh16{aData, true /*copy data*/};
         saveMeshAsTiff(filename, mesh16);
     }
 }
diff --git a/src/numerics/APRRaycaster.hpp b/src/numerics/APRRaycaster.hpp
index d1350148..6c986e50 100644
--- a/src/numerics/APRRaycaster.hpp
+++ b/src/numerics/APRRaycaster.hpp
@@ -42,10 +42,10 @@ class APRRaycaster {
     std::string name = "raycast";
 
     template<typename U,typename S,typename V,class BinaryOperation>
-    void perform_raycast(APR<U>& apr,ExtraParticleData<S>& particle_data,MeshData<V>& cast_views,BinaryOperation op);
+    void perform_raycast(APR<U>& apr,ExtraParticleData<S>& particle_data,PixelData<V>& cast_views,BinaryOperation op);
 
     template<typename S,typename U>
-    float perpsective_mesh_raycast(MeshData<S>& image,MeshData<U>& cast_views);
+    float perpsective_mesh_raycast(PixelData<S>& image,PixelData<U>& cast_views);
 
     // Stuff below is hiding implementation so eventually glm header files do not
     // have to be delivered with libAPR
@@ -58,7 +58,7 @@ class APRRaycaster {
 
 
 template<typename U,typename S,typename V,class BinaryOperation>
-void APRRaycaster::perform_raycast(APR<U>& apr,ExtraParticleData<S>& particle_data,MeshData<V>& cast_views,BinaryOperation op) {
+void APRRaycaster::perform_raycast(APR<U>& apr,ExtraParticleData<S>& particle_data,PixelData<V>& cast_views,BinaryOperation op) {
 
     //
     //  Bevan Cheeseman 2018
@@ -117,7 +117,7 @@ void APRRaycaster::perform_raycast(APR<U>& apr,ExtraParticleData<S>& particle_da
     ///
     /////////////////////////////////////
 
-    std::vector<MeshData<S>> depth_slice;
+    std::vector<PixelData<S>> depth_slice;
     depth_slice.resize(apr.level_max() + 1);
 
     std::vector<float> depth_vec;
@@ -294,7 +294,7 @@ void APRRaycaster::perform_raycast(APR<U>& apr,ExtraParticleData<S>& particle_da
 }
 
 template<typename S,typename U>
-float APRRaycaster::perpsective_mesh_raycast(MeshData<S>& image,MeshData<U>& cast_views) {
+float APRRaycaster::perpsective_mesh_raycast(PixelData<S>& image,PixelData<U>& cast_views) {
     //
     //  Bevan Cheeseman 2017
     //
@@ -343,7 +343,7 @@ float APRRaycaster::perpsective_mesh_raycast(MeshData<S>& image,MeshData<U>& cas
     for (float theta = theta_0; theta <= theta_f; theta += theta_delta) {
         initObjects(imageWidth, imageHeight, radius, theta, x0, y0, z0, x0f, y0f, z0f);
 
-        MeshData<S> proj_img;
+        PixelData<S> proj_img;
         proj_img.init(imageHeight, imageWidth, 1, 0);
 
         unsigned int z_, x_, j_;
diff --git a/src/numerics/APRReconstruction.hpp b/src/numerics/APRReconstruction.hpp
index 200663cc..669f65ca 100644
--- a/src/numerics/APRReconstruction.hpp
+++ b/src/numerics/APRReconstruction.hpp
@@ -22,7 +22,7 @@ class APRReconstruction {
 public:
 
     template<typename U,typename V,typename S>
-    void interp_img(APR<S>& apr, MeshData<U>& img,ExtraParticleData<V>& parts){
+    void interp_img(APR<S>& apr, PixelData<U>& img,ExtraParticleData<V>& parts){
         //
         //  Bevan Cheeseman 2016
         //
@@ -74,7 +74,7 @@ class APRReconstruction {
     }
 
     template<typename U,typename V,typename S>
-    void interp_image_patch(APR<S>& apr, MeshData<U>& img,ExtraParticleData<V>& parts,ReconPatch& reconPatch){
+    void interp_image_patch(APR<S>& apr, PixelData<U>& img,ExtraParticleData<V>& parts,ReconPatch& reconPatch){
         //
         //  Bevan Cheeseman 2016
         //
@@ -187,7 +187,7 @@ class APRReconstruction {
 
 
     template<typename U,typename S>
-    void interp_depth_ds(APR<S>& apr,MeshData<U>& img){
+    void interp_depth_ds(APR<S>& apr,PixelData<U>& img){
         //
         //  Returns an image of the depth, this is down-sampled by one, as the Particle Cell solution reflects this
         //
@@ -209,7 +209,7 @@ class APRReconstruction {
 
         }
 
-        MeshData<U> temp;
+        PixelData<U> temp;
 
         interp_img(apr,temp,depth_parts);
 
@@ -220,7 +220,7 @@ class APRReconstruction {
     }
 
     template<typename U,typename S>
-    void interp_level(APR<S> &apr, MeshData<U> &img){
+    void interp_level(APR<S> &apr, PixelData<U> &img){
         //
         //  Returns an image of the depth, this is down-sampled by one, as the Particle Cell solution reflects this
         //
@@ -250,7 +250,7 @@ class APRReconstruction {
     }
 
     template<typename U,typename S>
-    void interp_type(APR<S>& apr,MeshData<U>& img){
+    void interp_type(APR<S>& apr,PixelData<U>& img){
 
         //get depth
         ExtraParticleData<U> type_parts(apr);
@@ -279,7 +279,7 @@ class APRReconstruction {
     }
 
     template<typename T>
-    void calc_sat_adaptive_y(MeshData<T>& input,MeshData<uint8_t>& offset_img,float scale_in,unsigned int offset_max_in,const unsigned int d_max){
+    void calc_sat_adaptive_y(PixelData<T>& input,PixelData<uint8_t>& offset_img,float scale_in,unsigned int offset_max_in,const unsigned int d_max){
         //
         //  Bevan Cheeseman 2016
         //
@@ -403,7 +403,7 @@ class APRReconstruction {
 
 
     template<typename T>
-    void calc_sat_adaptive_x(MeshData<T>& input,MeshData<uint8_t>& offset_img,float scale_in,unsigned int offset_max_in,const unsigned int d_max){
+    void calc_sat_adaptive_x(PixelData<T>& input,PixelData<uint8_t>& offset_img,float scale_in,unsigned int offset_max_in,const unsigned int d_max){
         //
         //  Adaptive form of Matteusz' SAT code.
         //
@@ -521,7 +521,7 @@ class APRReconstruction {
 
 
     template<typename T>
-    void calc_sat_adaptive_z(MeshData<T>& input,MeshData<uint8_t>& offset_img,float scale_in,unsigned int offset_max_in,const unsigned int d_max ){
+    void calc_sat_adaptive_z(PixelData<T>& input,PixelData<uint8_t>& offset_img,float scale_in,unsigned int offset_max_in,const unsigned int d_max ){
 
         // The same, but in place
 
@@ -633,7 +633,7 @@ class APRReconstruction {
 
 
     template<typename U,typename V,typename S>
-    void interp_parts_smooth(APR<S>& apr,MeshData<U>& out_image,ExtraParticleData<V>& interp_data,std::vector<float> scale_d = {2,2,2}){
+    void interp_parts_smooth(APR<S>& apr,PixelData<U>& out_image,ExtraParticleData<V>& interp_data,std::vector<float> scale_d = {2,2,2}){
         //
         //  Performs a smooth interpolation, based on the depth (level l) in each direction.
         //
@@ -641,8 +641,8 @@ class APRReconstruction {
         APRTimer timer;
         timer.verbose_flag = false;
 
-        MeshData<U> pc_image;
-        MeshData<uint8_t> k_img;
+        PixelData<U> pc_image;
+        PixelData<uint8_t> k_img;
 
         unsigned int offset_max = 20;
 
diff --git a/src/numerics/MeshNumerics.hpp b/src/numerics/MeshNumerics.hpp
index c1f22a0e..bdd1919e 100644
--- a/src/numerics/MeshNumerics.hpp
+++ b/src/numerics/MeshNumerics.hpp
@@ -5,7 +5,7 @@
 #ifndef PARTPLAY_MESHNUMERICS_HPP
 #define PARTPLAY_MESHNUMERICS_HPP
 
-#include "../data_structures/Mesh/MeshData.hpp"
+#include "data_structures/Mesh/PixelData.hpp"
 #include "../misc/APRTimer.hpp"
 
 class MeshNumerics {
@@ -13,7 +13,7 @@ class MeshNumerics {
 public:
 
     template<typename U,typename V>
-    static float compute_gradient(const MeshData<U>& input_data,std::vector<MeshData<V>>& output_data,std::vector<float> delta = {1.0f,1.0f,1.0f}){
+    static float compute_gradient(const PixelData<U>& input_data,std::vector<PixelData<V>>& output_data,std::vector<float> delta = {1.0f,1.0f,1.0f}){
         //
         //  Computes gradient magnitude using finite differences
         //
diff --git a/test/APRTest.cpp b/test/APRTest.cpp
index b4a9a6dd..f60294af 100644
--- a/test/APRTest.cpp
+++ b/test/APRTest.cpp
@@ -4,7 +4,7 @@
 
 #include <gtest/gtest.h>
 #include "data_structures/APR/APR.hpp"
-#include "data_structures/Mesh/MeshData.hpp"
+#include "data_structures/Mesh/PixelData.hpp"
 #include "algorithm/APRConverter.hpp"
 #include <utility>
 #include <cmath>
@@ -12,13 +12,13 @@
 struct TestData{
 
     APR<uint16_t> apr;
-    MeshData<uint16_t> img_level;
-    MeshData<uint16_t> img_type;
-    MeshData<uint16_t> img_original;
-    MeshData<uint16_t> img_pc;
-    MeshData<uint16_t> img_x;
-    MeshData<uint16_t> img_y;
-    MeshData<uint16_t> img_z;
+    PixelData<uint16_t> img_level;
+    PixelData<uint16_t> img_type;
+    PixelData<uint16_t> img_original;
+    PixelData<uint16_t> img_pc;
+    PixelData<uint16_t> img_x;
+    PixelData<uint16_t> img_y;
+    PixelData<uint16_t> img_z;
 
     std::string filename;
     std::string output_name;
diff --git a/test/ComputeGradientTest.cpp b/test/ComputeGradientTest.cpp
index bcad417a..ebf2557a 100644
--- a/test/ComputeGradientTest.cpp
+++ b/test/ComputeGradientTest.cpp
@@ -3,7 +3,7 @@
  */
 
 #include <gtest/gtest.h>
-#include "data_structures/Mesh/MeshData.hpp"
+#include "data_structures/Mesh/PixelData.hpp"
 #include "algorithm/ComputeGradient.hpp"
 
 namespace {
@@ -14,7 +14,7 @@ namespace {
      * @return true if same
      */
     template <typename T>
-    bool compare(MeshData<T> &mesh, const float *data, const float epsilon) {
+    bool compare(PixelData<T> &mesh, const float *data, const float epsilon) {
         size_t dataIdx = 0;
         for (size_t z = 0; z < mesh.z_num; ++z) {
             for (size_t y = 0; y < mesh.y_num; ++y) {
@@ -33,7 +33,7 @@ namespace {
 
     TEST(ComputeGradientTest, 2D_XY) {
         {   // Corner points
-            MeshData<float> m(6, 6, 1, 0);
+            PixelData<float> m(6, 6, 1, 0);
             // expect gradient is 3x3 X/Y plane
             float expect[] = {1.41, 0, 4.24,
                               0, 0, 0,
@@ -43,33 +43,33 @@ namespace {
             m(5, 0, 0) = 4;
             m(0, 5, 0) = 6;
             m(5, 5, 0) = 8;
-            MeshData<float> grad;
+            PixelData<float> grad;
             grad.initDownsampled(m, 0);
             ComputeGradient cg;
             cg.calc_bspline_fd_ds_mag(m, grad, 1, 1, 1);
             ASSERT_TRUE(compare(grad, expect, 0.01));
         }
         {   // In the middle
-            MeshData<float> m(6, 6, 1, 0);
+            PixelData<float> m(6, 6, 1, 0);
             // expect gradient is 3x3 X/Y plane
             float expect[] = {1, 1, 0,
                               1, 0, 0,
                               0, 0, 0};
             // put values in corners
             m(1, 1, 0) = 2;
-            MeshData<float> grad;
+            PixelData<float> grad;
             grad.initDownsampled(m, 0);
             ComputeGradient cg;
             cg.calc_bspline_fd_ds_mag(m, grad, 1, 1, 1);
             ASSERT_TRUE(compare(grad, expect, 0.01));
         }
         {   // One pixel image 1x1x1
-            MeshData<float> m(1, 1, 1, 0);
+            PixelData<float> m(1, 1, 1, 0);
             // expect gradient is 3x3 X/Y plane
             float expect[] = {0};
             // put values in corners
             m(0, 0, 0) = 2;
-            MeshData<float> grad;
+            PixelData<float> grad;
             grad.initDownsampled(m, 0);
             ComputeGradient cg;
             cg.calc_bspline_fd_ds_mag(m, grad, 1, 1, 1);
@@ -79,7 +79,7 @@ namespace {
     }
 
     TEST(ComputeGradientTest, Corners3D) {
-        MeshData<float> m(6, 6, 4, 0);
+        PixelData<float> m(6, 6, 4, 0);
         // expect gradient is 3x3x2 X/Y/Z plane
         float expect[] = { 1.73, 0, 5.19,
                            0,    0,    0,
@@ -98,7 +98,7 @@ namespace {
         m(0, 5, 3) = 14;
         m(5, 5, 3) = 16;
 
-        MeshData<float> grad;
+        PixelData<float> grad;
         grad.initDownsampled(m, 0);
         ComputeGradient cg;
         cg.calc_bspline_fd_ds_mag(m, grad, 1, 1, 1);
diff --git a/test/MeshDataTest.cpp b/test/MeshDataTest.cpp
index 421d984e..861887c3 100644
--- a/test/MeshDataTest.cpp
+++ b/test/MeshDataTest.cpp
@@ -2,7 +2,7 @@
  * Created by Krzysztof Gonciarz 2018
  */
 #include <gtest/gtest.h>
-#include "data_structures/Mesh/MeshData.hpp"
+#include "data_structures/Mesh/PixelData.hpp"
 
 namespace {
     class MeshDataTest : public ::testing::Test {
@@ -25,7 +25,7 @@ namespace {
         const int xLen = 254;
         const int zLen = 123;
         const size_t sizeOfMesh = (size_t)yLen * xLen * zLen;
-        MeshData<MESH_TYPE> m{yLen, xLen, zLen};
+        PixelData<MESH_TYPE> m{yLen, xLen, zLen};
     };
 
     class MeshDataParameterTest : public ::testing::TestWithParam<int> {
@@ -48,13 +48,13 @@ namespace {
         const int xLen = 20;
         const int zLen = 30;
         const size_t sizeOfMesh = (size_t)yLen * xLen * zLen;
-        MeshData<MESH_TYPE> m{yLen, xLen, zLen};
+        PixelData<MESH_TYPE> m{yLen, xLen, zLen};
     };
 
     TEST(MeshDataSimpleTest, ConstructorTest) {
         // default
         {
-            MeshData<int> md;
+            PixelData<int> md;
             ASSERT_EQ(md.x_num, 0);
             ASSERT_EQ(md.y_num, 0);
             ASSERT_EQ(md.z_num, 0);
@@ -63,7 +63,7 @@ namespace {
 
         // size provided
         {
-            MeshData<int> md(100, 200, 300);
+            PixelData<int> md(100, 200, 300);
             ASSERT_EQ(md.x_num, 200);
             ASSERT_EQ(md.y_num, 100);
             ASSERT_EQ(md.z_num, 300);
@@ -73,11 +73,11 @@ namespace {
         // mesh provided
         {
             // generate some data
-            MeshData<int> md(3,4,5);
+            PixelData<int> md(3,4,5);
             for (size_t i = 0; i < md.mesh.size(); ++i) md.mesh[i] = i + 1;
 
             // test constructor
-            MeshData<char> testedMesh(md, true);
+            PixelData<char> testedMesh(md, true);
             ASSERT_EQ(md.mesh.size(), testedMesh.mesh.size());
             ASSERT_EQ(md.x_num, testedMesh.x_num);
             ASSERT_EQ(md.y_num, testedMesh.y_num);
@@ -91,7 +91,7 @@ namespace {
     TEST_F(MeshDataTest, ToTypeTest) {
         // Change type and compare if still OK
         typedef short NEW_TYPE;
-        MeshData<NEW_TYPE> mf = m.toType<NEW_TYPE>();
+        PixelData<NEW_TYPE> mf = m.toType<NEW_TYPE>();
         ASSERT_EQ(mf.x_num, xLen);
         ASSERT_EQ(mf.y_num, yLen);
         ASSERT_EQ(mf.z_num, zLen);
@@ -121,7 +121,7 @@ namespace {
     }
 
     TEST_P(MeshDataParameterTest, BlockCopyDataTest) {
-        MeshData<unsigned short> mNew(yLen, xLen, zLen);
+        PixelData<unsigned short> mNew(yLen, xLen, zLen);
 
         int numOfBlocks = GetParam();
         mNew.copyFromMesh(m, numOfBlocks);
@@ -142,7 +142,7 @@ namespace {
 
     TEST_F(MeshDataTest, InitializeTest) {
         {   // Size and initial value known
-            MeshData<int> md;
+            PixelData<int> md;
             md.init(3, 4, 5, 123);
             ASSERT_EQ(md.y_num, 3);
             ASSERT_EQ(md.x_num, 4);
@@ -152,7 +152,7 @@ namespace {
             for (int i = 0; i < size; ++i) ASSERT_EQ(md.mesh[i], 123);
         }
         {   // Use data from other mesh
-            MeshData<int> md;
+            PixelData<int> md;
             md.init(m);
             ASSERT_EQ(md.x_num, xLen);
             ASSERT_EQ(md.y_num, yLen);
@@ -160,7 +160,7 @@ namespace {
             ASSERT_EQ(md.mesh.size(), sizeOfMesh);
         }
         {   // Size and default value for type used
-            MeshData<int> md;
+            PixelData<int> md;
             md.init(3, 4, 5);
             ASSERT_EQ(md.y_num, 3);
             ASSERT_EQ(md.x_num, 4);
@@ -172,7 +172,7 @@ namespace {
 
     TEST_F(MeshDataTest, InitDownsampledTest) {
         {
-            MeshData<int> md;
+            PixelData<int> md;
             md.initDownsampled(3, 5, 7);
             ASSERT_EQ(md.y_num, 2);
             ASSERT_EQ(md.x_num, 3);
@@ -181,7 +181,7 @@ namespace {
             ASSERT_EQ(md.mesh.size(), size);
         }
         {
-            MeshData<int> md;
+            PixelData<int> md;
             md.initDownsampled(4, 6, 8);
             ASSERT_EQ(md.y_num, 2);
             ASSERT_EQ(md.x_num, 3);
@@ -190,8 +190,8 @@ namespace {
             ASSERT_EQ(md.mesh.size(), size);
         }
         {
-            MeshData<int> md;
-            md.initDownsampled(MeshData<char>(4, 6, 8));
+            PixelData<int> md;
+            md.initDownsampled(PixelData<char>(4, 6, 8));
             ASSERT_EQ(md.y_num, 2);
             ASSERT_EQ(md.x_num, 3);
             ASSERT_EQ(md.z_num, 4);
@@ -199,8 +199,8 @@ namespace {
             ASSERT_EQ(md.mesh.size(), size);
         }
         {
-            MeshData<int> md;
-            md.initDownsampled(MeshData<float>(3, 5, 7), 2);
+            PixelData<int> md;
+            md.initDownsampled(PixelData<float>(3, 5, 7), 2);
             ASSERT_EQ(md.y_num, 2);
             ASSERT_EQ(md.x_num, 3);
             ASSERT_EQ(md.z_num, 4);
@@ -211,10 +211,10 @@ namespace {
     }
 
     TEST(MeshDataSimpleTest, UnaryOpTest) {
-        MeshData<int> m(1, 5, 1, 1);
+        PixelData<int> m(1, 5, 1, 1);
         for (size_t i = 0; i < m.mesh.size(); ++i) m.mesh[i] = i + 1;
 
-        MeshData<int> m2(1, 5, 1);
+        PixelData<int> m2(1, 5, 1);
         m2.copyFromMeshWithUnaryOp(m, [](const int &a) { return a + 5; });
 
         for (size_t i = 0; i < m.mesh.size(); ++i) {
@@ -224,10 +224,10 @@ namespace {
 
     TEST(MeshDataSimpleTest, DownSample) {
         {   // reduce/constant_operator calculate maximum value when downsampling
-            MeshData<int> m(5, 6, 4);
+            PixelData<int> m(5, 6, 4);
             for (size_t i = 0; i < m.mesh.size(); ++i) m.mesh[i] = i + 1;
 
-            MeshData<int> m2;
+            PixelData<int> m2;
             downsample(m, m2,
                        [](float x, float y) { return std::max(x, y); },
                        [](float x) { return x; },
@@ -238,10 +238,10 @@ namespace {
             }
         }
         {   // reduce/constant_operator calculate maximum value when downsampling
-            MeshData<int> m(5, 6, 3);
+            PixelData<int> m(5, 6, 3);
             for (size_t i = 0; i < m.mesh.size(); ++i) m.mesh[i] = 5 * 6 * 3 - i;
 
-            MeshData<int> m2;
+            PixelData<int> m2;
             downsample(m, m2,
                        [](float x, float y) { return std::max(x, y); },
                        [](float x) { return x; },
@@ -254,10 +254,10 @@ namespace {
         }
         {
             // reduce/constant_operator calculate average value of pixels when downsampling
-            MeshData<uint16_t> m(2, 2, 2);
+            PixelData<uint16_t> m(2, 2, 2);
             for (size_t i = 0; i < m.mesh.size(); ++i) m.mesh[i] = 8 - i;
 
-            MeshData<float> m2;
+            PixelData<float> m2;
             downsample(m, m2,
                        [](const float &x, const float &y) -> float { return x + y; },
                        [](const float &x) -> float { return x / 8.0; },
@@ -270,10 +270,10 @@ namespace {
     }
 
     TEST(MeshDataSimpleTest, DownSamplePyramid) {
-        MeshData<float> m(4, 4, 4);
+        PixelData<float> m(4, 4, 4);
         for (size_t i = 0; i < m.mesh.size(); ++i) m.mesh[i] = i + 1;
 
-        std::vector<MeshData<float>> ds;
+        std::vector<PixelData<float>> ds;
         downsamplePyrmaid(m, ds, 3, 1);
 
         ASSERT_EQ(ds.size(), 4);
@@ -288,7 +288,7 @@ namespace {
     }
 
     TEST(MeshDataSimpleTest, GetIdx) {
-        MeshData<int> m(5, 6, 3);
+        PixelData<int> m(5, 6, 3);
 
         ASSERT_STREQ(m.getStrIndex(0).c_str(), "(0, 0, 0)");
         ASSERT_STREQ(m.getStrIndex(60).c_str(), "(0, 0, 2)");
diff --git a/test/TiffTest.cpp b/test/TiffTest.cpp
index f8dbcea0..ff4ea184 100644
--- a/test/TiffTest.cpp
+++ b/test/TiffTest.cpp
@@ -11,7 +11,7 @@ namespace {
     }
 
     TEST(TiffTest, LoadUint8) {
-        const MeshData<uint8_t> mesh = TiffUtils::getMesh<uint8_t>(testFilesDirectory() + "files/tiffTest/4x3x2x8bit.tif");
+        const PixelData<uint8_t> mesh = TiffUtils::getMesh<uint8_t>(testFilesDirectory() + "files/tiffTest/4x3x2x8bit.tif");
         for (int i = 0; i < 24; ++i) {
             ASSERT_EQ(mesh.mesh[i], i + 1);
         }
@@ -21,7 +21,7 @@ namespace {
         TiffUtils::TiffInfo t1(testFilesDirectory() + "files/tiffTest/4x3x2x8bit.tif");
         std::cout << t1 << std::endl;
 
-        MeshData<uint8_t> meshIn(t1.iImgHeight, t1.iImgWidth, t1.iNumberOfDirectories);
+        PixelData<uint8_t> meshIn(t1.iImgHeight, t1.iImgWidth, t1.iNumberOfDirectories);
         TiffUtils::getMesh<uint8_t>(t1, meshIn);
         for (int i = 0; i < 24; ++i) {
             ASSERT_EQ(meshIn.mesh[i], i + 1);
@@ -35,7 +35,7 @@ namespace {
         ASSERT_STREQ(t1.toString().c_str(), ("FileName: [" + fileName + "], Width/Height/Depth: 3/2/4, SamplesPerPixel: 1, Bits per sample: 16, ImageType: uint16, Photometric: 1, StripSize: 12").c_str());
         ASSERT_EQ(t1.isFileOpened(), true);
 
-        const MeshData<uint16_t> &mesh = TiffUtils::getMesh<uint16_t>(t1);
+        const PixelData<uint16_t> &mesh = TiffUtils::getMesh<uint16_t>(t1);
         ASSERT_EQ(mesh.x_num, 2);
         ASSERT_EQ(mesh.y_num, 3);
         ASSERT_EQ(mesh.z_num, 4);
@@ -49,7 +49,7 @@ namespace {
         TiffUtils::TiffInfo t1(testFilesDirectory() + "files/tiffTest/2x4x3xfloat.tif");
         std::cout << t1 << std::endl;
 
-        const MeshData<float> &mesh = TiffUtils::getMesh<float>(t1);
+        const PixelData<float> &mesh = TiffUtils::getMesh<float>(t1);
         ASSERT_EQ(mesh.x_num, 4);
         ASSERT_EQ(mesh.y_num, 2);
         ASSERT_EQ(mesh.z_num, 3);
@@ -70,7 +70,7 @@ namespace {
         // Then reads it again and compares input file and save file if same
         typedef uint16_t ImgType;
         TiffUtils::TiffInfo t(testFilesDirectory() + "files/tiffTest/3x2x4x16bit.tif");
-        MeshData<ImgType> mesh = TiffUtils::getMesh<ImgType>(t);
+        PixelData<ImgType> mesh = TiffUtils::getMesh<ImgType>(t);
         ASSERT_EQ(t.isFileOpened(), true);
 
         std::string fileName = "/tmp/testAprTiffSave" + std::to_string(time(nullptr)) + ".tif";
@@ -78,7 +78,7 @@ namespace {
 
         TiffUtils::TiffInfo t2(fileName);
         ASSERT_EQ(t2.isFileOpened(), true);
-        const MeshData<ImgType> &mesh2 = TiffUtils::getMesh<ImgType>(t2);
+        const PixelData<ImgType> &mesh2 = TiffUtils::getMesh<ImgType>(t2);
 
         ASSERT_EQ(mesh.mesh.size(), mesh2.mesh.size());
         for (size_t i = 0; i < mesh.mesh.size(); ++i)
@@ -94,14 +94,14 @@ namespace {
         // Then reads it again and compares input file and save file if same
         TiffUtils::TiffInfo t(testFilesDirectory() + "files/tiffTest/2x4x3xfloat.tif");
         ASSERT_EQ(t.isFileOpened(), true);
-        MeshData<float> mesh = TiffUtils::getMesh<float>(t);
+        PixelData<float> mesh = TiffUtils::getMesh<float>(t);
 
         std::string fileName = "/tmp/testAprTiffSave" + std::to_string(time(nullptr)) + ".tif";
         TiffUtils::saveMeshAsTiffUint16(fileName, mesh);
 
         TiffUtils::TiffInfo t2(fileName);
         ASSERT_EQ(t2.isFileOpened(), true);
-        const MeshData<uint16_t> &mesh2 = TiffUtils::getMesh<uint16_t>(t2);
+        const PixelData<uint16_t> &mesh2 = TiffUtils::getMesh<uint16_t>(t2);
 
         ASSERT_EQ(mesh.mesh.size(), mesh2.mesh.size());
         for (size_t i = 0; i < mesh.mesh.size(); ++i)