Adding example files

example/Numerics/Sussman_redistancing/example_sussman_circle/Makefile

+include ../../../example.mk
+
+CC=mpic++
+
+LDIR =
+
+OBJ = main.o
+
+%.o: %.cpp
+	$(CC) -O3 -c --std=c++11 -o $@ $< $(INCLUDE_PATH)
+
+example_sussman_circle: $(OBJ)
+	$(CC) -o $@ $^ $(CFLAGS) $(LIBS_PATH) $(LIBS)
+
+all: example_sussman_circle
+
+run: all
+	mpirun -np 2 ./example_sussman_circle
+
+.PHONY: clean all run
+
+clean:
+	rm -f *.o *~ core example_sussman_circle
+

example/Numerics/Sussman_redistancing/example_sussman_circle/config.cfg

+[pack]
+files = main.cpp Makefile

example/Numerics/Sussman_redistancing/example_sussman_images_2D/Makefile

+include ../../../example.mk
+
+CC=mpic++
+
+LDIR =
+
+OBJ = main.o
+
+%.o: %.cpp
+	$(CC) -O3 -c --std=c++11 -o $@ $< $(INCLUDE_PATH)
+
+example_sussman_images: $(OBJ)
+	$(CC) -o $@ $^ $(CFLAGS) $(LIBS_PATH) $(LIBS)
+
+all: example_sussman_images
+
+run: all
+	mpirun -np 2 ./example_sussman_images
+
+.PHONY: clean all run
+
+clean:
+	rm -f *.o *~ core example_sussman_images
+

example/Numerics/Sussman_redistancing/example_sussman_images_2D/config.cfg

+[pack]
+files = input/circle.bin input/dolphin.bin input/face.bin input/triangle.bin input/size_circle.csv input/size_dolphin.csv input/size_face.csv input/size_triangle.csv main.cpp Makefile

example/Numerics/Sussman_redistancing/example_sussman_images_2D/input/circle.bin

+








































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































\ No newline at end of file

example/Numerics/Sussman_redistancing/example_sussman_images_2D/input/face.bin

+





































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































































\ No newline at end of file

example/Numerics/Sussman_redistancing/example_sussman_images_2D/main.cpp

