Classes | |
class | ChanVeseForce |
This ForceGenerator creates a normal force field, that leads to a segmentation similar to the model proposed by Chan and Vese. More... | |
class | UpdatableImage |
Display the pixel distribution of the image as a graph (histogram).Graphical display of an Image, with basic user interaction functionality. More... | |
class | ConstantForce |
This ForceGenerator creates a constant normal force, which does not depend on the underlying image or its current segmentation. More... | |
class | ForceGenerator |
This is the abstract base class for the generators that define the normal force driving the evolution of a Region's LevelSetFunction. More... | |
class | Histogram |
Class representing a histogram. Its main use is in the ParzenDistributionForce class. More... | |
class | Image |
An 1, 2 or 3-dimensional image with an arbitrary number of channels. More... | |
class | LevelSetFunction |
A class representing a level-set function in 1, 2 or 3D. More... | |
class | MRBalloonForce |
This ForceGenerator creates a normal force field, which does not depend on the underlying image, only on its current segmentation. More... | |
class | NormalDistributionForce |
This ForceGenerator creates a normal force field derived from the minimization of the Mumford-Shah functional, in a setting where each region is modeled statistically using a normal distribution. More... | |
class | ParzenDistributionForce |
This ForceGenerator creates a normal force field derived from the minimization of the Mumford-Shah functional, in a setting where the pixel distribution of each region is modeled by a parzen estimate. More... | |
struct | Region |
This struct represents the information completely describing one particular segmented region of an image. It can be used alone in a two-region segmentation setting (the other region being its complement), or several Regions can be used in a multiregion segmentation setting. More... | |
Typedefs | |
typedef Image< char > | Mask |
The Mask is an Image used for defining active and inactive pixels of another, identically shaped image. | |
Enumerations | |
enum | BorderCondition { DIRICHLET = 0, NEUMANN = 1 } |
An enumeration of possible border conditions when solving a PDE. More... | |
enum | SEG_REGION { SEG_NEGATIVE = 0, SEG_POSITIVE } |
Functions | |
void | load_image (const char *name, Image< double > &target, bool convert_to_greyscale=false) |
load an Image<double> from file | |
void | load_image (const char *name, Image< int > &target, bool convert_to_greyscale=false) |
load an Image<int> from file | |
void | save_image (const char *name, const Image< double > &target) |
save an Image<double> to file | |
void | display_image (const Image< double > &img, int z=0) |
Display an image in a graphical window and pause the program until the window is closed. | |
void | permanent_display (const Image< double > &img) |
Display an Image in a graphical window an continue execution of program. | |
void | blur_image (Image< double > &img, double rho) |
Gaussian blur an Image<double> with a certain strength. | |
void | blur_1D (double *data, unsigned int data_size, double rho) |
Gaussian blur an 1D-array of double values with a certain strength. | |
void | gaussian_noise (Image< double > &img, double sigma=0) |
fill an Image<double> with Gaussian noise. | |
void | compute_squared_gradient_sum_2D (const Image< double > &image, Image< double > &result) |
Compute an estimate of the squared gradient norm, for each pixel in a two-dimensional image. | |
void | compute_squared_gradient_sum_3D (const Image< double > &image, Image< double > &result) |
Compute an estimate of the squared gradient norm, for each pixel in a three-dimensional image. | |
void | analytic_eigvals (double alpha1, double alpha2, double alpha3, double beta, double gamma, double delta, double &lambda1, double &lambda2, double &lambda3) |
Compute the eigenvalues of a 3x3 symmetric matrix by analytically solving the corresponding characteristic polynomial. | |
void | analytic_eigsys (double alpha1, double alpha2, double alpha3, double beta, double gamma, double delta, double &lambda1, double &lambda2, double &lambda3, double *v1, double *v2, double *v3) |
Compute the eigenvalues and eigenvectors of a 3x3 symmetric matrix by analytically solving the corresponding characteristic polynomial to find the eigenvalues, and then explicitly constructing the corresponding eigenvectors. | |
void | numeric_eigsys (double alpha1, double alpha2, double alpha3, double beta, double gamma, double delta, double &lambda1, double &lambda2, double &lambda3, double *v1, double *v2, double *v3) |
Compute the eigenvalues and eigenvectors of a 3x3 symmetric matrix by means of an iterative numerical algorithm using Householder transformations and Givens rotations. | |
void | compute_structure_tensor_2D (const Image< double > &img, Image< double > &G, bool square_root=false) |
Compute the structure tensor of a 2D-image, and writes the result to a 3-channeled 2D-image of the same resolution. | |
void | compute_structure_tensor_3D (const Image< double > &img, Image< double > &G, bool square_root=false) |
Compute the structure tensor of a 3D-image, and writes the result to a 6-channeled 2D-image of the same resolution. | |
void | compute_scale_factor_2D (const Image< double > &img, Image< double > &scale_accum, Image< double > &time_accum, double dt, double T) |
void | compute_scale_factor_3D (const Image< double > &img, Image< double > &scale_accum, Image< double > &time_accum, double dt, double T) |
void | nonlinear_gauss_filter_2D (const Image< double > &img_old, Image< double > &img_new, double dt, double sigma, double p) |
void | nonlinear_gauss_filter_3D (const Image< double > &img_old, Image< double > &img_new, double dt, double sigma, double p) |
void | anisotropic_smoothing (const Image< double > &img_old, Image< double > &img_new, double dt, double sigma, double p) |
void | compute_smoothing_geometry_3D (const Image< double > &G, double p1, double p2, double p3, Image< double > &T, bool take_square_root) |
void | compute_smoothing_geometry_2D (const Image< double > &G, double p1, double p2, Image< double > &T, bool take_square_root) |
void | normal_direction_flow_2D (const LevelSetFunction &phi, LevelSetFunction &advect, BorderCondition bcond, double &H1, double &H2, const Mask *const changes_allowed_mask, const Mask *const defined_region_mask) |
void | normal_direction_flow_3D (const LevelSetFunction &phi, LevelSetFunction &advect, BorderCondition bcond, double &H1, double &H2, double &H3, const Mask *const changes_allowed_mask, const Mask *const defined_region_mask) |
void | visualize_multisets (const LevelSetFunction **const images, int num_images, Image< double > &target, const int **const rgb_color) |
int | visualize_level_set (const LevelSetFunction &img, Image< double > &target, double threshold, double outside_intensity, double curve_intensity, double interior_intensity, const Mask *const mask, double mask_intensity) |
void | LIC_2D_FS (const Image< double > &src, const Image< double > &vec, Image< double > &target, double(*kernel_func)(double), const double dl, const double L) |
void | LIC_3D_FS (const Image< double > &src, const Image< double > &vec, Image< double > &target, double(*kernel_func)(double), const double dl, const double L) |
void | LIC_3D (const Image< double > &src, const Image< double > &vec, Image< double > &target, double(*kernel_func)(double), const double L) |
void | LIC_2D (const Image< double > &src, const Image< double > &vec, Image< double > &target, double(*kernel_func)(double), const double L) |
double | develop_multiregion_2D (Region *regs, int num_regs, int num_iter, int reinit_modulo, const Mask *geom_mask) |
double | develop_multiregion_3D (Region *regs, int num_regs, int num_iter, int reinit_modulo, const Mask *geom_mask) |
void | negate (Image< double > &img, double min=0, double max=255) |
make an Image<double> a negative of itself. | |
void | rescale (Image< double > &img, double cur_min, double cur_max, double to_min=0, double to_max=255) |
rescale the pixel range of an Image<double> so that its pixels spans a different range than before rescaling | |
void | rescale_channels (Image< double > &img, double to_min=0, double to_max=255) |
rescale the pixel ranges of each of an image's channels to a user-specified interval | |
void | clip (Image< double > &img, double min, double max) |
clip the pixel range of a image to a certain interval | |
void | transpose (Image< double > &img) |
transpose the x and y dimensions of an Image<double> | |
void | transpose (const Image< double > &img, Image< double > &target) |
Generate an Image<double> that is the transpose of another Image<double>. | |
void | horizontal_sinusoidal_bands (LevelSetFunction &img, int num_bands, double phase=0) |
Initializes a LevelSetFunction with horizontal bands created by a sine function. | |
void | rectangle (LevelSetFunction &img, double xmin_ratio, double xmax_ratio, double ymin_ratio, double ymax_ratio, double zmin_ratio=0, double zmax_ratio=1) |
Initializes a LevelSetFunction to become 1 everywhere, except for an interior rectangular domain where the LevelSetFunction has the value -1. | |
void | init_voronoi_regions (LevelSetFunction *regs, const double *center_coords,int num_regions, bool three_d=false) |
Initializes a set of LevelSetFunctions defined on a common domain so that their closed regions together constitute a voronoi subdivision of the domain. The points defining the voronoi subdivision is given by the user. | |
void | multiregion_bands (LevelSetFunction *regs, int num_regs, int pixel_bandwidth) |
Initializes a set of LevelSetFunctions defined on a common domain so that their closed regions together constitute a subdivision on the domain, and where each such region consists of a set of horizontal bands of a certain width along the y-direction. | |
void | random_scattered_voronoi (LevelSetFunction *regs, int num_regs, int num_fragments) |
Initializes a set of LevelSetFunctions defined on a common domain so that their closed regions together constitute a voronoi subdivision of the domain. The points defining the voronoi subdivision are randomly generated. | |
void | set_from_parzen (LevelSetFunction &img, const ParzenDistributionForce pf, const Mask *m=0) |
Initializes a LevelSetFunction to become a function with two values, 1 and -1, corresponding to the domains in which the force defined by the pf argument is positive or negative. | |
void | sphere (LevelSetFunction &img, double relrad,double xrelpos,double yrelpos,double zrelpos=0) |
Initializes a LevelSetFunction to become the signed distance function from a point in its domain minus some fixed value, thus its zero-set becomes a circle around . | |
void | make_border_mask (const LevelSetFunction &phi, Mask &target, int width=2, const Mask *geom_mask=0) |
Generate a Mask over a domain, where the active region of the Mask is specified as the region of a certain width (in pixels) around the zero-set of a LevelSetFunction defined on the same domain. | |
void | mask_from_segmentation (const LevelSetFunction &phi, Mask &target, SEG_REGION reg) |
Generate a Mask over a domain, where the active region of the Mask is specified as the region defined by the negative (or positive) region of a LevelSetFunction over the same domain. | |
void | read_image_sequence (istream &image_list, Image< double > &result, bool convert_to_grayscale) |
unsigned long int | nonzeroes (const Image< int > &img) |
Count the number of nonzero pixels of an Image<int>. | |
double | nonzero_ratio (const Image< int > &img) |
Calculate the ratio of the number pixels with a nonzero value to the total number of pixels in an Image<int>. | |
unsigned long int | positives (const Image< double > &img) |
Count the number of pixels with a nonnegative value in an Image<int>. | |
unsigned long int | negatives (const Image< double > &img) |
Count the number of pixels with a negative value in an Image<int>. | |
double | positive_ratio (const Image< double > &img) |
Calculate the ratio of the number of pixels with a nonnegative value to the total number of pixels in an Image<double>. | |
double | negative_ratio (const Image< double > &img) |
Calculate the ratio of the number of pixels with a negative value to the total number of pixels in an Image<double>. | |
double | develop_single_region_2D (Region ®, int num_iter, int reinit_modulo, const Mask *geom_mask) |
double | develop_single_region_3D (Region ®, int num_iter, int reinit_modulo, const Mask *geom_mask) |
template<typename T> | |
void | combine_channel_images (const Image< T > **channel_array, int num_input_images, Image< T > &result) |
template<typename T> | |
void | to_grayscale (Image< T > &img) |
convert a three-channel Image (RGB of double , int , etc.) to greyscale. | |
void | read_image_sequence (std::istream &image_list, Image< double > &result, bool convert_to_grayscale=false) |
Read a set of 2D-image files and write the result into a 3D Image<double>, where the read images will be stacked along the z-coordinate direction. | |
template<typename ImgType> | |
void | resample_into (const ImgType &src, ImgType &target, bool linear=true) |
Resample an image, using nearest-pixel or linear interpolation. | |
template<typename ImgType> | |
void | downsample_series (const ImgType &input, std::vector< ImgType > &result, int min_num_pixels, bool downscale_z=false, double downscale_factor=2, bool to_grayscale=false, bool linear=true) |
Create a sequence of progressively smaller, resampled copies of a reference image. | |
Variables | |
const int | RED [] = {255, 0, 0} |
3-tuple of int representing the color red in RGB format. | |
const int | BLUE [] = {0, 255,0} |
3-tuple of int representing the color blue in RGB format. | |
const int | GREEN [] = {0, 0, 255} |
3-tuple of int representing the color green in RGB format. | |
const int | YELLOW [] = {255,0, 255} |
3-tuple of int representing the color yellow in RGB format. | |
const int | CYAN [] = {0, 255, 255} |
3-tuple of int representing the color cyan in RGB format. | |
const int | MAGENTA [] = {255, 255, 0} |
3-tuple of int representing the color magenta in RGB format. | |
const int | WHITE [] = {255,255, 255} |
3-tuple of int representing the color white in RGB format. | |
const int | BLACK [] = {0, 0, 0} |
3-tuple of int representing the color black in RGB format. | |
const int | GREY [] = {200, 200, 200} |
3-tuple of int representing the color grey in RGB format. | |
const int | BROWN [] = {200, 100, 40} |
3-tuple of int representing the color brown in RGB format. |
The name lsseg
stands for "level-set segmentation", which was the initial aim of this code project. However, in addition to level-set segmentation code, based on [Brox05], this namespace also contains functionality for anisotropic smoothing of images, based on an algorithm presented by David Tschumperlé in [Tschumperle06], and with background material in [Tschumperle02].
typedef Image<char> lsseg::Mask |
The Mask is an Image used for defining active and inactive pixels of another, identically shaped image.
As such, the pixel values of a Mask as interpreted as boolean, i.e. as either true
or false
. (In the current implementation, these boolean values are stored as char
; the reason for this being that this was most efficient computationally speaking on the platform where this library was developed). The Mask is typically single-channelled, and used in conjunction with another Image of the same x
, y
and z
spatial resolution, so that there can be a one-to-one correspondence with the pixels in each channel of the Image and those of the Mask. Those image pixels who are associated with pixels in the Mask having false
as their value, are considered inactive ("masked out").
An enumeration of possible border conditions when solving a PDE.
Definition at line 65 of file BorderCondition.h.
enum lsseg::SEG_REGION |
void lsseg::load_image | ( | const char * | name, | |
Image< double > & | target, | |||
bool | convert_to_greyscale = false | |||
) |
load an Image<double> from file
name | - name of the file containing the image. The file extension, being a part of the name, specifies what format the image is stored in (png, jpg, etc.) | |
convert_to_greyscale | - if set to true , the loaded image will be converted to greyscale, and thus be guaranteed to contain exactly one image channel. If false , the number of image channels will be those found in the file. |
target | - Upon function completion, target will contain the image loaded from file |
Definition at line 104 of file cimg_dependent.C.
void lsseg::load_image | ( | const char * | name, | |
Image< int > & | target, | |||
bool | convert_to_greyscale = false | |||
) |
load an Image<int> from file
name | - name of the file containing the image. The file extension, being a part of the name, specifies what format the image is stored in (png, jpg, etc.) | |
convert_to_greyscale | - if set to true , the loaded image will be converted to greyscale, and thus be guaranteed to contain exactly one image channel. If false , the number of image channels will be those found in the file. |
target | - Upon function completion, target will contain the image loaded from file |
Definition at line 115 of file cimg_dependent.C.
void lsseg::save_image | ( | const char * | name, | |
const Image< double > & | target | |||
) |
save an Image<double> to file
name | the filename where the image should be stored. The storage format is specified by the filename extension (jpg, png, gif, etc..) | |
target | the image to be written to file. |
Definition at line 126 of file cimg_dependent.C.
void lsseg::display_image | ( | const Image< double > & | img, | |
int | z = 0 | |||
) |
Display an image in a graphical window and pause the program until the window is closed.
img | the image to be displayed | |
z | in case of a 3D image, z specifies which z-coordinate slice to display on screen (the image shown will always be in the (x, y)-plane). |
Definition at line 135 of file cimg_dependent.C.
void lsseg::permanent_display | ( | const Image< double > & | img | ) |
Display an Image in a graphical window an continue execution of program.
img | the image to display |
Definition at line 149 of file cimg_dependent.C.
void lsseg::blur_image | ( | Image< double > & | img, | |
double | rho | |||
) |
Gaussian blur an Image<double> with a certain strength.
img | the image to be blurred | |
rho | the parameter of the Gaussian kernel. Should be positive. The higher the value of rho , the more the picture will be blurred. |
Definition at line 162 of file cimg_dependent.C.
void lsseg::blur_1D | ( | double * | data, | |
unsigned int | data_size, | |||
double | rho | |||
) |
Gaussian blur an 1D-array of double
values with a certain strength.
data | pointer to the memory array to be blurred. | |
data_size | size of the memory array to be blurred | |
rho | the parameter of the Gaussian kernel. Should be positive. The higher the value of rho , the more the picture will be blurred. |
Definition at line 172 of file cimg_dependent.C.
void lsseg::gaussian_noise | ( | Image< double > & | img, | |
double | sigma = 0 | |||
) |
fill an Image<double> with Gaussian noise.
img | the image to be filled with Gaussian noise | |
sigma | the parameter of the Gaussian kernel. Should be nonnegative. |
Definition at line 182 of file cimg_dependent.C.
void lsseg::compute_squared_gradient_sum_2D | ( | const Image< double > & | image, | |
Image< double > & | result | |||
) |
Compute an estimate of the squared gradient norm, for each pixel in a two-dimensional image.
If the image contains more than one channel, the sum of the squared gradient norms of each channel will be computed for each image pixel.
[in] | image | The input image. It is supposed to be 2D (only 1 pixel along z -direction). It can have an arbitrary number of channels. |
[out] | result | Upon completion of the function, this image will have been resized to the same resolution as image , but will contain only one channel. Each pixel in this channel corresponds to the estimated squared gradient norm of the corresponding pixel in image (or, in the case that image is multichanneled, it will correspond to the sum of the squared gradient norms of this pixel for each image channel in image . |
Definition at line 55 of file DiscreteApproximations.C.
void lsseg::compute_squared_gradient_sum_3D | ( | const Image< double > & | image, | |
Image< double > & | result | |||
) |
Compute an estimate of the squared gradient norm, for each pixel in a three-dimensional image.
If the image contains more than one channel, the sum of the squared gradient norms of each channel will be computed for each image pixel.
[in] | image | The input image. It is supposed to be 3D (more than 1 pixel along z -direction). image can in theory also be a 2D image, but in that case it is better to call the slightly faster function compute_squared_gradient_sum_2D(). image can have an arbitrary number of channels. |
[out] | result | Upon completion of the function, this image will have been resized to the same resolution as image , but will contain only one channel. Each pixel in this channel corresponds to the estimated squared gradient norm of the corresponding pixel in image (or, in the case that image is multichanneled, it will correspond to the sum of the squared gradient norms of this pixel for each image channel in image . |
Definition at line 91 of file DiscreteApproximations.C.
void lsseg::analytic_eigvals | ( | double | alpha1, | |
double | alpha2, | |||
double | alpha3, | |||
double | beta, | |||
double | gamma, | |||
double | delta, | |||
double & | lambda1, | |||
double & | lambda2, | |||
double & | lambda3 | |||
) |
Compute the eigenvalues of a 3x3 symmetric matrix by analytically solving the corresponding characteristic polynomial.
The terms in the symmetric matrix are given like this:
[in] | alpha1 | in the matrix above (first diagonal term) |
[in] | alpha2 | in the matrix above (second diagonal term) |
[in] | alpha3 | in the matrix above (third diagonal term) |
[in] | beta | in the matrix above |
[in] | gamma | in the matrix above |
[in] | delta | in the matrix above |
[out] | lambda1 | the first of the computed eigenvalues (in arbitrary order) |
[out] | lambda2 | the second of the computed eigenvalues (in arbitrary order) |
[out] | lambda3 | the third of the computed eigenvalues (in arbitrary order) |
Definition at line 177 of file EigValComp3x3.C.
void lsseg::analytic_eigsys | ( | double | alpha1, | |
double | alpha2, | |||
double | alpha3, | |||
double | beta, | |||
double | gamma, | |||
double | delta, | |||
double & | lambda1, | |||
double & | lambda2, | |||
double & | lambda3, | |||
double * | v1, | |||
double * | v2, | |||
double * | v3 | |||
) |
Compute the eigenvalues and eigenvectors of a 3x3 symmetric matrix by analytically solving the corresponding characteristic polynomial to find the eigenvalues, and then explicitly constructing the corresponding eigenvectors.
The terms in the symmetric matrix are given like this:
[in] | alpha1 | in the matrix above (first diagonal term) |
[in] | alpha2 | in the matrix above (second diagonal term) |
[in] | alpha3 | in the matrix above (third diagonal term) |
[in] | beta | in the matrix above |
[in] | gamma | in the matrix above |
[in] | delta | in the matrix above |
[out] | lambda1 | the first of the computed eigenvalues (the largest one) |
[out] | lambda2 | the second of the computed eigenvalues (the middle one) |
[out] | lambda3 | the third of the computed eigenvalues (the smallest one) |
[out] | v1 | pointer to an array where the eigenvector corresponding to lambda1 has been stored. (The array contains 3 elements). |
[out] | v2 | pointer to an array where the eigenvector corresponding to lambda2 has been stored. (The array contains 3 elements). |
[out] | v3 | pointer to an array where the eigenvector corresponding to lambda3 has been stored. (The array contains 3 elements). |
Definition at line 232 of file EigValComp3x3.C.
void lsseg::numeric_eigsys | ( | double | alpha1, | |
double | alpha2, | |||
double | alpha3, | |||
double | beta, | |||
double | gamma, | |||
double | delta, | |||
double & | lambda1, | |||
double & | lambda2, | |||
double & | lambda3, | |||
double * | v1, | |||
double * | v2, | |||
double * | v3 | |||
) |
Compute the eigenvalues and eigenvectors of a 3x3 symmetric matrix by means of an iterative numerical algorithm using Householder transformations and Givens rotations.
The terms in the symmetric matrix are given like this:
[in] | alpha1 | in the matrix above (first diagonal term) |
[in] | alpha2 | in the matrix above (second diagonal term) |
[in] | alpha3 | in the matrix above (third diagonal term) |
[in] | beta | in the matrix above |
[in] | gamma | in the matrix above |
[in] | delta | in the matrix above |
[out] | lambda1 | the first of the computed eigenvalues (arbitrary order) |
[out] | lambda2 | the second of the computed eigenvalues (arbitrary order) |
[out] | lambda3 | the third of the computed eigenvalues (arbitrary order) |
[out] | v1 | pointer to an array where the eigenvector corresponding to lambda1 has been stored. (The array contains 3 elements). |
[out] | v2 | pointer to an array where the eigenvector corresponding to lambda2 has been stored. (The array contains 3 elements). |
[out] | v3 | pointer to an array where the eigenvector corresponding to lambda3 has been stored. (The array contains 3 elements). |
Definition at line 245 of file EigValComp3x3.C.
void lsseg::compute_structure_tensor_2D | ( | const Image< double > & | img, | |
Image< double > & | G, | |||
bool | square_root = false | |||
) |
Compute the structure tensor of a 2D-image, and writes the result to a 3-channeled 2D-image of the same resolution.
[in] | img | the image we want to compute the structure tensor for. It should have two dimensions (x and y ) and an arbitrary number of channels. |
[out] | G | this image will contain the three components of the structure tensor, each stored in a separate channel. The first channel will be the square of the estimated partial derivative in x , the second channel will be the estimated partial derivative in x multiplied by the estimated partial derivative in y , and the third channel will be the square of the estimated partial derivative in y . |
[in] | square_root | if this is set to 'true', then the square root of the tensor will be calculated in each pixel, rather than the tensor itself. |
void lsseg::compute_structure_tensor_3D | ( | const Image< double > & | img, | |
Image< double > & | G, | |||
bool | square_root = false | |||
) |
Compute the structure tensor of a 3D-image, and writes the result to a 6-channeled 2D-image of the same resolution.
[in] | img | the image we want to compute the structure tensor for. It should have three dimensions (x , y and z ) and an arbitrary number of channels. |
[out] | G | this image will contain the six components of the structure tensor, each stored in a separate channel. The first channel will be the square of the estimated partial derivative in x , the second channel will be the estimated partial derivative in x multiplied by the estimated partial derivative in y . The third channel will be the square of the estimated partial derivative in y . The fourth cahnnel will be the estimated partial derivative in x multiplied by the estimated partial derivative in z . The fifth channel will be the estimated partial derivative in y multiplied by the estimated partial derivative in z . The sixth channel will be the square of the estimated partial derivative in z . |
[in] | square_root | if this is set to 'true', then the square root of the tensor will be calculated in each pixel, rather than the tensor itself. |
void lsseg::compute_scale_factor_2D | ( | const Image< double > & | img, | |
Image< double > & | scale_accum, | |||
Image< double > & | time_accum, | |||
double | dt, | |||
double | T | |||
) |
compute the 'scale factor' of an image. The scale factor gives some information about the global scale of the structure that each pixel belongs to. The scale factor and how it is computed is explained by Thomas Brox. Basically, it is determined by how much each pixel changes during an anisotropic smoothing process defined by a time-dependent partial differential equation.
img | the image for which we want to compute the scale factor. | |
scale_accum | an image containing the accumulated scale of each pixel in img over the total 'simulation time'. | |
time_accum | an image containing, for each pixel, the accumulated 'simulation time' during which its accumulated scale value has changed. | |
dt | simulation time step used when computing the scale factor. | |
T | total simulation time used. |
void lsseg::compute_scale_factor_3D | ( | const Image< double > & | img, | |
Image< double > & | scale_accum, | |||
Image< double > & | time_accum, | |||
double | dt, | |||
double | T | |||
) |
compute the 'scale factor' of an image. The scale factor gives some information about the global scale of the structure that each pixel belongs to. The scale factor and how it is computed is explained by Thomas Brox. Basically, it is determined by how much each pixel changes during an anisotropic smoothing process defined by a time-dependent partial differential equation.
[in] | img | the image for which we want to compute the scale factor. |
[out] | scale_accum | an image containing the accumulated scale of each pixel in img over the total 'simulation time'. |
[out] | time_accum | an image containing, for each pixel, the accumulated 'simulation time' during which its accumulated scale value has changed. |
[in] | dt | simulation time step used when computing the scale factor. |
[in] | T | total simulation time used. |
void lsseg::nonlinear_gauss_filter_2D | ( | const Image< double > & | img_old, | |
Image< double > & | img_new, | |||
double | dt, | |||
double | sigma, | |||
double | p | |||
) |
Makes a smoothed version of a given image. The smoothing is nonlinear, meaning that the edges of the images will be smoothed differently than regular areas. The process is based on the time-dependent partial differntial equation , where represent the image pixel values and is a fixed scalar value. An explanation can be found in section 2.1 of Thomas Brox' thesis, particularly under 2.1.2.
[in] | img_old | the original image that we want to smooth. |
[out] | img_new | the resulting image after smoothing |
[in] | dt | the time step to be used in the smoothing process |
[in] | sigma | strength of (nonlinear) gaussian pre-smoothing of the image before applying the nonlinear smoothing process. (Recommended: a small, nonzero value). |
[in] | p | this is the user-given of the equation above. It can be seen as a "nonlinearity" parameter. A value of 0 yields pure gaussian smoothing. A value of 1 yields the TVD (total-variation-diminishing) smoothing. Values beyond 1 yields backward diffusion, meaning that edges will actually be enhanced (edge-enhancing flow). |
void lsseg::nonlinear_gauss_filter_3D | ( | const Image< double > & | img_old, | |
Image< double > & | img_new, | |||
double | dt, | |||
double | sigma, | |||
double | p | |||
) |
Makes a smoothed version of a given image. The smoothing is nonlinear, meaning that the edges of the images will be smoothed differently than regular areas. The process is based on the time-dependent partial differntial equation , where represent the image pixel values and is a fixed scalar value. An explanation can be found in section 2.1 of Thomas Brox' thesis, particularly under 2.1.2.
[in] | img_old | the original image that we want to smooth. |
[out] | img_new | the resulting image after smoothing |
[in] | dt | the time step to be used in the smoothing process |
[in] | sigma | strength of (nonlinear) gaussian pre-smoothing of the image before applying the nonlinear smoothing process. (Recommended: a small, nonzero value). |
[in] | p | this is the user-given of the equation above. It can be seen as a "nonlinearity" parameter. A value of 0 yields pure gaussian smoothing. A value of 1 yields the TVD (total-variation-diminishing) smoothing. Values beyond 1 yields backward diffusion, meaning that edges will actually be enhanced (edge-enhancing flow). |
void lsseg::anisotropic_smoothing | ( | const Image< double > & | img_old, | |
Image< double > & | img_new, | |||
double | dt, | |||
double | sigma = 0 , |
|||
double | p = 1 | |||
) |
A wrapper function for nonlinear_gauss_filter_2D and nonlinear_gauss_filter_3D. It checks whether the input image is 2D or 3D, and calls the respective function accordingly.
[in] | img_old | the original image that we want to smooth. |
[out] | img_new | the resulting image after smoothing |
[in] | dt | the time step to be used in the smoothing process |
[in] | sigma | strength of (nonlinear) gaussian pre-smoothing of the image before applying the nonlinear smoothing process. (Recommended: a small, nonzero value). |
[in] | p | this is the user-given of the equation above. It can be seen as a "nonlinearity" parameter. A value of 0 yields pure gaussian smoothing. A value of 1 yields the TVD (total-variation-diminishing) smoothing. Values beyond 1 yields backward diffusion, meaning that edges will actually be enhanced (edge-enhancing flow). |
void lsseg::compute_smoothing_geometry_3D | ( | const Image< double > & | G, | |
double | p1, | |||
double | p2, | |||
double | p3, | |||
Image< double > & | T, | |||
bool | take_square_root = false | |||
) |
Compute the smoothing geometry tensor field T from a (blured) structure tensor field G, i.e., , where , and are the three eigenvectors of the structure tensor in each point. correspond here to the eigenvector with the smallest eigenvalue and to the one with the largest. The multiplicative factors , and are computed from the eigenvalues of the structure tensor, , and in the following way:
For explanation of use, read section 2.1.5 of D. Tschumperle's thesis, as well as 2.1 of [Tschumperle06].
G | the input structure tensor field | |
p1 | Power of the inverse of the multiplicative factor for the direction with the smallest eigenvalue for G. Should generally be smaller than p2. | |
p2 | Power of the inverse of the multiplicative factor for the direction with the intermediate eigenvalue for G. Should generally be larger than p1. | |
p3 | Power of the inverse of the multiplicative factor for the direction with the largest eigenvalue for G. Should generally be larger than p1 and p2. | |
T | The resulting structure tensor field (6 components) | |
take_square_root | If 'true', the generated structure tensor field will be sqrt(T) rather than T. |
void lsseg::compute_smoothing_geometry_2D | ( | const Image< double > & | G, | |
double | p1, | |||
double | p2, | |||
Image< double > & | T, | |||
bool | take_square_root = false | |||
) |
Compute the smoothing geometry tensor field T from a (blured) structure tensor field G, i.e., , where and are the eigenvectors of the structure tensor in each point. correspond to the direction with the smallest eigenvalue and to the one with the largest. The multiplicative factors and are computed from the eigenvalues of the structure tensor, and in the following way:
For explanation of use, read section 2.1.5 of D. Tschumperle's thesis, as well as 2.1 of [Tschumperle06].
G | the input structure tensor field | |
p1 | Power of the inverse of the multiplicative factor for the direction with the smallest eigenvalue. Should generally be smaller than p2. | |
p2 | Power of the inverse of the multiplicative factor for the direction with the largest eigenvalue. Should generally be larger than p1. | |
T | The resulting structure tensor field | |
take_square_root | If 'true', the generated structure tensor field will be sqrt(T) rather than T. |
void lsseg::normal_direction_flow_2D | ( | const LevelSetFunction & | phi, | |
LevelSetFunction & | advect, | |||
BorderCondition | bcond, | |||
double & | H1, | |||
double & | H2, | |||
const Mask *const | changes_allowed_mask = 0 , |
|||
const Mask *const | defined_region_mask = 0 | |||
) |
Compute the change that must be made to a LevelSetFunction when it undergoes one timestep of a normal direction flow, under an externally imposed force field. The length of the timestep is not specified, but is taken to be the maximum stable one. As we are in fact solving a Hamilton-Jacobi equation, the maximum allowable timestep can be computed by satisfying the CFL condition: . The book Level set methods and dynamic implicit surfaces has a good primer on Hamilton-Jacobi equations and their numerical discretisation, and the reader is encouraged to consult it for details.
[in] | phi | the level-set function that should undergo the normal-direction flow |
[in,out] | advect | on input, this represent the normal force field, which should have the same dimensions as 'phi' (one value per pixel). Upon return of the function, it will contain the "modification function" that should be multiplied by the timestep and added to 'phi'. |
[in] | bcond | specify the boundary condition to use at the edges of the domain (or the active parts of the domain, in case we use a mask). |
[out] | H1 | returns the max partial derivative of the Hamiltonian with respect to x . This is used in order to compute the maximum allowable timestep, according to the CFL condition presented above. |
[out] | H2 | returns the max partial derivative of the Hamiltonian with respect to y . This is used in order to compute the maximum allowable timestep, according to the CFL condition presented above. |
[in] | changes_allowed_mask | If nonzero, a mask defining which parts of the LevelSetFunction 'phi' that are subject to change for this timestep. Values corresponding to inactive parts of this mask will not get any update value in 'advect' upon return of the function. If zero, the all of (the defined part of) 'phi' is subject to change in this timestep. |
[in] | defined_region_mask | If nonzero, a mask defining which parts of (the domain of) the LevelSetFunction 'phi' are really defined. Only these parts will be used when computing update values. The boundary of the defined parts will be treated using the boundary condition specified in 'bcond'. If zero, then the whole domain of 'phi' is considered as defined. |
changes_allowed_mask
and the defined_region_mask
that deserves to be pointed out. The defined_region_mask
really specifies which pixels of the LevelSetFunction that contain any meaningful value. I.e. when estimating derivatives, etc., the algorithm must take care never to use any pixel that lies outside this domain. Moreover, it is of no use to compute updates for pixels outside this domain, as these pixels do not make sense anyway. On the other hand, the changes_allowed_mask
expresses which parts of (the defined region of) the LevelSetFunction that are subject to change during this timestep. I.e. updates will only be computed for pixels lying inside this region. However, pixels outside this region can be used for the algorithm's computational needs (estimation of derivatives, etc.). Definition at line 165 of file level_set.C.
void lsseg::normal_direction_flow_3D | ( | const LevelSetFunction & | phi, | |
LevelSetFunction & | advect, | |||
BorderCondition | bcond, | |||
double & | H1, | |||
double & | H2, | |||
double & | H3, | |||
const Mask *const | changes_allowed_mask = 0 , |
|||
const Mask *const | defined_region_mask = 0 | |||
) |
Compute the change that must be made to a LevelSetFunction when it undergoes one timestep of a normal direction flow, under an externally imposed force field. The length of the timestep is not specified, but is taken to be the maximum stable one. As we are in fact solving a Hamilton-Jacobi equation, the maximum allowable timestep can be computed by satisfying the CFL condition: . The book Level set methods and dynamic implicit surfaces has a good primer on Hamilton-Jacobi equations and their numerical discretisation, and the reader is encouraged to consult it for details.
[in] | phi | the level-set function that should undergo the normal-direction flow |
[in,out] | advect | on input, this represent the normal force field, which should have the same dimensions as 'phi' (one value per pixel). Upon return of the function, it will contain the "modification function" that should be multiplied by the timestep and added to 'phi'. |
[in] | bcond | specify the boundary condition to use at the edges of the domain (or the active parts of the domain, in case we use a mask). |
[out] | H1 | returns the max partial derivative of the Hamiltonian with respect to x . This is used in order to compute the maximum allowable timestep, according to the CFL condition presented above. |
[out] | H2 | returns the max partial derivative of the Hamiltonian with respect to y . This is used in order to compute the maximum allowable timestep, according to the CFL condition presented above. |
[out] | H3 | returns the max partial derivative of the Hamiltonian with respect to z . This is used in order to compute the maximum allowable timestep, according to the CFL condition presented above. |
[in] | changes_allowed_mask | If nonzero, a mask defining which parts of the LevelSetFunction 'phi' that are subject to change for this timestep. Values corresponding to inactive parts of this mask will not get any update value in 'advect' upon return of the function. If zero, the all of (the defined part of) 'phi' is subject to change in this timestep. |
[in] | defined_region_mask | If nonzero, a mask defining which parts of (the domain of) the LevelSetFunction 'phi' are really defined. Only these parts will be used when computing update values. The boundary of the defined parts will be treated using the boundary condition specified in 'bcond'. If zero, then the whole domain of 'phi' is considered as defined. |
changes_allowed_mask
and the defined_region_mask
that deserves to be pointed out. The defined_region_mask
really specifies which pixels of the LevelSetFunction that contain any meaningful value. I.e. when estimating derivatives, etc., the algorithm must take care never to use any pixel that lies outside this domain. Moreover, it is of no use to compute updates for pixels outside this domain, as these pixels do not make sense anyway. On the other hand, the changes_allowed_mask
expresses which parts of (the defined region of) the LevelSetFunction that are subject to change during this timestep. I.e. updates will only be computed for pixels lying inside this region. However, pixels outside this region can be used for the algorithm's computational needs (estimation of derivatives, etc.). Definition at line 231 of file level_set.C.
void lsseg::visualize_multisets | ( | const LevelSetFunction **const | phis, | |
int | num_phi, | |||
Image< double > & | target, | |||
const int **const | rgb_color | |||
) |
Creates a color image that visualizes a set of LevelSetFuncions, using a different color for the inside of each LevelSetFunction, and another color for regions outside of all LevelSetFunctions.
[in] | phis | pointer to array of LevelSetFunctions to visualize. |
[in] | num_phi | number of elements in the array pointed to by 'phis'. (I.e., number of LevelSetFunctions to visualize). |
[out] | target | the produced image |
[in] | rgb_color | pointer to an array of colors, where each color is defined as a 3-array of RGB values. There should be (num_phi + 1) colors in this array; the last one being the color for regions that fall outside of all given level-set functions. |
Definition at line 305 of file level_set.C.
int lsseg::visualize_level_set | ( | const LevelSetFunction & | phi, | |
Image< double > & | target, | |||
double | threshold, | |||
double | outside_intensity = 0 , |
|||
double | curve_intensity = 255 , |
|||
double | interior_intensity = 50 , |
|||
const Mask *const | mask = 0 , |
|||
double | mask_intensity = 255 | |||
) |
Creates a black-and-white image that visualizes a LevelSetFunction, using one intensity for the inside region, another for the outside region, and a third intensity for the pixels on the boundary between the two regions.
[in] | phi | the level-set function to visualize |
[out] | target | the produced image |
[in] | threshold | the threshold between what is considered 'inside' and 'outside' of the level-set function (usually zero). |
[in] | outside_intensity | image intensity of the region considered 'outside' of 'phi' |
[in] | curve_intensity | image intensity for pixels lying on the boundary between the 'inside' and 'outside' regions. |
[in] | interior_intensity | image intensity of the region considered 'inside' of 'phi'. |
[in] | mask | if nonzero, this mask specifies what is considered the defined part of the domain 'phi'. If zero, then the whole domain of 'phi' is considered defined. |
[in] | mask_intensity | image intensity of the undefined region of 'phi', as specified by 'mask'. |
Definition at line 319 of file level_set.C.
void lsseg::LIC_2D_FS | ( | const Image< double > & | src, | |
const Image< double > & | vec, | |||
Image< double > & | target, | |||
double(*)(double) | kernel_func, | |||
const double | dl = 0.8 , |
|||
const double | L = 10 | |||
) |
Carry out line integral convolution of a 2D-image with a kernel function, along streamlines given by a user-specified vector field. The result is an image smoothed along the streamlines, with a smoothing defined by the kernel function. For more on line integral convolution, refer to this paper [Cabral93], or (more adapted for the curvature-preserving smoothing algorithm): this paper [Tschumperle06].
[in] | src | the Image that will participate in the line integral convolution |
[in] | vec | a vector field defining the streamlines along which convolution will take place. The vector field is specified as an Image with three channels; the two first containing the x- and y-components of a normalized vector field, the third containing the magnitudes of the vectors in the field. |
[out] | target | the Image obtained when convoluting 'src' with 'kernel_func' along the streamlines of 'vec' (the result of the line integral convolution. |
[in] | kernel_func | pointer to the kernel function (a function from to . This is typically a symmetrical function with compact support, or with infinite support but with a value that goes off to zero relatively quickly. |
[in] | dl | steplength used for integration purposes |
[in] | L | total length along which to integrate (in each direction from zero). The important part of the support of 'kernel_func' should fall within this range. |
void lsseg::LIC_3D_FS | ( | const Image< double > & | src, | |
const Image< double > & | vec, | |||
Image< double > & | target, | |||
double(*)(double) | kernel_func, | |||
const double | dl = 0.8 , |
|||
const double | L = 10 | |||
) |
Carry out line integral convolution of a 3D-image with a kernel function, along streamlines given by a user-specified vector field. The result is an image smoothed along the streamlines, with a smoothing defined by the kernel function. For more on line integral convolution, refer to this paper [Cabral93], or (more adapted for the curvature-preserving smoothing algorithm): this paper [Tschumperle06].
[in] | src | the Image that will participate in the line integral convolution |
[in] | vec | a vector field defining the streamlines along which convolution will take place. The vector field is specified as an Image with four channels; the three first containing the x- , y- and z-components of a normalized vector field, the third containing the magnitudes of the vectors in the field. |
[out] | target | the Image obtained when convoluting 'src' with 'kernel_func' along the streamlines of 'vec' (the result of the line integral convolution. |
[in] | kernel_func | pointer to the kernel function (a function from to . This is typically a symmetrical function with compact support, or with infinite support but with a value that goes off to zero relatively quickly. |
[in] | dl | steplength used for integration purposes |
[in] | L | total length along which to integrate (in each direction from zero). The important part of the support of 'kernel_func' should fall within this range. |
void lsseg::LIC_3D | ( | const Image< double > & | src, | |
const Image< double > & | vec, | |||
Image< double > & | target, | |||
double(*)(double) | kernel_func, | |||
const double | L = 10 | |||
) |
Carry out line integral convolution of a 3D-image with a kernel function, along streamlines given by a user-specified vector field. The result is an image smoothed along the streamlines, with a smoothing defined by the kernel function. During integration, this function uses an adaptive steplength. For more on line integral convolution, refer to this paper [Cabral93], or (more adapted for the curvature-preserving smoothing algorithm): this paper [Tschumperle06].
[in] | src | the Image that will participate in the line integral convolution |
[in] | vec | a vector field defining the streamlines along which convolution will take place. The vector field is specified as an Image with three channels, expressing the x-, y- and z-components of the vector field. |
[out] | target | the Image obtained when convoluting 'src' with 'kernel_func' along the streamlines of 'vec' (the result of the line integral convolution. |
[in] | kernel_func | pointer to the kernel function (a function from to . This is typically a symmetrical function with compact support, or with infinite support but with a value that goes off to zero relatively quickly. |
[in] | L | total length along which to integrate (in each direction from zero). The important part of the support of 'kernel_func' should fall within this range. |
void lsseg::LIC_2D | ( | const Image< double > & | src, | |
const Image< double > & | vec, | |||
Image< double > & | target, | |||
double(*)(double) | kernel_func, | |||
const double | L = 10 | |||
) |
Carry out line integral convolution of a 2D-image with a kernel function, along streamlines given by a user-specified vector field. The result is an image smoothed along the streamlines, with a smoothing defined by the kernel function. During integration, this function uses an adaptive steplength. For more on line integral convolution, refer to this paper [Cabral93], or (more adapted for the curvature-preserving smoothing algorithm): this paper [Tschumperle06].
[in] | src | the Image that will participate in the line integral convolution |
[in] | vec | a vector field defining the streamlines along which convolution will take place. The vector field is specified as an Image with two channels, expressing the x- and y-components of the vector field. |
[out] | target | the Image obtained when convoluting 'src' with 'kernel_func' along the streamlines of 'vec' (the result of the line integral convolution. |
[in] | kernel_func | pointer to the kernel function (a function from to . This is typically a symmetrical function with compact support, or with infinite support but with a value that goes off to zero relatively quickly. |
[in] | L | total length along which to integrate (in each direction from zero). The important part of the support of 'kernel_func' should fall within this range. |
double lsseg::develop_multiregion_2D | ( | Region * | regs, | |
int | num_regs, | |||
int | num_iter, | |||
int | reinit_modulo, | |||
const Mask * | geom_mask | |||
) |
This function is the bread-and-butter of the level-set based image segmentation progess. It runs a multi-region segmentation on an Image according to rules specified by a specific set of ForceGenerators, and with a specified set of boundary smoothness criteria (scalar values). The segmentation is carried out by the development of a time-dependent partial differential equation. A certain, user-specified, number of iterations will be carried out - there is no stop-criterion or other convergence checks (such checks would be quite computationally expensive). The Image with its initial segmentation, as well as the ForceGenerators and the smoothness criteria are given to the function through a set of Region objects, each representing a particular segmented region of the image. The resulting segmentation will be returned through the same Region objects. For an overview of the segmentation procedure, read our text on how segmentation works, which also explains the general details on multi-region segmentation. This function is optimized for two-dimensional Images.
[in,out] | regs | 'regs' points to an array of Regions, each representing one of the regions that will be developing in the image. The information contained in the Region's structure gives the driving forces behind the segmentation for that Region, its initial segmentation and smoothness criterion. After calling this function, each Region will have its LevelSetFunction changed to represent its resulting segmentation. |
[in] | num_regs | Number of Regions in the array pointed to by 'regs' (ie. number of "competing regions" in the segmentation process). |
[in] | num_iter | number of iterations (timesteps) of applying the PDE-based process for developing the segmentation. Please note that the length of each timestep can vary and is not given here. The timestep is chosen by the algorithm to be the largest one that is guaranteed to be stable. |
[in] | reinit_modulo | Ever so often, the LevelSetFunction defining the image segmentation should be reinitialized. This stabilizes the segmentation process. Reinitialization is quite expensive, however, and it is not in general necessary to run it after each iteration of applying the PDE. This parameter lets the user decide how many iterations of applying the PDE should be carried out between each reinitialization. A value of 0 to this argument means that no reinitalization should ever be carried out. |
[in] | geom_mask | if nonzero, the Mask pointed to by this argument specifies which parts of the original Image that should participate in the segmentation process, the rest of the image will not be taken into account. If this pointer is zero, then the whole image participates in the segmentation process. |
Some ForceGenerators have different modes for normal, two-region segmentation and for multi-region segmentation. Make sure they are constructed with the multi-region mode if they are to be used in the context of this function.
Definition at line 98 of file MultiRegionAlgorithm.C.
double lsseg::develop_multiregion_3D | ( | Region * | regs, | |
int | num_regs, | |||
int | num_iter, | |||
int | reinit_modulo, | |||
const Mask * | geom_mask | |||
) |
This function is the bread-and-butter of the level-set based image segmentation progess. It runs a multi-region segmentation on an Image according to rules specified by a specific set of ForceGenerators, and with a specified set of boundary smoothness criteria (scalar values). The segmentation is carried out by the development of a time-dependent partial differential equation. A certain, user-specified, number of iterations will be carried out - there is no stop-criterion or other convergence checks (such checks would be quite computationally expensive). The Image with its initial segmentation, as well as the ForceGenerators and the smoothness criteria are given to the function through a set of Region objects, each representing a particular segmented region of the image. The resulting segmentation will be returned through the same Region objects. For an overview of the segmentation procedure, read our text on how segmentation works, which also explains the general details on multi-region segmentation. This function is optimized for three-dimensional Images.
[in,out] | regs | 'regs' points to an array of Regions, each representing one of the regions that will be developing in the image. The information contained in the Region's structure gives the driving forces behind the segmentation for that Region, its initial segmentation and smoothness criterion. After calling this function, each Region will have its LevelSetFunction changed to represent its resulting segmentation. |
[in] | num_regs | Number of Regions in the array pointed to by 'regs' (ie. number of "competing regions" in the segmentation process). |
[in] | num_iter | number of iterations (timesteps) of applying the PDE-based process for developing the segmentation. Please note that the length of each timestep can vary and is not given here. The timestep is chosen by the algorithm to be the largest one that is guaranteed to be stable. |
[in] | reinit_modulo | Ever so often, the LevelSetFunction defining the image segmentation should be reinitialized. This stabilizes the segmentation process. Reinitialization is quite expensive, however, and it is not in general necessary to run it after each iteration of applying the PDE. This parameter lets the user decide how many iterations of applying the PDE should be carried out between each reinitialization. A value of 0 to this argument means that no reinitalization should ever be carried out. |
[in] | geom_mask | if nonzero, the Mask pointed to by this argument specifies which parts of the original Image that should participate in the segmentation process, the rest of the image will not be taken into account. If this pointer is zero, then the whole image participates in the segmentation process. |
Some ForceGenerators have different modes for normal, two-region segmentation and for multi-region segmentation. Make sure they are constructed with the multi-region mode if they are to be used in the context of this function.
Definition at line 211 of file MultiRegionAlgorithm.C.
void lsseg::negate | ( | Image< double > & | img, | |
double | min = 0 , |
|||
double | max = 255 | |||
) |
make an Image<double> a negative of itself.
[in] | min | the min. image pixel value |
[in] | max | the max. image pixel value |
[in,out] | img | the Image<double> which will be turned into its own negative |
min
and max
Definition at line 87 of file simple_tools.C.
void lsseg::rescale | ( | Image< double > & | img, | |
double | cur_min, | |||
double | cur_max, | |||
double | to_min = 0 , |
|||
double | to_max = 255 | |||
) |
rescale the pixel range of an Image<double> so that its pixels spans a different range than before rescaling
[in] | cur_min | the current min. of the image's pixel range |
[in] | cur_max | the current max. of the image's pixel range |
[in] | to_min | the min. of the image's pixel range after rescaling |
[in] | to_max | the max. of the image's pixel range after rescaling |
[out] | img | the image whose pixels are to be rescaled |
cur_min
and cur_max
parameters could possibly be directly determined from the image. However, this is not always desirable. Firstly, it is not sure that the image spans the complete possible range (for instance, many pictures have a formal pixel value range from 0 to 255, but where no pixel in the image attains these extreme values). Secondly, computing the range directly from the image is costly. Definition at line 97 of file simple_tools.C.
void lsseg::rescale_channels | ( | Image< double > & | img, | |
double | to_min = 0 , |
|||
double | to_max = 255 | |||
) |
rescale the pixel ranges of each of an image's channels to a user-specified interval
NB Upon function completion, each channel of the image will have a pixel range that spans the full interval specified by the user.
[in] | to_min | the min. value of each channel's pixel range after rescaling |
[in] | to_max | the max. value of each channel's pixel range after rescaling |
[in,out] | img | the image whose pixel ranges are to be rescaled |
Definition at line 111 of file simple_tools.C.
void lsseg::clip | ( | Image< double > & | img, | |
double | min, | |||
double | max | |||
) |
clip the pixel range of a image to a certain interval
Pixel values below/above the user-specified interval will be set to the interval's min/max value respectively.
[in,out] | img | image whose pixels will be clipped to the specified range |
[in] | min | min. bound of the user-specified range |
[in] | max | max. bound of the user-specified range |
Definition at line 136 of file simple_tools.C.
void lsseg::transpose | ( | Image< double > & | img | ) |
transpose the x
and y
dimensions of an Image<double>
The z-dimension
will remain unchanged. The transpose operation is of course carried out on all image channels.
[in,out] | img | the image to be transposed |
Definition at line 151 of file simple_tools.C.
void lsseg::transpose | ( | const Image< double > & | img, | |
Image< double > & | target | |||
) |
Generate an Image<double> that is the transpose of another Image<double>.
The transpose is taking place in the x
and y
dimensions, with the z-dimension
remaining unchanged.
[in] | img | the original image |
[out] | target | upon function completion, target will contain the transposed version of img . |
Definition at line 160 of file simple_tools.C.
void lsseg::horizontal_sinusoidal_bands | ( | LevelSetFunction & | img, | |
int | num_bands, | |||
double | phase = 0 | |||
) |
Initializes a LevelSetFunction with horizontal bands created by a sine function.
This is a simple way of creating a smooth LevelSetFunction whose two regions consist of a number of horizontal bands (the positive vs. the negative regions of the sine function). This can be used as an initialization of a LevelSetFunction used for segmentation purposes, since it starts out with a "segmentation" where the two regions are somewhat equally distributed over the image.
[in,out] | img | the LevelSetFunction we want to fill with horizontal sinusoidal bands |
[in] | num_bands | number of distinct bands in the image |
[in] | phase | use this variable to change the phase of the sine function and thus specify exactly the vertical position of the bands. The phase is specified in multiples of , so that phase = 0.5 gives a shift of and phase = 1 gives the same image as phase = 0 |
Definition at line 176 of file simple_tools.C.
void lsseg::rectangle | ( | LevelSetFunction & | img, | |
double | xmin_ratio, | |||
double | xmax_ratio, | |||
double | ymin_ratio, | |||
double | ymax_ratio, | |||
double | zmin_ratio = 0 , |
|||
double | zmax_ratio = 1 | |||
) |
Initializes a LevelSetFunction to become 1 everywhere, except for an interior rectangular domain where the LevelSetFunction has the value -1.
This is a simple way of creating a LevelSetFunction where the closed region is a rectangle positioned by the user somewhere inside the image. As the resulting function is not smooth, it may be recommendable to reinitialize it before using it for segmentation purposes.
[in,out] | img | the LevelSetFunction we want to initialize |
[in] | xmin_ratio | position of the "leftmost" edge of the inner rectangular domain, as a fraction of total domain width of the LevelSetFunction; i.e. 0 means the leftmost edge of the domain. |
[in] | xmax_ratio | position of the "rightmost" edge of the inner rectangular domain, as a fraction of total domain width of the LevelSetFunction; i.e. 1 means the rightmost edge of the domain. |
[in] | ymin_ratio | As xmin_ratio , but concern the y -direction. |
[in] | ymax_ratio | As ymax_ratio , but concern the y -direction. |
[in] | zmin_ratio | As zmin_ratio , but concern the z -direction. |
[in] | zmax_ratio | As zmax_ratio , but concern the z -direction. |
xmin_ratio
< xmax_ratio
, and similar requirements for ymin_ratio
and zmin_ratio
. Definition at line 192 of file simple_tools.C.
void lsseg::init_voronoi_regions | ( | LevelSetFunction * | regs, | |
const double * | center_coords, | |||
int | num_regions, | |||
bool | three_d = false | |||
) |
Initializes a set of LevelSetFunctions defined on a common domain so that their closed regions together constitute a voronoi subdivision of the domain. The points defining the voronoi subdivision is given by the user.
Defined this way, the union of the closed regions will exactly cover the total domain, and the intersection of the closed regions will be empty. This can be a useful initialization for multiregion- segmentations where you have some idea about the position of each object.
[in,out] | regs | pointer to an array of LevelSetFunctions that will be initialized so that each of them represent one of the voronoi regions in the subdivision of the domain. NB:All the LevelSetFunctions in this array must have exactly the same shape, and there must be as many of them as the desired number of regions. |
[in] | center_coords | pointer to an array of 2D or 3D coordinates, which represent the points used as the basis for the voronoi subdivision of the domain. The number of points must be equal to the number of regions. The coordinates are stored in format (2D) or format (3D). |
[in] | num_regions | the number of voronoi regions into which the domain should be subdivided. |
[in] | three_d | specify whether the domain is 2D or 3D. |
Definition at line 224 of file simple_tools.C.
void lsseg::multiregion_bands | ( | LevelSetFunction * | regs, | |
int | num_regs, | |||
int | pixel_bandwidth | |||
) |
Initializes a set of LevelSetFunctions defined on a common domain so that their closed regions together constitute a subdivision on the domain, and where each such region consists of a set of horizontal bands of a certain width along the y-direction.
The union of the closed regions will exactly cover the whole domain, and the intersections of the closed regions will be empty.
[in,out] | regs | pointer to an array of LevelSetFunctions that will be initialized so that their closed domains consists of horizontal stripes, and form a disjoint partition of the domain. NB: All the LevelSetFunctions in this array must have exactly the same shape (the shape of the domain), and there must be as many of them as the desired number of regions, as expressed with num_regs . |
[in] | num_regs | The number of regions into which the domain should be subdivided. |
[in] | pixel_bandwidth | The number of pixels in the y -direction across each horizontal stripe that makes up a region. |
Definition at line 283 of file simple_tools.C.
void lsseg::random_scattered_voronoi | ( | LevelSetFunction * | regs, | |
int | num_regs, | |||
int | num_fragments | |||
) |
Initializes a set of LevelSetFunctions defined on a common domain so that their closed regions together constitute a voronoi subdivision of the domain. The points defining the voronoi subdivision are randomly generated.
Defined this way, the union of the closed regions will exactly cover the whole domain, and the intersection of the closed regions will be empty. The number of disjoint fragments of each region can be specified by the user (each such fragment constitutes a voronoi region in the total domain subdivision).
[in,out] | regs | pointer to an array of LevelSetFunctions that will be initialized so that each of them represent a certain number of disjoint voronoi regions. Taken together, the interior regions defined by regs will cover the whole domain. NB: All the LevelSetFunctions in this array must have exactly the same shape, and there must be as many of them as the desired number of regions, as expressed with num_regs |
[in] | num_regs | The number of regions into which the domain should be subdivided. |
[in] | num_fragments | The number of disjointed fragments (separate voronoi regions) that each region could be made up of. Should be at least 1. |
Definition at line 325 of file simple_tools.C.
void lsseg::set_from_parzen | ( | LevelSetFunction & | img, | |
const ParzenDistributionForce | pf, | |||
const Mask * | m = 0 | |||
) |
Initializes a LevelSetFunction to become a function with two values, 1 and -1, corresponding to the domains in which the force defined by the pf
argument is positive or negative.
This way of initalizing a LevelSetFunction before a segmentation may be useful in some cases, where the segmentation is strongly dependent on a fixed pixel distribution that can be specified beforehand (fixed in the ParzenDistributionForce).
[in,out] | img | the LevelSetFunction to initialize |
pf | the ForceGenerator that will directly decide which region each of the LevelSetFuction's pixels belong to. |
img
. m | An optional Mask defining which part of img that will be filled in, and which part will remain unchanged. If m is a zero-pointer, then no Mask will be used and the complete domain of img will be initialized. |
Definition at line 390 of file simple_tools.C.
void lsseg::sphere | ( | LevelSetFunction & | img, | |
double | relrad, | |||
double | xrelpos, | |||
double | yrelpos, | |||
double | zrelpos = 0 | |||
) |
Initializes a LevelSetFunction to become the signed distance function from a point in its domain minus some fixed value, thus its zero-set becomes a circle around .
[in,out] | img | the LevelSetFunction to initialize |
[in] | relrad | the relative radius of the sphere (circle in 2D) defined by img 's zero-set after initialization. 1 correspond to the length of the longest dimension of img , so usually you will want to give this a value between 0 and 1. |
[in] | xrelpos | the relative position of the x-coordinate of the center of the point as defined above. A value of 0 corresponds to the LevelSetFunction's domain's leftmost edge, and 1 to its rightmost. |
[in] | yrelpos | the relative position of the y-coordinate of the center of the point as defined above. Analoguous to xrelpos but concerning the y-dimension. |
[in] | zrelpos | the relative position of the z-coordinate of the center of the point as defined above. Analoguous to xrelpos but concerning the z-dimension. |
Definition at line 417 of file simple_tools.C.
void lsseg::make_border_mask | ( | const LevelSetFunction & | phi, | |
Mask & | target, | |||
int | width = 2 , |
|||
const Mask * | geom_mask = 0 | |||
) |
Generate a Mask over a domain, where the active region of the Mask is specified as the region of a certain width (in pixels) around the zero-set of a LevelSetFunction defined on the same domain.
[in] | phi | the LevelSetFunction defined over the domain in question. Its zero-set will be used in the definition of the active region of the Mask to create |
[out] | target | the Mask that will be defined. It will be resized so that its domain coincides with that of phi , and its active region will be set according to the description above. |
[in] | width | Specify the width of the active region. The width is specified in number of pixels around the zero-set of the LevelSetFunction given by phi . |
[in] | geom_mask | An optional pointer to a mask specifying the exact geometry of the domain to consider. If this pointer is 0, the domain to consider is the whole domain covered by phi . |
Definition at line 498 of file simple_tools.C.
void lsseg::mask_from_segmentation | ( | const LevelSetFunction & | phi, | |
Mask & | target, | |||
SEG_REGION | reg | |||
) |
Generate a Mask over a domain, where the active region of the Mask is specified as the region defined by the negative (or positive) region of a LevelSetFunction over the same domain.
[in] | phi | the LevelSetFunction specifying which part of the domain that shall constitute the Mask's active region. |
[out] | target | the Mask that will be generated |
[in] | reg | set this to SEG_NEGATIVE if you want the active region of the Mask to be the part of the domain where phi is negative. Set it to SEG_POSITIVE if you want the active region of the Mask to be the part of the domain where phi is positive. |
Definition at line 559 of file simple_tools.C.
void lsseg::read_image_sequence | ( | istream & | image_list, | |
Image< double > & | result, | |||
bool | convert_to_grayscale | |||
) |
Definition at line 586 of file simple_tools.C.
unsigned long int lsseg::nonzeroes | ( | const Image< int > & | img | ) |
Count the number of nonzero pixels of an Image<int>.
[in] | img | the image for which you want to count the number of nonzero pixels. |
returns | the number of nonzero pixels in img . |
Definition at line 667 of file simple_tools.C.
double lsseg::nonzero_ratio | ( | const Image< int > & | img | ) |
Calculate the ratio of the number pixels with a nonzero value to the total number of pixels in an Image<int>.
[in] | img | the image for which we want to compute this ratio. |
returns | the ratio of nonzero pixels to the total number of image pixels in img . |
Definition at line 680 of file simple_tools.C.
unsigned long int lsseg::positives | ( | const Image< double > & | img | ) |
Count the number of pixels with a nonnegative value in an Image<int>.
[in] | img | the image for which we want to compute the number of nonnegative pixels. |
returns | the number of nonnegative pixels in img . |
Definition at line 689 of file simple_tools.C.
unsigned long int lsseg::negatives | ( | const Image< double > & | img | ) |
Count the number of pixels with a negative value in an Image<int>.
[in] | img | the image for which we want to compute the number of negative pixels. |
returns | the number of negative pixels in img . |
Definition at line 702 of file simple_tools.C.
double lsseg::positive_ratio | ( | const Image< double > & | img | ) |
Calculate the ratio of the number of pixels with a nonnegative value to the total number of pixels in an Image<double>.
[in] | img | the image for which we want to compute this ratio |
returns | the ratio of nonnegative pixels to the total number of image pixels in img . |
Definition at line 709 of file simple_tools.C.
double lsseg::negative_ratio | ( | const Image< double > & | img | ) |
Calculate the ratio of the number of pixels with a negative value to the total number of pixels in an Image<double>.
[in] | img | the image for which we want to compute this ratio |
returns | the ratio of negative pixels to the total number of image pixels in img . |
Definition at line 716 of file simple_tools.C.
double lsseg::develop_single_region_2D | ( | Region & | reg, | |
int | num_iter, | |||
int | reinit_modulo = 0 , |
|||
const Mask * | geom_mask = 0 | |||
) |
This function is the bread-and-butter of the level-set based image segmentation progess. It runs a segmentation on an Image according to rules specified by a specific ForceGenerator, and with a specified boundary smoothness criterion (scalar value). The segmentation is carried out by the development of a time-dependent partial differential equation. A certain, user-specified, number of iterations will be carried out - there is no stop-criterion or other convergence checks (such checks would be quite computationally expensive). The Image with its initial segmentation, as well as the ForceGenerator and the smoothness criterion in question are all given to the function through a Region object, and the resulting segmentation will be returned as part of the same Region object. For an overview of the segmentation procedure, read our text on how segmentation works. This function is optimized for two-dimensional Images.
[in,out] | reg | Through this argument, the image, the initial segmentation, the driving force behind the segmentation and the smoothness criterion are specified. After calling this function, the LevelSetFunction defining the segmentation will have changed to the resulting one. |
[in] | num_iter | number of iterations (timesteps) of applying the PDE-based process for developing the segmentation. Please note that the length of each timestep can vary and is not given here. The timestep is chosen by the algorithm to be the largest one that is guaranteed to be stable. |
[in] | reinit_modulo | Ever so often, the LevelSetFunction defining the image segmentation should be reinitialized. This stabilizes the segmentation process. Reinitialization is quite expensive, however, and it is not in general necessary to run it after each iteration of applying the PDE. This parameter lets the user decide how many iterations of applying the PDE should be carried out between each reinitialization. A value of 0 to this argument means that no reinitalization should ever be carried out. |
[in] | geom_mask | if nonzero, the Mask pointed to by this argument specifies which parts of the original Image that should participate in the segmentation process, the rest of the image will not be taken into account. If this pointer is zero, then the whole image participates in the segmentation process. |
Definition at line 65 of file SingleRegionAlgorithm.C.
double lsseg::develop_single_region_3D | ( | Region & | reg, | |
int | num_iter, | |||
int | reinit_modulo = 0 , |
|||
const Mask * | geom_mask = 0 | |||
) |
This function is the bread-and-butter of the level-set based image segmentation progess. It runs a segmentation on an Image according to rules specified by a specific ForceGenerator, and with a specified boundary smoothness criterion (scalar value). The segmentation is carried out by the development of a time-dependent partial differential equation. A certain, user-specified, number of iterations will be carried out - there is no stop-criterion or other convergence checks (such checks would be quite computationally expensive). The Image with its initial segmentation, as well as the ForceGenerator and the smoothness criterion in question are all given to the function through a Region object, and the resulting segmentation will be returned as part of the same Region object. For an overview of the segmentation procedure, read our text on how segmentation works. This function is optimized for three-dimensional Images.
[in,out] | reg | Through this argument, the image, the initial segmentation, the driving force behind the segmentation and the smoothness criterion are specified. After calling this function, the LevelSetFunction defining the segmentation will have changed to the resulting one. |
[in] | num_iter | number of iterations (timesteps) of applying the PDE-based process for developing the segmentation. Please note that the length of each timestep can vary and is not given here. The timestep is chosen by the algorithm to be the largest one that is guaranteed to be stable. |
[in] | reinit_modulo | Ever so often, the LevelSetFunction defining the image segmentation should be reinitialized. This stabilizes the segmentation process. Reinitialization is quite expensive, however, and it is not in general necessary to run it after each iteration of applying the PDE. This parameter lets the user decide how many iterations of applying the PDE should be carried out between each reinitialization. A value of 0 to this argument means that no reinitalization should ever be carried out. |
[in] | geom_mask | if nonzero, the Mask pointed to by this argument specifies which parts of the original Image that should participate in the segmentation process, the rest of the image will not be taken into account. If this pointer is zero, then the whole image participates in the segmentation process. |
Definition at line 133 of file SingleRegionAlgorithm.C.
void lsseg::combine_channel_images | ( | const Image< T > ** | channel_array, | |
int | num_input_images, | |||
Image< T > & | result | |||
) |
Function that takes each image channel of a given number of images, and combines them all into a new image.
[in] | channel_array | array of pointers to the images whose channels will participate in the new image. It is of course required that all the images refered to are spatially compatible. |
[in] | num_input_images | specify the size of the array of pointers in channel_array . |
[out] | result | the resulting image, created from combining all of the channels in the input images. |
Definition at line 66 of file image_utils.h.
void lsseg::to_grayscale | ( | Image< T > & | img | ) |
convert a three-channel Image (RGB of double
, int
, etc.) to greyscale.
[in,out] | img | the image to convert to greyscale. Should have exactly three channels. Upon function completion, it will have only one channel, which will be a weighed average of those three. |
Definition at line 167 of file simple_tools_templates.h.
void lsseg::read_image_sequence | ( | std::istream & | image_list, | |
Image< double > & | result, | |||
bool | convert_to_grayscale = false | |||
) |
Read a set of 2D-image files and write the result into a 3D Image<double>, where the read images will be stacked along the z-coordinate direction.
[in] | image_list | a stream containing the list of names of the image files to read. The names should be in ASCII, separated by whitespace. NB: All the image must have the same shape (resolution) and number of channels. The fileformats are specified by the filename suffixes (.png, .gif, .jpeg, etc..) |
[out] | result | The images will be written to this Image<double>, stacked along the z -direction. |
[in] | convert_to_grayscale | set this to true if you want to convert the read images into grayscale before writing them to result . If you do this, result will only contain one channel upon return. |
void lsseg::resample_into | ( | const ImgType & | src, | |
ImgType & | target, | |||
bool | linear = true | |||
) |
Resample an image, using nearest-pixel or linear interpolation.
The function is a template on the image type, ImgType
, which should be an Image<> of some kind (LevelSetFunction, Mask, etc..).
template
argument here and not a reference to a base class... At the time of writing of this documentation, I don't remember any more. [in] | src | The image to resample |
[out] | target | Before calling this function, set target to the desired shape of the resampled image (in order to communicate it to this function). On function completion, it will contain the resampled image. |
[in] | linear | set this to true for linear pixel interpolation, otherwise, the resampling will make use of the "nearest pixel" in the image before resampling when determining the value of a pixel in the resampled image. |
Definition at line 130 of file simple_tools_templates.h.
void lsseg::downsample_series | ( | const ImgType & | input, | |
std::vector< ImgType > & | result, | |||
int | min_num_pixels, | |||
bool | downscale_z = false , |
|||
double | downscale_factor = 2 , |
|||
bool | to_grayscale = false , |
|||
bool | linear = true | |||
) |
Create a sequence of progressively smaller, resampled copies of a reference image.
The function is a template on the image type, ImgType
, which should be an Image<> of some kind (LevelSetFunction, Mask, etc..).
template
argument here and not a reference to a base class... At the time of writing of this documentation, I don't remember any more. [in] | input | the image we will use as reference image for generating the sequence. |
[out] | result | When calling the function, this should be an empty vector. Upon function completion, it will contain a number of progressively smaller resampled versions of the image input . |
[in] | min_num_pixels | This is a criterion for determining how many progressively smaller images will be generated. When the resampled image would have a coordinate direction with a resolution less than min_num_pixel , the function stops, and this image will not be generated. NB: The z -coordinate direction will only be taken into consideration here if downscale_z = true . |
[in] | downscale_z | If 'true' , downscaling will also take place in the z -direction, otherwise the image will only be downscaled in the x - and y -directions. Do remember that if the image is only 2D (ie. with a z -dimension of length 1), setting this parameter to true would result in the functions immediate return (no resampled images generated at all), due to the stop criterion explained in the explanation of the parameter min_num_pixels |
[in] | downscale_factor | Specify the downscale factor. If downscale_factor is equal to 2, the resolution along the x - and y - (and possibly z -) directions will be halved from one image in the sequence to the next (divided by 2). If downscale_factor = 1.5 , the resolution along each direction will decrease with 50% from one image in the sequence to the next (divided by 1.5), etc. |
[in] | to_grayscale | Specify whether to convert the images to grayscale before adding them to the result vector. |
[in] | linear | set this to true for linear pixel interpolation, otherwise, the resampling will make use of the "nearest pixel" in the image before resampling when determining the value of a pixel in the resampled image. |
Definition at line 75 of file simple_tools_templates.h.
const int lsseg::RED[] = {255, 0, 0} |
const int lsseg::BLUE[] = {0, 255,0} |
3-tuple of int
representing the color blue in RGB format.
Definition at line 59 of file colordefs.h.
const int lsseg::GREEN[] = {0, 0, 255} |
3-tuple of int
representing the color green in RGB format.
Definition at line 61 of file colordefs.h.
const int lsseg::YELLOW[] = {255,0, 255} |
3-tuple of int
representing the color yellow in RGB format.
Definition at line 63 of file colordefs.h.
const int lsseg::CYAN[] = {0, 255, 255} |
3-tuple of int
representing the color cyan in RGB format.
Definition at line 65 of file colordefs.h.
const int lsseg::MAGENTA[] = {255, 255, 0} |
3-tuple of int
representing the color magenta in RGB format.
Definition at line 67 of file colordefs.h.
const int lsseg::WHITE[] = {255,255, 255} |
3-tuple of int
representing the color white in RGB format.
Definition at line 69 of file colordefs.h.
const int lsseg::BLACK[] = {0, 0, 0} |
3-tuple of int
representing the color black in RGB format.
Definition at line 71 of file colordefs.h.
const int lsseg::GREY[] = {200, 200, 200} |
3-tuple of int
representing the color grey in RGB format.
Definition at line 73 of file colordefs.h.
const int lsseg::BROWN[] = {200, 100, 40} |
3-tuple of int
representing the color brown in RGB format.
Definition at line 75 of file colordefs.h.