+//
+// Created by jstark on 2020-05-18.
+//
+/**
+ * @file example_sussman_images_2D/main.cpp
+ * @page example_sussman_images_2D Images 2D
+ *
+ * # Load the geometrical object from a binary 2D image #
+ * In this example, we will:
+ * * Read a 2D binary image from a binary file (for a 2D image volume see @ref example_sussman_images_2D)
+ * * Build a 2D cartesian OpenFPM grid with same dimensions as the image (1 particle for each pixel in x and y) or
+ *   refined by arbitrary factor in dimension of choice (e.g. to get a isotropic grid)
+ * * Assign pixel value to a grid node property
+ * * Run the Sussman Redistancing and get the narrow band around the interface (see also: @ref example_sussman_circle
+ *   and example_sussman_sphere)
+ *
+ * Output:
+ * print on promt : Iteration, Change, Residual (see: #DistFromSol::change, #DistFromSol::residual)
+ * writes vtk and hdf5 files of:
+ * 1.) 2D grid with geometrical object pre-redistancing and post-redistancing (Phi_0 and Phi_SDF, respectively)
+ * 2.) particles on narrow band around interface.
+ **/
+
+/**
+ * @page example_sussman_images_2D Images 2D
+ *
+ * ## Include ## {#e2d_img_include}
+ *
+ * These are the header files that we need to include:
+ *
+ * @snippet examples/example_sussman_images_2D/main.cpp Include
+ *
+ */
+//! @cond [Include] @endcond
+// Include standard library header files
+#include <iostream>
+#include <typeinfo>
+#include <cmath>
+#include <cstdio>
+
+// Include header files from other libraries
+#include <boost/log/core.hpp>
+#include <boost/log/trivial.hpp>
+#include <boost/log/expressions.hpp>
+#include <boost/log/utility/setup/file.hpp>
+#include <boost/log/utility/setup/common_attributes.hpp>
+
+// Include OpenFPM header files
+#include "Vector/vector_dist.hpp"
+#include "Grid/grid_dist_id.hpp"
+#include "data_type/aggregate.hpp"
+#include "Decomposition/CartDecomposition.hpp"
+
+// Include redistancing header files
+#include "util/PathsAndFiles.hpp"
+#include "level_set/redistancing_Sussman/RedistancingSussman.hpp"
+#include "level_set/redistancing_Sussman/NarrowBand.hpp"
+#include "RawReader/InitGridWithPixel.hpp"
+//! @cond [Include] @endcond
+
+
+/**
+ * @page example_sussman_images_2D Images 2D
+ *
+ * ## Initialization, indices and output folder ## {#e2d_img_init}
+ *
+ * Again we start with
+ * * Initializing OpenFPM
+ * * Setting the output path and creating an output folder
+ * This time, we also set the input path and name of the binary image that we want to load onto the grid. For this
+ * example we provide 3 simple example images. The binary images have been converted into -1 / +1 values.
+ * A jupyter notebook how to do this can be found here:
+ * Optionally, we can define the grid dimensionality and some indices for better code readability later on.
+ * * \p x: First dimension
+ * * \p y: Second dimension
+ * * \p Phi_0_grid: Index of property that stores the initial level-set-function
+ * * \p Phi_SDF_grid: Index of property where the redistancing result should be written to
+ *
+ * @snippet examples/example_sussman_images_2D/main.cpp Initialization
+ *
+ */
+//! @cond [Initialization] @endcond
+
+int main(int argc, char* argv[])
+{
+//	std::string image_name = "circle";
+//  std::string image_name = "triangle";
+//  std::string image_name = "face";
+	std::string image_name = "dolphin";
+
+
+	//	initialize library
+	openfpm_init(&argc, &argv);
+	///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+	// Set current working directory, define output paths and create folders where output will be saved
+	std::string cwd                     = get_cwd();
+	const std::string path_output       = cwd + "/output_" + image_name;
+//	const std::string path_output       = cwd + "/output_face";
+
+	create_directory_if_not_exist(path_output);
+	
+	//////////////////////////////////////////////////////////////////////////////////////////////
+	// Now we set the input paths. We need a binary file with the pixel values and a csv file with the
+	// size of the stack (in #pixels / dimension)
+	const std::string path_input        ="../input/";
+	const std::string path_to_image     = path_input + image_name + ".bin";
+	const std::string path_to_size      = path_input + "size_" + image_name + ".csv";
+
+	/*
+	 * in case of 2D (single image):
+	 */
+	const unsigned int grid_dim         = 2;
+	// some indices
+	const size_t x                      = 0;
+	const size_t y                      = 1;
+	
+	const size_t Phi_0_grid             = 0;
+	const size_t Phi_SDF_grid           = 1;
+	
+	const size_t Phi_SDF_vd             = 0;
+	const size_t Phi_grad_vd            = 1;
+	const size_t Phi_magnOfGrad_vd      = 2;
+	//! @cond [Initialization] @endcond
+	
+	
+	/**
+	 * @page example_sussman_images_2D Images 2D
+	 *
+	 * ## Set refinement factor ## {#e2d_img_refine}
+	 *
+	 * If we want that the grid has a different resolution as the image, we can set a refinement factor. In the
+	 * refinement array, we define by which factor the grid resolution should be changed w.r.t. the image resolution.
+	 * This can be useful for example when you want to have an isotropic grid but the underlying image is
+	 * anisotropic (as it often happens for microcsopic z-stacks rather than 2D images).
+	 *
+	 * @snippet examples/example_sussman_images_2D/main.cpp Refinement
+	 *
+	 */
+	//! @cond [Refinement] @endcond
+	const double refinement []          = {0.5, 0.5};   // (2D) factor by which grid should be finer / coarser as
+	// underlying image in each dimension (e.g. to get isotropic grid from anisotropic image resolution)
+	//! @cond [Refinement] @endcond
+	///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+	//	read the stack size (number of pixel values per dimension) from a binary file
+	// alternatively, you can directly define the stack size as: std::vector<size_t> stack_size {#pixels in x, #pixels in y}
+	/**
+	 * @page example_sussman_images_2D Images 2D
+	 *
+	 * ## Get size from image_size.csv ## {#e2d_img_size}
+	 *
+	 * Here we read the image size (number of pixel values per dimension) from a csv file. Alternatively, you can
+	 * directly define the image size as: \p std::vector<size_t> \p stack_size \p {#pixels \p in \p x, \p #pixels \p in
+	 * \p y}.
+	 * We use this image size and the refinement factor to set the grid size \p sz.
+	 *
+	 *
+	 * @snippet examples/example_sussman_images_2D/main.cpp Size
+	 *
+	 */
+	//! @cond [Size] @endcond
+	std::vector<size_t> stack_size = get_size(path_to_size);
+	auto & v_cl = create_vcluster();
+	if (v_cl.rank() == 0)
+	{
+		for(std::vector<int>::size_type i = 0; i != stack_size.size(); i++)
+		{
+			std::cout << "#Pixel in dimension " << i << " = " << stack_size[i] << std::endl;
+			std::cout << "original refinement in dimension " << i << " = " << refinement[i] << std::endl;
+		}
+	}
+	
+	//  Array with grid dimensions after refinement. This size-array will be used for the grid creation.
+	const size_t sz[grid_dim] = {(size_t)std::round(stack_size[x] * refinement[x]), (size_t)std::round(stack_size[y] * refinement[y])}; // 2D
+	//! @cond [Size] @endcond
+	// Here we create a 2D grid that stores 2 properties:
+	// Prop1: store the initial Phi;
+	// Prop2: here the re-initialized Phi (signed distance function) will be written to in the re-distancing step
+	/**
+	 * @page example_sussman_images_2D Images 2D
+	 *
+	 * ## Run Sussman redistancing and get narrow band ## {#e2d_img_redist}
+	 *
+	 * Once we have loaded the geometrical object from the 2D image onto the grid, we can perform Sussman
+	 * redistancing and get the narrow band the same way as it is explained in detail here: @ref
+	 * example_sussman_circle and here: @ref example_sussman_sphere.
+	 *
+	 *
+	 * @snippet examples/example_sussman_images_2D/main.cpp Redistancing
+	 *
+	 */
+	//! @cond [Redistancing] @endcond
+	Box<grid_dim, double> box({0.0, 0.0}, {5.0, 5.0}); // 2D
+	
+	Ghost<grid_dim, long int> ghost(0);
+	typedef aggregate<double, double> props;
+	typedef grid_dist_id<grid_dim, double, props > grid_in_type;
+	grid_in_type g_dist(sz, box, ghost);
+	g_dist.setPropNames({"Phi_0", "Phi_SDF"});
+	
+	// Now we can initialize the grid with the pixel values from the image stack
+	load_pixel_onto_grid<Phi_0_grid>(g_dist, path_to_image, stack_size);
+	g_dist.write(path_output + "/grid_from_images_initial", FORMAT_BINARY); // Save the initial grid as vtk file
+
+
+	///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+	// Now we want to convert the initial Phi into a signed distance function (SDF) with magnitude of gradient = 1
+	// For the initial re-distancing we use the Sussman method
+	// 1.) Set some redistancing options
+	Redist_options redist_options;
+	redist_options.sigma                                = 0;        // if the initial gradient of phi at the interface is too large, gauss smoothing will automatically be applied until phi gradient magnitude <= 12, regardless of the given sigma in the main
+	redist_options.min_iter                             = 100;      // min. number of iterations before steady state in narrow band will be checked (default: 100)
+	redist_options.max_iter                             = 10000;    // max. number of iterations you want to run the
+																	// redistancing, even if steady state might not yet have been reached (default: 1e6)
+	
+	redist_options.convTolChange.value                  = 1e-6;     // convolution tolerance for the normalized total change of Phi in the narrow band between two consecutive iteratins (default: 1e-6)
+	redist_options.convTolChange.check                  = true;     // define here which of the convergence criteria above should be used. If both are true, termination only occurs when both are fulfilled or when iter > max_iter
+	redist_options.convTolResidual.value                = 1e-1;     // convolution tolerance for the normalized total residual of Phi versus the SDF (aka the error). Don't choose too small to not run forever. (default: 1e-1)
+	redist_options.convTolResidual.check                = false;    // (default: false)
+	
+	redist_options.interval_check_convergence           = 1;        // interval of #iterations at which convergence is checked (default: 100)
+	redist_options.width_NB_in_grid_points              = 6;        // width of narrow band in number of grid points.
+	// Must be at least 4, in order to have at least 2 grid points on each side of the interface. (default: 4)
+	redist_options.print_current_iterChangeResidual     = true;     // if true, prints out every current iteration + corresponding change from the previous iteration + residual from SDF (default: false)
+	redist_options.print_steadyState_iter               = true;     // if true, prints out the final iteration number when steady state was reached + final change + residual (default: true)
+	
+	RedistancingSussman<grid_in_type> redist_obj(g_dist, redist_options);   // Instantiation of Sussman-redistancing class
+//	std::cout << "dt = " << redist_obj.get_time_step() << std::endl;
+	// Run the redistancing. in the <> brackets provide property-index where 1.) your initial Phi is stored and 2.) where the resulting SDF should be written to.
+	redist_obj.run_redistancing<Phi_0_grid, Phi_SDF_grid>();
+	
+	g_dist.write(path_output + "/grid_images_post_redistancing", FORMAT_BINARY); // Save the result as vtk file
+	
+	///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+	//	Get narrow band: Place particles on interface (narrow band width e.g. 2 grid points on each side of the interface)
+	size_t bc[grid_dim] = {PERIODIC, PERIODIC};
+	// Create an empty vector to which narrow-band particles will be added. You can choose, how many properties you want.
+	// Minimum is 1 property, to which the Phi_SDF can be written
+	// In this example we chose 3 properties. The 1st for the Phi_SDF, the 2nd for the gradient of phi and the 3rd for
+	// the magnitude of the gradient
+	typedef aggregate<double, double[grid_dim], double> props_nb;
+	typedef vector_dist<grid_dim, double, props_nb> vd_type;
+	vd_type vd_narrow_band(0, box, bc, ghost);
+	vd_narrow_band.setPropNames({"Phi_SDF", "Phi_grad", "Phi_magnOfGrad"});
+	
+	NarrowBand<grid_in_type> narrowBand(g_dist, redist_options.width_NB_in_grid_points); // Instantiation of NarrowBand class
+	
+	// Get the narrow band. You can decide, if you only want the Phi_SDF saved to your particles or
+	// if you also want the gradients or gradients and magnitude of gradient.
+	// The function will know depending on how many property-indices you provide in the <> brackets.
+	// First property-index must always be the index of the SDF on the grid!
+	// E.g.: The following line would only write only the Phi_SDF from the grid to your narrow-band particles
+	// narrowBand.get_narrow_band<Phi_SDF_grid, Phi_SDF_vd>(g_dist, vd_narrow_band);
+	// Whereas this will give you the gradients and magnOfGrad of Phi as well:
+	narrowBand.get_narrow_band<Phi_SDF_grid, Phi_SDF_vd, Phi_grad_vd, Phi_magnOfGrad_vd>(g_dist, vd_narrow_band);
+	
+	vd_narrow_band.write(path_output + "/vd_narrow_band_images", FORMAT_BINARY); // Save particles as vtk file
+
+	openfpm_finalize(); // Finalize openFPM library
+	return 0;
+}
+//! @cond [Redistancing] @endcond
+
+/**
+ * @page example_sussman_images_2D Images 2D
+ *
+ * ## Full code ## {#e2d_img_full}
+ *
+ * @include examples/example_sussman_images_2D/main.cpp
+ */
+
+
+
+
+
+
+
+
+

example/Numerics/Sussman_redistancing/example_sussman_images_3D/Makefile

+include ../../../example.mk
+
+CC=mpic++
+
+LDIR =
+
+OBJ = main.o
+
+%.o: %.cpp
+	$(CC) -O3 -c --std=c++11 -o $@ $< $(INCLUDE_PATH)
+
+example_sussman_images: $(OBJ)
+	$(CC) -o $@ $^ $(CFLAGS) $(LIBS_PATH) $(LIBS)
+
+all: example_sussman_images
+
+run: all
+	mpirun -np 2 ./example_sussman_images
+
+.PHONY: clean all run
+
+clean:
+	rm -f *.o *~ core example_sussman_images
+

example/Numerics/Sussman_redistancing/example_sussman_images_3D/config.cfg

+[pack]
+files = input/sphere.bin  input/size_sphere.csv main.cpp Makefile

example/Numerics/Sussman_redistancing/example_sussman_images_3D/main.cpp

+//
+// Created by jstark on 2020-05-18.
+//
+/**
+ * @file example_sussman_images_3D/main.cpp
+ * @page example_sussman_images_3D Images 3D
+ *
+ * [TOC]
+ *
+ * # Example for loading a 3D object from an image stack (binary) onto a grid and applying Sussman redistancing #
+ *
+ * In this example the image stack is read from a binary file. A 3D cartesian grid with same dimensions as image
+ * stack is constructed. The grid resolution can be either 1 grid node for each pixel in x and y) or the resolution
+ * can be higher/lower as the image stack. This can be achieved by setting the refinement factor to a value of choice in
+ * dimension of choice (e.g. to get a isotropic grid). The pixel value is stored in a property of the grid.
+ *
+ * The initial, image-derived level set function phi (indicator function) is then converted into a signed distance
+ * function by using the Sussman redistancing (see @ref RedistancingSussman.hpp).
+ *
+ * When the grid contains phi as a signed distance function, particles can be placed on narrow band around the
+ * interface.
+ *
+ * Output:
+ * print on promt : Iteration, Change, Residual (see: #DistFromSol::change, #DistFromSol::residual)
+ * writes vtk and hdf5 files of:
+ * 1.) 3D grid with geometrical object pre-redistancing and post-redistancing (Phi_0 and Phi_SDF, respectively)
+ * 2.) particles on narrow band around interface.
+ **/
+
+/**
+ * @page example_sussman_images_3D Images 3D
+ *
+ * ## Include ## {#e2d_img_include}
+ *
+ * These are the header files that we need to include:
+ *
+ * @snippet examples/example_sussman_images_3D/main.cpp Include
+ *
+ */
+//! @cond [Include] @endcond
+// Include standard library header files
+#include <iostream>
+#include <typeinfo>
+#include <cmath>
+#include <cstdio>
+
+// Include header files from other libraries
+#include <boost/log/core.hpp>
+#include <boost/log/trivial.hpp>
+#include <boost/log/expressions.hpp>
+#include <boost/log/utility/setup/file.hpp>
+#include <boost/log/utility/setup/common_attributes.hpp>
+
+// Include OpenFPM header files
+#include "Vector/vector_dist.hpp"
+#include "Grid/grid_dist_id.hpp"
+#include "data_type/aggregate.hpp"
+#include "Decomposition/CartDecomposition.hpp"
+
+// Include redistancing header files
+#include "util/PathsAndFiles.hpp"
+#include "level_set/redistancing_Sussman/RedistancingSussman.hpp"
+#include "level_set/redistancing_Sussman/NarrowBand.hpp"
+#include "RawReader/InitGridWithPixel.hpp"
+//! @cond [Include] @endcond
+
+
+/**
+ * @page example_sussman_images_3D Images 3D
+ *
+ * ## Initialization, indices and output folder ## {#e2d_img_init}
+ *
+ * Again we start with
+ * * Initializing OpenFPM
+ * * Setting the output path and creating an output folder
+ * This time, we also set the input path and name of the binary image stack that we want to load onto the grid. For this
+ * example we provide 1 simple example image stack. The binary image stack have been converted into -1 / +1 values.
+ * A jupyter notebook how to do this can be found here:
+ * Optionally, we can define the grid dimensionality and some indices for better code readability later on.
+ * * \p x: First dimension
+ * * \p y: Second dimension
+ * * \p Phi_0_grid: Index of property that stores the initial level-set-function
+ * * \p Phi_SDF_grid: Index of property where the redistancing result should be written to
+ *
+ * @snippet examples/example_sussman_images_3D/main.cpp Initialization
+ *
+ */
+//! @cond [Initialization] @endcond
+int main(int argc, char* argv[])
+{
+	//	initialize library
+	openfpm_init(&argc, &argv);
+	///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+	// Set current working directory, define output paths and create folders where output will be saved
+	std::string cwd                     = get_cwd();
+	const std::string path_output       = cwd + "/output_3D_sphere";
+	create_directory_if_not_exist(path_output);
+	
+	//////////////////////////////////////////////////////////////////////////////////////////////
+	// Now we set the input paths. We need a binary file with the pixel values and a csv file with the
+	// size of the stack (in #pixels / dimension)
+	const std::string path_input        ="../input";
+	const std::string path_to_zstack    = path_input + "/sphere.bin";
+	const std::string path_to_size      = path_input + "/size_sphere.csv";
+	/*
+	 * in case of 3D (image volume):
+	 */
+	const unsigned int grid_dim         = 3;
+	const size_t x                      = 0;
+	const size_t y                      = 1;
+	const size_t z                      = 2;
+
+	const size_t Phi_0_grid             = 0;
+	const size_t Phi_SDF_grid           = 1;
+
+	const size_t Phi_SDF_vd             = 0;
+	const size_t Phi_grad_vd            = 1;
+	const size_t Phi_magnOfGrad_vd      = 2;
+	//! @cond [Initialization] @endcond
+	
+	
+	/**
+	 * @page example_sussman_images_3D Images 3D
+	 *
+	 * ## Set refinement factor ## {#e2d_img_refine}
+	 *
+	 * If we want that the grid has a different resolution as the image stack, we can set a refinement factor. In the
+	 * refinement array, we define by which factor the grid resolution should be changed w.r.t. the image stack
+	 * resolution.
+	 * This can be useful for example when you want to have an isotropic grid but the underlying image stack is
+	 * anisotropic (as it often happens for the z-resolution of volumetric microscopy image data).
+	 *
+	 * @snippet examples/example_sussman_images_3D/main.cpp Refinement
+	 *
+	 */
+	//! @cond [Refinement] @endcond
+	const double refinement []          = {1.0, 1.0, 1.0};      // without refinement
+//	const double refinement []          = {0.8, 1.5, 2.0};      // factors by which grid should be finer as underlying image stack in each dimension (e.g. to get isotropic grid from anisotropic stack resolution)
+	//! @cond [Refinement] @endcond
+	///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+	//	read the stack size (number of pixel values per dimension) from a binary file
+	// alternatively, you can directly define the stack size as: std::vector<size_t> stack_size {#pixels in x, #pixels in y, #pixels in z}
+	/**
+	 * @page example_sussman_images_3D Images 3D
+	 *
+	 * ## Get size from image_size.csv ## {#e2d_img_size}
+	 *
+	 * Here we read the stack size (number of pixel values per dimension) from a csv file. Alternatively, you can
+	 * directly define the stack size as: \p std::vector<size_t> \p stack_size \p {#pixels \p in \p x, \p #pixels \p in
+	 * \p y}.
+	 * We use this stack size and the refinement factor to set the grid size \p sz.
+	 *
+	 *
+	 * @snippet examples/example_sussman_images_3D/main.cpp Size
+	 *
+	 */
+	//! @cond [Size] @endcond
+	std::vector<size_t> stack_size = get_size(path_to_size);
+	auto & v_cl = create_vcluster();
+	if (v_cl.rank() == 0)
+	{
+ 		for(std::vector<int>::size_type i = 0; i != stack_size.size(); i++)
+		{
+			std::cout << "#Pixel in dimension " << i << " = " << stack_size[i] << std::endl;
+			std::cout << "original refinement in dimension " << i << " = " << refinement[i] << std::endl;
+		}
+	}
+	
+	
+	//  Array with grid dimensions after refinement. This size-array will be used for the grid creation.
+	const size_t sz[grid_dim] = {(size_t)std::round(stack_size[x] * refinement[x]), (size_t)std::round(stack_size[y] * refinement[y]), (size_t)std::round(stack_size[z] * refinement[z])}; // 3D
+	//! @cond [Size] @endcond
+	// Here we create a 3D grid that stores 2 properties:
+	// Prop1: store the initial Phi;
+	// Prop2: here the re-initialized Phi (signed distance function) will be written to in the re-distancing step
+	/**
+	 * @page example_sussman_images_3D Images 3D
+	 *
+	 * ## Run Sussman redistancing and get narrow band ## {#e2d_img_redist}
+	 *
+	 * Once we have loaded the geometrical object from the 3D stack onto the grid, we can perform Sussman
+	 * redistancing and get the narrow band the same way as it is explained in detail here: @ref
+	 * example_sussman_circle and here: @ref example_sussman_sphere.
+	 *
+	 *
+	 * @snippet examples/example_sussman_images_3D/main.cpp Redistancing
+	 *
+	 */
+	//! @cond [Redistancing] @endcond
+	Box<grid_dim, double> box({0.0, 0.0, 0.0}, {5.0, 5.0, 5.0}); // 3D
+
+	Ghost<grid_dim, long int> ghost(0);
+	typedef aggregate<double, double> props;
+	typedef grid_dist_id<grid_dim, double, props > grid_in_type;
+	grid_in_type g_dist(sz, box, ghost);
+	g_dist.setPropNames({"Phi_0", "Phi_SDF"});
+	
+	// Now we can initialize the grid with the pixel values from the image stack
+	load_pixel_onto_grid<Phi_0_grid>(g_dist, path_to_zstack, stack_size);
+	g_dist.write(path_output + "/grid_from_images_initial", FORMAT_BINARY); // Save the initial grid as vtk file
+
+
+	///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+	// Now we want to convert the initial Phi into a signed distance function (SDF) with magnitude of gradient = 1
+	// For the initial re-distancing we use the Sussman method
+	// 1.) Set some redistancing options
+	Redist_options redist_options;
+	redist_options.sigma                                = 0;        // if the initial gradient of phi at the interface is too large, gauss smoothing will automatically be applied until phi gradient magnitude <= 12, regardless of the given sigma in the main
+	redist_options.min_iter                             = 100;      // min. number of iterations before steady state in narrow band will be checked (default: 100)
+	redist_options.max_iter                             = 10000;    // max. number of iterations you want to run the
+															        // redistancing, even if steady state might not yet
+															        // have been reached (default: 1e6)
+	redist_options.convTolChange.value                  = 1e-6;     // convolution tolerance for the normalized total change of Phi in the narrow band between two consecutive iteratins (default: 1e-6)
+	redist_options.convTolChange.check                  = true;     // define here which of the convergence criteria above should be used. If both are true, termination only occurs when both are fulfilled or when iter > max_iter
+	redist_options.convTolResidual.value                = 1e-1;     // convolution tolerance for the normalized total residual of Phi versus the SDF (aka the error). Don't choose too small to not run forever. (default: 1e-1)
+	redist_options.convTolResidual.check                = false;    // (default: false)
+	
+	redist_options.interval_check_convergence           = 1;        // interval of #iterations at which convergence is checked (default: 100)
+	redist_options.width_NB_in_grid_points              = 6;        // width of narrow band in number of grid points.
+															        // Must be at least 4, in order to have at least 2
+															        // grid points on each side of the interface.
+															        // (default: 4)
+	redist_options.print_current_iterChangeResidual     = true;     // if true, prints out every current iteration + corresponding change from the previous iteration + residual from SDF (default: false)
+	redist_options.print_steadyState_iter               = true;     // if true, prints out the final iteration number when steady state was reached + final change + residual (default: true)
+	
+	RedistancingSussman<grid_in_type> redist_obj(g_dist, redist_options);   // Instantiation of Sussman-redistancing class
+//	std::cout << "dt = " << redist_obj.get_time_step() << std::endl;
+	// Run the redistancing. in the <> brackets provide property-index where 1.) your initial Phi is stored and 2.) where the resulting SDF should be written to.
+	redist_obj.run_redistancing<Phi_0_grid, Phi_SDF_grid>();
+	
+	g_dist.write(path_output + "/grid_images_post_redistancing", FORMAT_BINARY); // Save the result as vtk file
+	
+	///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+	//	Get narrow band: Place particles on interface (narrow band width e.g. 2 grid points on each side of the interface)
+	size_t bc[grid_dim] = {PERIODIC, PERIODIC, PERIODIC};
+	// Create an empty vector to which narrow-band particles will be added. You can choose, how many properties you want.
+	// Minimum is 1 property, to which the Phi_SDF can be written
+	// In this example we chose 3 properties. The 1st for the Phi_SDF, the 2nd for the gradient of phi and the 3rd for
+	// the magnitude of the gradient
+	typedef aggregate<double, double[grid_dim], double> props_nb;
+	typedef vector_dist<grid_dim, double, props_nb> vd_type;
+	vd_type vd_narrow_band(0, box, bc, ghost);
+	vd_narrow_band.setPropNames({"Phi_SDF", "Phi_grad", "Phi_magnOfGrad"});
+	
+	NarrowBand<grid_in_type> narrowBand(g_dist, redist_options.width_NB_in_grid_points); // Instantiation of NarrowBand class
+	
+	// Get the narrow band. You can decide, if you only want the Phi_SDF saved to your particles or
+	// if you also want the gradients or gradients and magnitude of gradient.
+	// The function will know depending on how many property-indices you provide in the <> brackets.
+	// First property-index must always be the index of the SDF on the grid!
+	// E.g.: The following line would only write only the Phi_SDF from the grid to your narrow-band particles
+	// narrowBand.get_narrow_band<Phi_SDF_grid, Phi_SDF_vd>(g_dist, vd_narrow_band);
+	// Whereas this will give you the gradients and magnOfGrad of Phi as well:
+	narrowBand.get_narrow_band<Phi_SDF_grid, Phi_SDF_vd, Phi_grad_vd, Phi_magnOfGrad_vd>(g_dist, vd_narrow_band);
+	
+	vd_narrow_band.write(path_output + "/vd_narrow_band_images", FORMAT_BINARY); // Save particles as vtk file
+
+	openfpm_finalize(); // Finalize openFPM library
+	return 0;
+}
+//! @cond [Redistancing] @endcond
+
+/**
+ * @page example_sussman_images_3D Images 3D
+ *
+ * ## Full code ## {#e2d_img_full}
+ *
+ * @include examples/example_sussman_images_3D/main.cpp
+ */
+
+
+
+
+
+
+
+
+
+

example/Numerics/Sussman_redistancing/example_sussman_sphere/Makefile

+include ../../../example.mk
+
+CC=mpic++
+
+LDIR =
+
+OBJ = main.o
+
+%.o: %.cpp
+	$(CC) -O3 -c --std=c++11 -o $@ $< $(INCLUDE_PATH)
+
+example_sussman_sphere: $(OBJ)
+	$(CC) -o $@ $^ $(CFLAGS) $(LIBS_PATH) $(LIBS)
+
+all: example_sussman_sphere
+
+run: all
+	mpirun -np 2 ./example_sussman_sphere
+
+.PHONY: clean all run
+
+clean:
+	rm -f *.o *~ core example_sussman_sphere
+

example/Numerics/Sussman_redistancing/example_sussman_sphere/config.cfg

+[pack]
+files = main.cpp Makefile

example/Numerics/Sussman_redistancing/image_binary_conversion/config.cfg

+[pack]
+files = image2binary_circle.ipynb image2binary_create_sphere.ipynb image2binary_dolphin.ipynb lib/image_process_and_save.py lib/morphological_operations.py

example/Numerics/Sussman_redistancing/image_binary_conversion/image2binary_circle.ipynb

+{
+ "cells": [
+  {
+   "cell_type": "code",
+   "execution_count": 1,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from tifffile import imread, imsave\n",
+    "import numpy as np\n",
+    "import os\n",
+    "import matplotlib.pyplot as plt\n",
+    "import cv2\n",
+    "import sys\n",
+    "sys.path.append(\"./lib\")\n",
+    "import morphological_operations as mo\n",
+    "import image_process_and_save as ips"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Open and plot image"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 2,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "filename = '/Users/jstark/Desktop/image_examples/circle.tiff'\n",
+    "x = imread('%s' % filename)\n",
+    "x = x[:, :, 0]"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 3,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "(255, 0, (64, 64))"
+      ]
+     },
+     "execution_count": 3,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "np.max(x), np.min(x), x.shape"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 4,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# fig, ax = plt.subplots(figsize=(7, 7))\n",
+    "# ax.imshow(x, cmap='gray')"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Convert binary image: >1 -> -1, 0 -> +1 (indicator function)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 5,
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "format of stack x_bin =  (64, 64) \n",
+      "\n",
+      "minimum pixel value after conversion =  -1 \n",
+      "\n",
+      "maximum pixel value after conversion =  1 \n",
+      "\n"
+     ]
+    },
+    {
+     "data": {
+      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAaIAAAGfCAYAAAAd79YcAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAEj5JREFUeJzt3V+sZWd53/HfEw8WKaGxDenI8uDaCAvERTBoREFBleuIyE1R8AWyiBJpEqHOTSpRJVUwuUkTCQl6EeKLKpJlSOYiCVgOxBYXSSzHaXPlMMZEGA8Ehxoxlu1JZFtNckFleHpxlttjZ2bOnjl/npm9Px9pdPZaZ5+939dnj7+z1nrPPtXdAYApPzQ9AAA2mxABMEqIABglRACMEiIARgkRAKOECIBRQgTAKCECYNSh3XxxVd2W5K4kVyS5p7s/ucP9vY0DwIbo7lrlfnWxb/FTVVck+Zsk709yOsmXk/xsdz9xnq8RIoANsWqIdnNq7t1Jnuzub3f3/0nyuSQf3MXjAbCBdhOi65J8d9v26WXfK1TV8ao6WVUnd/FcAKypXV0jWkV3353k7sSpOQD+ud0cET2d5E3bto8s+wBgZbsJ0ZeT3FRVN1bVlUk+nOSBvRkWAJviok/NdfdLVfWfkvxptpZvf7a7v75nIwNgI1z08u2LejLXiAA2xkEs3waAXRMiAEYJEQCjhAiAUUIEwCghAmCUEAEwSogAGCVEAIwSIgBGCREAo4QIgFFCBMAoIQJglBABMEqIABglRACMEiIARgkRAKOECIBRQgTAKCECYJQQATBKiAAYJUQAjBIiAEYJEQCjhAiAUUIEwCghAmCUEAEwSogAGCVEAIwSIgBGCREAo4QIgFFCBMAoIQJglBABMEqIABglRACMEiIARgkRAKOECIBRQgTAKCECYJQQATBKiAAYJUQAjBIiAEbtGKKq+mxVnamqx7ftu6aqHqyqby0fr97fYQKwrlY5Ivq9JLe9at+dSR7q7puSPLRsA8AF2zFE3f0/kzz/qt0fTHJiuX0iye17PC4ANsTFXiM63N3PLLefTXJ4j8YDwIY5tNsH6O6uqj7X56vqeJLju30eANbTxR4RPVdV1ybJ8vHMue7Y3Xd399HuPnqRzwXAGrvYED2Q5Nhy+1iS+/dmOABsmuo+51m1rTtU/WGSW5K8MclzSX49yR8nuTfJ9Um+k+SO7n71goazPdb5nwyAtdHdtcr9dgzRXhIigM2xaoi8swIAo4QIgFFCBMAoIQJglBABMEqIABglRACMEiIARgkRAKN2/e7bcDk5yHcSmVa10g+1wzhHRACMEiIARgkRAKOECIBRQgTAKKvmuCxs0mq3vbJX/82svmO/OSICYJQQATBKiAAYJUQAjBIiAEZZNccIq+AuHxf6vbLKjgvliAiAUUIEwCghAmCUEAEwSogAGGXVHPvK6rjNc67vudV0nIsjIgBGCREAo4QIgFFCBMAoIQJglFVzXBCr4LhY3rOOc3FEBMAoIQJglBABMEqIABglRACMsmqOs7I6jmnes25zOCICYJQQATBKiAAYJUQAjBIiAEZZNbfhrI7jcmM13fpxRATAKCECYJQQATBKiAAYJUQAjNoxRFX1pqp6uKqeqKqvV9VHl/3XVNWDVfWt5ePV+z9cANZN7bR8t6quTXJtd3+lql6f5NEktyf5hSTPd/cnq+rOJFd398d2eCxrhS8xlm+zLizfvvR090rflB2PiLr7me7+ynL7H5KcSnJdkg8mObHc7US24gQAF+SCfqC1qm5I8s4kjyQ53N3PLJ96Nsnhc3zN8STHL36IAKyzHU/N/b87Vv1Ikv+R5BPd/YWqerG7r9r2+Re6+7zXiZyau/Q4Nce6cGru0rNnp+aSpKpek+SPkvx+d39h2f3ccv3o5etIZy5moABstlVWzVWSzyQ51d2/te1TDyQ5ttw+luT+vR8ee6W7z/oH1oXX+OVrlVVz70vyl0m+luQHy+5fy9Z1onuTXJ/kO0nu6O7nd3gsr4oh/kKyqZyym7PqqbmVrxHtBSGaI0RsKiGas6fXiABgvwgRAKOECIBRfkPrmnEtCF7Jb3S99DkiAmCUEAEwSogAGCVEAIwSIgBGCREAo4QIgFFCBMAoIQJglBABMEqIABjlveYuU95TDnbHe9BdOhwRATBKiAAYJUQAjBIiAEYJEQCjhAiAUUIEwCghAmCUEAEwSogAGCVEAIwSIgBGCREAo4QIgFFCBMAoIQJglBABMMpvaL3E+U2scLD85taD54gIgFFCBMAoIQJglBABMEqIABglRACMEiIARgkRAKOECIBRQgTAKCECYJQQATBKiAAYJUQAjBIiAEYJEQCjdgxRVb22qv6qqv66qr5eVb+x7L+xqh6pqier6vNVdeX+DxeAdbPKEdH3ktza3e9IcnOS26rqPUk+leTT3f2WJC8k+cj+DROAdbVjiHrLPy6br1n+dJJbk9y37D+R5PZ9GSEAa22la0RVdUVVfTXJmSQPJvnbJC9290vLXU4nuW5/hgjAOlspRN39/e6+OcmRJO9O8rZVn6CqjlfVyao6eZFjBGCNXdCque5+McnDSd6b5KqqOrR86kiSp8/xNXd399HuPrqrkQKwllZZNfdjVXXVcvuHk7w/yalsBelDy92OJbl/vwYJwPqq7j7/Hap+PFuLEa7IVrju7e7frKo3J/lckmuSPJbk57v7ezs81vmfjH9mp+8PcDCqanoIl53uXuk/2o4h2ktCdOGECC4NQnThVg2Rd1YAYJQQATBKiAAYJUQAjBIiAEYJEQCjhAiAUUIEwCghAmCUEAEwSogAGCVEAIwSIgBGCREAo4QIgFFCBMAoIQJglBABMEqIABglRACMEiIARgkRAKOECIBRQgTAKCECYJQQATBKiAAYJUQAjBIiAEYJEQCjhAiAUUIEwCghAmDUoekBcH5Vddb93X3AI4HNcK6/c+wfR0QAjBIiAEYJEQCjhAiAUUIEwCghAmCUEAEwSogAGCVEAIwSIgBGCREAo4QIgFFCBMAoIQJglBABMEqIABi1coiq6oqqeqyqvrRs31hVj1TVk1X1+aq6cv+GCcC6upAjoo8mObVt+1NJPt3db0nyQpKP7OXAOL+qOusfYDX+Dl06VgpRVR1J8h+S3LNsV5Jbk9y33OVEktv3Y4AArLdVj4h+O8mvJvnBsv2GJC9290vL9ukk1+3x2ADYADuGqKo+kORMdz96MU9QVcer6mRVnbyYrwdgvR1a4T4/keRnquqnk7w2yb9McleSq6rq0HJUdCTJ02f74u6+O8ndSVJVvSejBmBt7HhE1N0f7+4j3X1Dkg8n+fPu/rkkDyf50HK3Y0nu37dRArC2dvNzRB9L8stV9WS2rhl9Zm+GBMAmqe6DO1vm1Nz+O8jvJ1zOLNXef9290n9k76wAwCghAmCUEAEwSogAGLXKzxFxGTnXBViLGNhUFiVc+hwRATBKiAAYJUQAjBIiAEYJEQCjrJrbEFbTse6sjrt8OSICYJQQATBKiAAYJUQAjBIiAEYJEQCjhAiAUUIEwCghAmCUEAEwSogAGOW95jac96DjcuM95daPIyIARgkRAKOECIBRQgTAKCECYJRVc5yV1XRMszpuczgiAmCUEAEwSogAGCVEAIwSIgBGWTXHBbnQlUxW2fEyq+A4F0dEAIwSIgBGCREAo4QIgFFCBMAoq+bYV96zbvNYHceFckQEwCghAmCUEAEwSogAGCVEAIyyao4R3rPu8mEVHPvNEREAo4QIgFFCBMAoIQJg1EqLFarqqST/kOT7SV7q7qNVdU2Szye5IclTSe7o7hf2Z5gArKtaZTXSEqKj3f332/b9tyTPd/cnq+rOJFd398d2eBxLnxi1SavvrHZjWnev9CLczam5DyY5sdw+keT2XTwWABtq1RB1kj+rqker6viy73B3P7PcfjbJ4T0fHQBrb9UfaH1fdz9dVf8qyYNV9Y3tn+zuPtdptyVcx8/2OQBY6RrRK76g6r8m+cck/zHJLd39TFVdm+QvuvutO3zt5pyg55LkGhEcnD27RlRVr6uq1798O8lPJXk8yQNJji13O5bk/osbKgCbbMcjoqp6c5IvLpuHkvxBd3+iqt6Q5N4k1yf5TraWbz+/w2Ntzj9HATbcqkdEF3xqbjeECGBzHMTybQDYNSECYJQQATBKiAAYJUQAjBIiAEYJEQCjhAiAUUIEwCghAmCUEAEwSogAGCVEAIwSIgBGCREAo4QIgFFCBMAoIQJglBABMEqIABglRACMEiIARgkRAKOECIBRQgTAKCECYJQQATBKiAAYJUQAjBIiAEYJEQCjhAiAUUIEwCghAmCUEAEwSogAGCVEAIwSIgBGCREAo4QIgFFCBMAoIQJglBABMEqIABglRACMEiIARgkRAKOECIBRQgTAKCECYNRKIaqqq6rqvqr6RlWdqqr3VtU1VfVgVX1r+Xj1fg8WgPWz6hHRXUn+pLvfluQdSU4luTPJQ919U5KHlm0AuCDV3ee/Q9WPJvlqkjf3tjtX1TeT3NLdz1TVtUn+orvfusNjnf/JAFgb3V2r3G+VI6Ibk/xdkt+tqseq6p6qel2Sw939zHKfZ5McvrihArDJVgnRoSTvSvI73f3OJP+UV52GW46Uznq0U1XHq+pkVZ3c7WABWD+rhOh0ktPd/ciyfV+2wvTcckouy8czZ/vi7r67u49299G9GDAA62XHEHX3s0m+W1UvX//5ySRPJHkgybFl37Ek9+/LCAFYazsuVkiSqro5yT1Jrkzy7SS/mK2I3Zvk+iTfSXJHdz+/w+NYrACwIVZdrLBSiPaKEAFsjr1cNQcA+0aIABglRACMEiIARgkRAKOECIBRQgTAKCECYJQQATBKiAAYJUQAjBIiAEYJEQCjhAiAUUIEwCghAmDUoQN+vr/P1m9zTZI3LtubYJPmmmzWfDdprslmzXeT5prs/Xz/9ap3PNDf0PqKJ6462d1HR578gG3SXJPNmu8mzTXZrPlu0lyT2fk6NQfAKCECYNRkiO4efO6DtklzTTZrvps012Sz5rtJc00G5zt2jQgAEqfmABh24CGqqtuq6ptV9WRV3XnQz7/fquqzVXWmqh7ftu+aqnqwqr61fLx6cox7pareVFUPV9UTVfX1qvrosn9d5/vaqvqrqvrrZb6/sey/saoeWV7Tn6+qK6fHuleq6oqqeqyqvrRsr/Ncn6qqr1XVV6vq5LJvXV/LV1XVfVX1jao6VVXvnZzrgYaoqq5I8t+T/Pskb0/ys1X19oMcwwH4vSS3vWrfnUke6u6bkjy0bK+Dl5L8Sne/Pcl7kvzS8v1c1/l+L8mt3f2OJDcnua2q3pPkU0k+3d1vSfJCko8MjnGvfTTJqW3b6zzXJPl33X3ztmXM6/pavivJn3T325K8I1vf47m5dveB/Uny3iR/um3740k+fpBjOKB53pDk8W3b30xy7XL72iTfnB7jPs37/iTv34T5JvkXSb6S5N9k64cADy37X/Eav5z/JDmSrf8h3ZrkS0lqXee6zOepJG981b61ey0n+dEk/yvLGoFLYa4HfWruuiTf3bZ9etm37g539zPL7WeTHJ4czH6oqhuSvDPJI1nj+S6nqr6a5EySB5P8bZIXu/ul5S7r9Jr+7SS/muQHy/Ybsr5zTZJO8mdV9WhVHV/2reNr+cYkf5fkd5fTrvdU1esyOFeLFQ5Yb/1zY62WKlbVjyT5oyT/ubv/9/bPrdt8u/v73X1zto4W3p3kbcND2hdV9YEkZ7r70emxHKD3dfe7snXp4Jeq6t9u/+QavZYPJXlXkt/p7ncm+ae86jTcQc/1oEP0dJI3bds+suxbd89V1bVJsnw8MzyePVNVr8lWhH6/u7+w7F7b+b6su19M8nC2Tk9dVVUvv2/jurymfyLJz1TVU0k+l63Tc3dlPeeaJOnup5ePZ5J8MVv/0FjH1/LpJKe7+5Fl+75shWlsrgcdoi8nuWlZeXNlkg8neeCAxzDhgSTHltvHsnUt5bJXVZXkM0lOdfdvbfvUus73x6rqquX2D2fretipbAXpQ8vd1mK+3f3x7j7S3Tdk6+/pn3f3z2UN55okVfW6qnr9y7eT/FSSx7OGr+XufjbJd6vqrcuun0zyRAbneuA/0FpVP52tc89XJPlsd3/iQAewz6rqD5Pckq13sn0uya8n+eMk9ya5PlvvPn5Hdz8/Nca9UlXvS/KXSb6W/38d4deydZ1oHef740lOZOu1+0NJ7u3u36yqN2frqOGaJI8l+fnu/t7cSPdWVd2S5L909wfWda7LvL64bB5K8gfd/YmqekPW87V8c5J7klyZ5NtJfjHLazoDc/XOCgCMslgBgFFCBMAoIQJglBABMEqIABglRACMEiIARgkRAKP+L2Y93tJ+ejV5AAAAAElFTkSuQmCC\n",
+      "text/plain": [
+       "<Figure size 504x504 with 1 Axes>"
+      ]
+     },
+     "metadata": {},
+     "output_type": "display_data"
+    }
+   ],
+   "source": [
+    "x_bin = np.where(x > 1, -1, 1)\n",
+    "print(\"format of stack x_bin = \", x_bin.shape, \"\\n\"), \n",
+    "print(\"minimum pixel value after conversion = \", np.amin(x_bin), \"\\n\"), \n",
+    "print(\"maximum pixel value after conversion = \", np.amax(x_bin), \"\\n\")\n",
+    "# fig, ax = plt.subplots(figsize=(7, 7))\n",
+    "# ax.imshow(x_bin[30], cmap = 'gray')\n",
+    "fig, ax = plt.subplots(figsize=(7, 7))\n",
+    "ax.imshow(x_bin, cmap='gray')\n",
+    "None\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 6,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# save flat image stack as binary file\n",
+    "path = '/Users/jstark/Desktop/image_examples/'\n",
+    "filename = 'circle'\n",
+    "ips.save_array_toBin(x_bin, path, filename)\n",
+    "\n",
+    "# Save image size in #pixels / axis as csv file\n",
+    "dim = np.asarray(x_bin.shape)\n",
+    "np.savetxt('%s/size_%s.csv' % (path, filename), dim, fmt='%i', delimiter=',')"
+   ]
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Python 3",
+   "language": "python",
+   "name": "python3"
+  },
+  "language_info": {
+   "codemirror_mode": {
+    "name": "ipython",
+    "version": 3
+   },
+   "file_extension": ".py",
+   "mimetype": "text/x-python",
+   "name": "python",
+   "nbconvert_exporter": "python",
+   "pygments_lexer": "ipython3",
+   "version": "3.6.5"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 2
+}

example/Numerics/Sussman_redistancing/image_binary_conversion/image2binary_create_sphere.ipynb

+{
+ "cells": [
+  {
+   "cell_type": "code",
+   "execution_count": 1,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from tifffile import imread, imsave\n",
+    "import numpy as np\n",
+    "import os\n",
+    "import matplotlib.pyplot as plt\n",
+    "import cv2\n",
+    "import sys\n",
+    "sys.path.append(\"./lib\")\n",
+    "import morphological_operations as mo\n",
+    "import image_process_and_save as ips"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 2,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import numpy as np\n",
+    "from copy import deepcopy\n",
+    "\n",
+    "''' size : size of original 3D numpy matrix A.\n",
+    "    radius : radius of circle inside A which will be filled with ones. \n",
+    "'''\n",
+    "size, radius = 64, 5\n",
+    "\n",
+    "''' A : numpy.ndarray of shape size*size*size. '''\n",
+    "A = np.zeros((size,size, size)) \n",
+    "\n",
+    "''' AA : copy of A (you don't want the original copy of A to be overwritten.) '''\n",
+    "AA = deepcopy(A) \n",
+    "\n",
+    "''' (x0, y0, z0) : coordinates of center of circle inside A. '''\n",
+    "x0, y0, z0 = int(np.floor(A.shape[0]/2)), \\\n",
+    "        int(np.floor(A.shape[1]/2)), int(np.floor(A.shape[2]/2))\n",
+    "\n",
+    "\n",
+    "for x in range(AA.shape[0]):\n",
+    "    for y in range(AA.shape[1]):\n",
+    "        for z in range(AA.shape[2]):   \n",
+    "            d = radius*radius*radius - (x0-x)**2 - (y0-y)**2 - (z0-z)**2 \n",
+    "            if d <= 0: AA[x,y,z] = -1\n",
+    "            elif d >= 0: AA[x,y,z] = 1"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 3,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "(64, 64, 64)"
+      ]
+     },
+     "execution_count": 3,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "AA.shape"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 4,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAaIAAAGfCAYAAAAd79YcAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAEdVJREFUeJzt3VGs5ndd5/HP1w4NiqxtwZ00HbotoYFwIYVMWIhkU2swXZdIL0iD0WTWkJ0bTdjsbrR442pCIjdiLzYmTUHnQoWmgm24UJtaXa8qU4qhdGCpWMI0bcdN26x6gSl89+L8ux66nZ5n5jznfDvP83olJ+f5/8//nOf3m3na9/z/z+88T3V3AGDKD0wPAIDtJkQAjBIiAEYJEQCjhAiAUUIEwCghAmCUEAEwSogAGHVkP99cVbckuSPJZUnu6u7f3ON4L+MAsCW6u1Y5ri72JX6q6rIk/yvJ+5OcTfLFJD/b3Y+9wvcIEcCWWDVE+7k09+4kj3f3N7v7n5N8JskH9/HzANhC+wnRNUm+vWv77LLv+1TVyao6XVWn93FfAGyofT1HtIruvjPJnYlLcwD8//ZzRvRkkjft2j627AOAle0nRF9MckNVXV9Vlyf5cJL71jMsALbFRV+a6+4XquqXkvxpdpZvf7q7v7q2kQGwFS56+fZF3ZnniAC2xmEs3waAfRMiAEYJEQCjhAiAUUIEwCghAmCUEAEwSogAGCVEAIwSIgBGCREAo4QIgFFCBMAoIQJglBABMEqIABglRACMEiIARgkRAKOECIBRQgTAKCECYJQQATBKiAAYJUQAjBIiAEYJEQCjhAiAUUIEwCghAmCUEAEwSogAGCVEAIwSIgBGCREAo4QIgFFCBMAoIQJglBABMEqIABglRACMEiIARgkRAKOECIBRQgTAKCECYJQQATBKiAAYJUQAjBIiAEbtGaKq+nRVnauqR3ftu6qq7q+qbyyfrzzYYQKwqVY5I/q9JLe8ZN/tSR7o7huSPLBsA8AF2zNE3f0/kzz7kt0fTHJquX0qya1rHhcAW+JinyM62t1PLbefTnJ0TeMBYMsc2e8P6O6uqj7f16vqZJKT+70fADbTxZ4RPVNVVyfJ8vnc+Q7s7ju7+3h3H7/I+wJgg11siO5LcmK5fSLJvesZDgDbprrPe1Vt54CqP0xyU5I3Jnkmya8l+eMkdye5Nsm3ktzW3S9d0PByP+uV7wyAjdHdtcpxe4ZonYQIYHusGiKvrADAKCECYJQQATBKiAAYJUQAjBIiAEYJEQCjhAiAUUIEwCghAmCUEAEwSogAGCVEAIwSIgBGCREAo4QIgFFCBMAoIQJglBABMEqIABglRACMEiIARgkRAKOECIBRQgTAKCECYJQQATBKiAAYJUQAjBIiAEYJEQCjhAiAUUIEwCghAmCUEAEwSogAGCVEAIwSIgBGCREAo4QIgFFCBMAoIQJglBABMEqIABglRACMEiIARgkRAKOECIBRQgTAKCECYJQQATBqzxBV1Zuq6sGqeqyqvlpVH132X1VV91fVN5bPVx78cAHYNNXdr3xA1dVJru7uL1XV65M8nOTWJP8xybPd/ZtVdXuSK7v7V/b4Wa98ZwBsjO6uVY7b84you5/q7i8tt/8hyZkk1yT5YJJTy2GnshMnALggRy7k4Kq6Lsk7kzyU5Gh3P7V86ekkR8/zPSeTnLz4IQKwyfa8NPf/Dqz64SR/meTj3f25qnq+u6/Y9fXnuvsVnydyaQ5ge6zt0lySVNVrkvxRkt/v7s8tu59Znj968XmkcxczUAC22yqr5irJp5Kc6e7f2vWl+5KcWG6fSHLv+ocHl5bu3vcHbJtVVs29L8lfJflKku8tu381O88T3Z3k2iTfSnJbdz+7x8/yXxkbbR0h2fm3H1z6Vr00t/JzROsgRGw6IYJ/sdbniADgoAgRAKOECIBRF/QLrbCpXk2r1dY1Fs81calwRgTAKCECYJQQATBKiAAYJUQAjLJqjo30aloFN+VC/wyssmOKMyIARgkRAKOECIBRQgTAKCECYJRVc1zSrI5bn/P9WVpNx0FzRgTAKCECYJQQATBKiAAYJUQAjLJqjkuC1XFzrKbjoDkjAmCUEAEwSogAGCVEAIwSIgBGWTXHq4rVcZcOq+lYF2dEAIwSIgBGCREAo4QIgFFCBMAoIQJglBABMEqIABglRACMEiIARgkRAKO81hwjvKbc5vIadFwoZ0QAjBIiAEYJEQCjhAiAUUIEwCghAmCUEAEwSogAGCVEAIwSIgBG7RmiqnptVf11Vf1NVX21qn592X99VT1UVY9X1Wer6vKDHy4Am2aVM6LvJLm5u9+R5MYkt1TVe5J8Isknu/stSZ5L8pGDGyYAm2rPEPWOf1w2X7N8dJKbk9yz7D+V5NYDGSEAG22l54iq6rKq+nKSc0nuT/K3SZ7v7heWQ84mueZghgjAJlspRN393e6+McmxJO9O8rZV76CqTlbV6ao6fZFjBGCDXdCque5+PsmDSd6b5IqqevH9jI4lefI833Nndx/v7uP7GikAG2mVVXM/WlVXLLd/MMn7k5zJTpA+tBx2Ism9BzVIADZX7fVOmVX1Y9lZjHBZdsJ1d3f/RlW9OclnklyV5JEkP9/d39njZ3lbzi3jnVjZi3du3VzdvdJf7p4hWich2j5CxF6EaHOtGiKvrADAKCECYJQQATBKiAAYJUQAjBIiAEYJEQCjhAiAUUIEwCghAmCUEAEwSogAGCVEAIwSIgBGCREAo47sfQhcvPO914z3Kdo+3neI83FGBMAoIQJglBABMEqIABglRACMEiIARgkRAKOECIBRQgTAKCECYJQQATBKiAAYJUQAjBIiAEYJEQCjhAiAUUIEwCjv0MoI79y6ubwTKxfKGREAo4QIgFFCBMAoIQJglBABMEqIABglRACMEiIARgkRAKOECIBRQgTAKK81x6uK16C7dHhNOdbFGREAo4QIgFFCBMAoIQJg1MohqqrLquqRqvrCsn19VT1UVY9X1Wer6vKDGyYAm+pCzog+muTMru1PJPlkd78lyXNJPrLOgcFuVfWyHxw8f/YctJVCVFXHkvyHJHct25Xk5iT3LIecSnLrQQwQgM226hnRbyf55STfW7bfkOT57n5h2T6b5Jo1jw2ALbBniKrqA0nOdffDF3MHVXWyqk5X1emL+X4ANtsqr6zw40l+pqp+Oslrk/yrJHckuaKqjixnRceSPPly39zddya5M0mqyq/HA/B99jwj6u6Pdfex7r4uyYeT/Hl3/1ySB5N8aDnsRJJ7D2yUAGys/fwe0a8k+S9V9Xh2njP61HqGBKuzomt9/FkypQ7zxSRdmuOweJHUCyc6rFt3r/Sg8soKAIwSIgBGCREAo4QIgFHeoZWNdKFPvG/i4gaLD7hUOCMCYJQQATBKiAAYJUQAjBIiAEZZNQdZ3wqzday+s9qNbeOMCIBRQgTAKCECYJQQATBKiAAYZdUcrJEVb3DhnBEBMEqIABglRACMEiIARgkRAKOECIBRQgTAKCECYJQQATBKiAAYJUQAjBIiAEYJEQCjhAiAUUIEwCghAmCUEAEwSogAGCVEAIwSIgBGCREAo4QIgFFCBMAoIQJglBABMEqIABglRACMEiIARgkRAKOECIBRQgTAKCECYJQQATBKiAAYdWSVg6rqiST/kOS7SV7o7uNVdVWSzya5LskTSW7r7ucOZpgAbKoLOSP6ie6+sbuPL9u3J3mgu29I8sCyDQAXZD+X5j6Y5NRy+1SSW/c/HAC2zaoh6iR/VlUPV9XJZd/R7n5quf10kqNrHx0AG2+l54iSvK+7n6yqf53k/qr62u4vdndXVb/cNy7hOvlyXwOA6n7Zfpz/G6r+e5J/TPKfktzU3U9V1dVJ/qK737rH917YnQFwyeruWuW4PS/NVdXrqur1L95O8lNJHk1yX5ITy2Enktx7cUMFYJvteUZUVW9O8vll80iSP+juj1fVG5LcneTaJN/KzvLtZ/f4Wc6IALbEqmdEF3xpbj+ECGB7rO3SHAAcJCECYJQQATBKiAAYJUQAjBIiAEYJEQCjhAiAUUIEwCghAmCUEAEwSogAGCVEAIwSIgBGCREAo4QIgFFCBMAoIQJglBABMEqIABglRACMEiIARgkRAKOECIBRQgTAKCECYJQQATBKiAAYJUQAjBIiAEYJEQCjhAiAUUIEwCghAmCUEAEwSogAGCVEAIwSIgBGCREAo4QIgFFCBMAoIQJglBABMEqIABglRACMEiIARgkRAKOECIBRQgTAKCECYNRKIaqqK6rqnqr6WlWdqar3VtVVVXV/VX1j+XzlQQ8WgM2z6hnRHUn+pLvfluQdSc4kuT3JA919Q5IHlm0AuCDV3a98QNWPJPlykjf3roOr6utJburup6rq6iR/0d1v3eNnvfKdAbAxurtWOW6VM6Lrk/x9kt+tqkeq6q6qel2So9391HLM00mOXtxQAdhmq4ToSJJ3Jfmd7n5nkn/KSy7DLWdKL3u2U1Unq+p0VZ3e72AB2DyrhOhskrPd/dCyfU92wvTMckkuy+dzL/fN3X1ndx/v7uPrGDAAm2XPEHX300m+XVUvPv/zk0keS3JfkhPLvhNJ7j2QEQKw0fZcrJAkVXVjkruSXJ7km0l+ITsRuzvJtUm+leS27n52j59jsQLAllh1scJKIVoXIQLYHutcNQcAB0aIABglRACMEiIARgkRAKOECIBRQgTAKCECYJQQATBKiAAYJUQAjBIiAEYJEQCjhAiAUUIEwCghAmDUkUO+v/+dnXdzTZI3LtvbYJvmmmzXfLdprsl2zXeb5pqsf77/ZtUDD/UdWr/vjqtOd/fxkTs/ZNs012S75rtNc022a77bNNdkdr4uzQEwSogAGDUZojsH7/uwbdNck+2a7zbNNdmu+W7TXJPB+Y49RwQAiUtzAAw79BBV1S1V9fWqeryqbj/s+z9oVfXpqjpXVY/u2ndVVd1fVd9YPl85OcZ1qao3VdWDVfVYVX21qj667N/U+b62qv66qv5mme+vL/uvr6qHlsf0Z6vq8umxrktVXVZVj1TVF5btTZ7rE1X1lar6clWdXvZt6mP5iqq6p6q+VlVnquq9k3M91BBV1WVJ/keSf5/k7Ul+tqrefphjOAS/l+SWl+y7PckD3X1DkgeW7U3wQpL/2t1vT/KeJL+4/H1u6ny/k+Tm7n5HkhuT3FJV70nyiSSf7O63JHkuyUcGx7huH01yZtf2Js81SX6iu2/ctYx5Ux/LdyT5k+5+W5J3ZOfveG6u3X1oH0nem+RPd21/LMnHDnMMhzTP65I8umv760muXm5fneTr02M8oHnfm+T92zDfJD+U5EtJ/m12fgnwyLL/+x7jl/JHkmPZ+R/SzUm+kKQ2da7LfJ5I8saX7Nu4x3KSH0nyd1nWCLwa5nrYl+auSfLtXdtnl32b7mh3P7XcfjrJ0cnBHISqui7JO5M8lA2e73Kp6stJziW5P8nfJnm+u19YDtmkx/RvJ/nlJN9btt+QzZ1rknSSP6uqh6vq5LJvEx/L1yf5+yS/u1x2vauqXpfBuVqscMh6558bG7VUsap+OMkfJfnP3f1/dn9t0+bb3d/t7huzc7bw7iRvGx7SgaiqDyQ5190PT4/lEL2vu9+VnacOfrGq/t3uL27QY/lIkncl+Z3ufmeSf8pLLsMd9lwPO0RPJnnTru1jy75N90xVXZ0ky+dzw+NZm6p6TXYi9Pvd/bll98bO90Xd/XySB7NzeeqKqnrxdRs35TH940l+pqqeSPKZ7FyeuyObOdckSXc/uXw+l+Tz2fmHxiY+ls8mOdvdDy3b92QnTGNzPewQfTHJDcvKm8uTfDjJfYc8hgn3JTmx3D6RnedSLnlVVUk+leRMd//Wri9t6nx/tKquWG7/YHaeDzuTnSB9aDlsI+bb3R/r7mPdfV12/jv98+7+uWzgXJOkql5XVa9/8XaSn0ryaDbwsdzdTyf5dlW9ddn1k0key+BcD/0XWqvqp7Nz7fmyJJ/u7o8f6gAOWFX9YZKbsvNKts8k+bUkf5zk7iTXZufVx2/r7menxrguVfW+JH+V5Cv5l+cRfjU7zxNt4nx/LMmp7Dx2fyDJ3d39G1X15uycNVyV5JEkP9/d35kb6XpV1U1J/lt3f2BT57rM6/PL5pEkf9DdH6+qN2QzH8s3JrkryeVJvpnkF7I8pjMwV6+sAMAoixUAGCVEAIwSIgBGCREAo4QIgFFCBMAoIQJglBABMOr/AvCduw1wUNEzAAAAAElFTkSuQmCC\n",
+      "text/plain": [
+       "<Figure size 504x504 with 1 Axes>"
+      ]
+     },
+     "metadata": {},
+     "output_type": "display_data"
+    }
+   ],
+   "source": [
+    "fig, ax = plt.subplots(figsize=(7, 7))\n",
+    "ax.imshow(AA[int(AA.shape[0]/2)], cmap='gray')\n",
+    "None"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 5,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# flatten stack to get a 1D array\n",
+    "x_bin = AA.flatten();"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 6,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "(262144,)"
+      ]
+     },
+     "execution_count": 6,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "x_bin.shape"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 7,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# save flat image stack as binary file\n",
+    "path = '/Users/jstark/Desktop/image_examples/'\n",
+    "filename = 'sphere'\n",
+    "ips.save_array_toBin(x_bin, path, filename)\n",
+    "\n",
+    "# Save image size in #pixels / axis as csv file\n",
+    "dim = np.asarray(AA.shape)\n",
+    "np.savetxt('%s/size_%s.csv' % (path, filename), dim, fmt='%i', delimiter=',')"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": []
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Python 3",
+   "language": "python",
+   "name": "python3"
+  },
+  "language_info": {
+   "codemirror_mode": {
+    "name": "ipython",
+    "version": 3
+   },
+   "file_extension": ".py",
+   "mimetype": "text/x-python",
+   "name": "python",
+   "nbconvert_exporter": "python",
+   "pygments_lexer": "ipython3",
+   "version": "3.6.5"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 2
+}