Commit b89b9b6c authored by Fabien Freling's avatar Fabien Freling
Browse files

Add source for dicom2dump tool.

	* fabien/TODO: Simple TODO reminder.
	* fabien/bin/Makefile: Makefile to compile dicom binaries.
	* fabien/bin/dicom2dump.cc: Binary to convert dicom to dump.
	* fabien/bin/dicom_mask.cc: Binary to extract simple mask
	  from dicom images.
	* fabien/igr/Makefile: Update.
	* fabien/igr/launch2d.sh: Update.
	* fabien/igr/launch3d.sh: Update.
	* fabien/igr/seg2d.cc: Update.
	* fabien/igr/seg3d.cc: Update.
	* fabien/igr/seg_vol_irm.hh: Implement different threshold
	  techniques: double and deviation.

git-svn-id: https://svn.lrde.epita.fr/svn/oln/trunk@3482 4aad255d-cdde-0310-9447-f3009e2ae8c0
parent 2cca7666
2009-03-05 Fabien Freling <fabien.freling@lrde.epita.fr>
Add source for dicom2dump tool.
* fabien/TODO: Simple TODO reminder.
* fabien/bin/Makefile: Makefile to compile dicom binaries.
* fabien/bin/dicom2dump.cc: Binary to convert dicom to dump.
* fabien/bin/dicom_mask.cc: Binary to extract simple mask
from dicom images.
* fabien/igr/Makefile: Update.
* fabien/igr/launch2d.sh: Update.
* fabien/igr/launch3d.sh: Update.
* fabien/igr/seg2d.cc: Update.
* fabien/igr/seg3d.cc: Update.
* fabien/igr/seg_vol_irm.hh: Implement different threshold
techniques: double and deviation.
2009-03-05 Thierry Geraud <thierry.geraud@lrde.epita.fr>
Fix fastest distance computation.
......
......@@ -9,5 +9,14 @@
[X] Create binaries
[X] Create scripts shell
[ ] Check standard deviation
[X] Generate color images for regions
[X] Generate histograms (normal, bg, obj, p_bg, p_obj)
[ ] Test processing chain on US
[X] Create README file for special images
[ ] Implement watershed
[X] Create tool for dicom mask
[X] Check conversion for mask (everything except int_u12(0))
[ ] Extract ROI (projection or fill holes)
[ ] Create routine for region colors
[ ] Create dicom2dump (2d, 3d, int_u8, int_u12)
[ ] After threshold, take biggest object and then erode, dilate
GDCM_SRC_DIR = /Users/HiSoKa/Downloads/gdcm-2.0.10
GDCM_BIN_DIR = /Users/HiSoKa/Downloads/gdcmbin
DICOM_INC = -I${GDCM_SRC_DIR}/Source/Common/ \
-I${GDCM_BIN_DIR}/Source/Common/ \
-I${GDCM_SRC_DIR}/Source/DataDictionary/ \
-I${GDCM_SRC_DIR}/Source/MediaStorageAndFileFormat/ \
-I${GDCM_SRC_DIR}/Source/DataStructureAndEncodingDefinition/
# "-framework CoreFoundation" is a Mac OS X specific flag
DICOM_LIB = -L${GDCM_BIN_DIR}/bin \
-lgdcmCommon -lgdcmDICT -lgdcmDSED -lgdcmIOD -lgdcmMSFF -lgdcmexpat -lgdcmjpeg12 -lgdcmjpeg16 -lgdcmjpeg8 -lgdcmopenjpeg -lgdcmuuid -lgdcmzlib \
-framework CoreFoundation
CXXFLAGS = -DNDEBUG -O1
all: mask dump
mask: dicom_mask.cc
g++ -I../../../ ${DICOM_INC} ${DICOM_LIB} ${CXXFLAGS} $^ -o dicom_mask
dump: dicom2dump.cc
g++ -I../../../ ${DICOM_INC} ${DICOM_LIB} ${CXXFLAGS} $^ -o dicom2dump
#include <mln/core/concept/image.hh>
#include <mln/core/image/image3d.hh>
#include <mln/value/int_u12.hh>
#include <mln/io/dicom/load.hh>
#include <mln/io/dump/save.hh>
int usage(char* argv[])
{
std::cerr << "usage: " << argv[0] << " input.dcm output.dump" << std::endl;
std::cerr << "\t work for 3D images encoded in int_u12" << std::endl;
return 1;
}
int main(int argc, char* argv[])
{
using namespace mln;
using value::int_u12;
if (argc != 3)
return usage(argv);
image3d<int_u12> ima;
io::dicom::load(ima, argv[1]);
io::dump::save(ima, argv[2]);
return 0;
}
#include <mln/core/concept/image.hh>
#include <mln/core/image/image2d.hh>
#include <mln/value/int_u8.hh>
#include <mln/io/dicom/load.hh>
#include <mln/io/pbm/save.hh>
#include <mln/literal/colors.hh>
#include <mln/level/transform.hh>
#include <mln/fun/v2b/threshold.hh>
#include <mln/data/fill.hh>
#include <mln/pw/all.hh>
void usage(char* argv[])
{
std::cerr << "usage: " << argv[0] << " input.dcm output.pbm" << std::endl;
abort();
}
int main(int argc, char* argv[])
{
using namespace mln;
using value::int_u8;
if (argc != 3)
usage(argv);
image2d<int_u8> input;
io::dicom::load(input, argv[1]);
image2d<bool> ima = level::transform(input, fun::v2b::threshold<int_u8>(1));
io::pbm::save(ima, argv[2]);
}
......@@ -16,5 +16,11 @@ all: 2d 3d
3d: seg_vol_irm.hh seg3d.cc
g++ -I../../../ ${DICOM_INC} ${DICOM_LIB} ${CXXFLAGS} seg3d.cc -o seg3d
grad: grad_clo_and_wshd.cc
g++ -I../../../ ${DICOM_INC} ${DICOM_LIB} ${CXXFLAGS} $^ -o grad_clo
wst: wst_rag.cc
g++ -I../../../ ${DICOM_INC} ${DICOM_LIB} ${CXXFLAGS} $^ -o wst_rag
clean:
rm -rf *.dump *.p?m *.plot *.log *.csv
......@@ -3,8 +3,9 @@
process_file ()
{
./seg2d $1
if [ -f result.pbm ]; then
mv result.pbm results/${2}_06_result.pbm
if [ -f result_double.pbm ]; then
mv result_double.pbm results/${2}_06_result_double.pbm
mv result_deviation.pbm results/${2}_07_result_deviation.pbm
fi
if [ -f bg.pbm ]; then
......@@ -14,10 +15,15 @@ process_file ()
mv bg_closed.pbm results/${2}_02_bg_closed.pbm
mv obj_closed.pbm results/${2}_04_obj_closed.pbm
mv regions_color.ppm results/${2}_05_regions_color.ppm
if [ -f regions_color.ppm ]; then
mv regions_color.ppm results/${2}_05_regions_color.ppm
fi
mv bg_histo_norm.plot results/${2}_bg.plot
mv obj_histo_norm.plot results/${2}_obj.plot
mv histo.plot results/${2}.plot
mv bg_histo.plot results/${2}_bg.plot
mv obj_histo.plot results/${2}_obj.plot
mv bg_histo_norm.plot results/${2}_p_bg.plot
mv obj_histo_norm.plot results/${2}_p_obj.plot
fi
}
......
......@@ -3,7 +3,8 @@
process_file ()
{
./seg3d $1
../bin/dump2pbm result.dump results/${2}_06_result.pbm
../bin/dump2pbm result_double.dump results/${2}_06_result_double.pbm
../bin/dump2pbm result_deviation.dump results/${2}_07_result_deviation.pbm
../bin/dump2pbm bg.dump results/${2}_01_bg.pbm
../bin/dump2pbm obj.dump results/${2}_03_obj.pbm
......@@ -11,12 +12,15 @@ process_file ()
../bin/dump2pbm bg_closed.dump results/${2}_02_bg_closed.pbm
../bin/dump2pbm obj_closed.dump results/${2}_04_obj_closed.pbm
#../bin/dump2ppm regions_color.dump results/${2}_05_colors.ppm
../bin/dump2ppm regions_color.dump results/${2}_05_regions_colors.ppm
rm *.dump
mv bg_histo_norm.plot results/${2}_bg.plot
mv obj_histo_norm.plot results/${2}_obj.plot
mv histo.plot results/${2}.plot
mv bg_histo.plot results/${2}_bg.plot
mv obj_histo.plot results/${2}_obj.plot
mv bg_histo_norm.plot results/${2}_p_bg.plot
mv obj_histo_norm.plot results/${2}_p_obj.plot
}
process_file "/Users/HiSoKa/Work/IGR/souris18/irm/IM_0052.dcm" "52"
......
......@@ -34,7 +34,22 @@ int main(int argc, char* argv[])
label_16 nlabels;
image2d<int_u12> ima;
io::dicom::load(ima, argv[1]);
io::pbm::save(igr_seg(ima, c4(), nlabels), "result.pbm");
image2d<bool> ima_thres;
initialize(ima_thres, ima);
data::fill(ima_thres, false);
unsigned threshold_value = find_threshold_value(ima, c4());
std::cout << "double threshold value = " << threshold_value << std::endl;
data::fill((ima_thres | pw::value(ima) < pw::cst(threshold_value)).rw(), true);
io::pbm::save(ima_thres, "result_double.pbm");
data::fill(ima_thres, false);
threshold_value = find_threshold_mean(ima, c4());
std::cout << " deviation threshold value = " << threshold_value << std::endl;
data::fill((ima_thres | pw::value(ima) < pw::cst(threshold_value)).rw(), true);
io::pbm::save(ima_thres, "result_deviation.pbm");
return 0;
}
......@@ -36,7 +36,22 @@ int main(int argc, char* argv[])
label_16 nlabels;
image3d<int_u12> ima;
io::dicom::load(ima, argv[1]);
io::dump::save(igr_seg(ima, c6(), nlabels), "result.dump");
image3d<bool> ima_thres;
initialize(ima_thres, ima);
data::fill(ima_thres, false);
unsigned threshold_value = find_threshold_value(ima, c6());
std::cout << "double threshold value = " << threshold_value << std::endl;
data::fill((ima_thres | pw::value(ima) < pw::cst(threshold_value)).rw(), true);
io::dump::save(ima_thres, "result_double.dump");
data::fill(ima_thres, false);
threshold_value = find_threshold_mean(ima, c6());
std::cout << " deviation threshold value = " << threshold_value << std::endl;
data::fill((ima_thres | pw::value(ima) < pw::cst(threshold_value)).rw(), true);
io::dump::save(ima_thres, "result_deviation.dump");
return 0;
}
......@@ -39,6 +39,7 @@
#include <mln/value/int_u12.hh>
#include <mln/value/int_u16.hh>
#include <mln/value/rgb8.hh>
#include <mln/value/label_16.hh>
#include <mln/io/pgm/all.hh>
#include <mln/io/pbm/all.hh>
......@@ -56,6 +57,7 @@
#include <mln/labeling/blobs.hh>
#include <mln/labeling/compute.hh>
#include <mln/labeling/foreground.hh>
#include <mln/labeling/fill_holes.hh>
#include <mln/labeling/n_max.hh>
......@@ -85,7 +87,6 @@
#include <mln/pw/all.hh>
#include <mln/morpho/elementary/gradient_internal.hh>
#include <mln/morpho/closing.hh>
#include <mln/morpho/dilation.hh>
#include <mln/morpho/erosion.hh>
......@@ -149,12 +150,30 @@ close_threshold(const Image<I>& input, metal::int_<3>, int dil, int ero)
template <typename I>
void
save_regions_color(const Image<I>& ima, metal::int_<2>)
{
io::ppm::save(ima, "regions_color.ppm");
}
template <typename I>
void
save_regions_color(const Image<I>& ima, metal::int_<3>)
{
io::dump::save(ima, "regions_color.dump");
}
template <typename I, typename N>
unsigned
find_threshold_value(const Image<I>& input, const Neighborhood<N>& nbh)
{
int bg_thres = 30;
int obj_thres = 15;
int bg_thres = 20;
int obj_thres = 10;
mln_ch_value(I, bool) ima_bg;
initialize(ima_bg, input);
......@@ -166,15 +185,16 @@ find_threshold_value(const Image<I>& input, const Neighborhood<N>& nbh)
unsigned result = 0;
histo::array<mln_value(I)> arr_histo = histo::compute(input);
// We remove the 0 value because it is not part of the image.
histo::array<mln_value(I)> arr_histo = histo::compute(input | pw::value(input) != 0);
image1d<unsigned> ima_histo;
convert::from_to(arr_histo, ima_histo);
// We remove the 0 value because it is not part of the image.
ima_histo(point1d(0)) = 0;
std::ofstream fout("histo.plot");
fout << "0 0" << std::endl;
for (unsigned int i = 1; i < ima_histo.nelements(); ++i)
{
fout << i << " " << ima_histo(point1d(i)) << std::endl;
ima_histo(point1d(i)) += ima_histo(point1d(i - 1));
}
accu::max<unsigned> max_accu;
......@@ -207,26 +227,41 @@ find_threshold_value(const Image<I>& input, const Neighborhood<N>& nbh)
io::dump::save(ima_obj, "obj.dump");
}
ima_bg = close_threshold(ima_bg, metal::int_<I::site::dim>(), 3, 5); // 5, 7?
ima_obj = close_threshold(ima_obj, metal::int_<I::site::dim>(), 3, 5);
ima_bg = close_threshold(ima_bg, metal::int_<I::site::dim>(), 9, 15);
ima_obj = close_threshold(ima_obj, metal::int_<I::site::dim>(), 9, 11);
// Debug output images
mln_ch_value(I, rgb8) out = level::convert(rgb8(), level::stretch(int_u8(), input));
data::fill((out | pw::value(morpho::elementary::gradient_internal(ima_bg, nbh)) == true).rw(), literal::red);
data::fill((out | pw::value(morpho::elementary::gradient_internal(ima_obj, nbh)) == true).rw(), literal::green);
save_regions_color(out, metal::int_<I::site::dim>());
if (I::site::dim == 2)
{
io::pbm::save(ima_bg, "bg_closed.pbm");
io::pbm::save(ima_obj, "obj_closed.pbm");
io::ppm::save(out, "regions_color.ppm");
}
if (I::site::dim == 3)
{
io::dump::save(ima_bg, "bg_closed.dump");
io::dump::save(ima_obj, "obj_closed.dump");
io::dump::save(out, "regions_color.dump");
}
// Labeling
/*label_16 nlabels = 0;
mln_ch_value(I, label_16) bg_labels = labeling::foreground(ima_bg, nbh, nlabels);
mln_ch_value(I, label_16) obj_labels = labeling::foreground(ima_obj, nbh, nlabels);
accu::count<int_u8> a_;
util::array<unsigned> arr_label = labeling::compute(a_, ima_bg, bg_labels, nlabels);
util::array<label_16> arr_big = labeling::n_max<label_16>(arr_label, 1);
data::fill((ima_bg | (pw::value(bg_labels) != pw::cst(arr_big[1]))).rw(), false);
arr_label = labeling::compute(a_, ima_obj, obj_labels, nlabels);
arr_big = labeling::n_max<label_16>(arr_label, 1);
data::fill((ima_obj | (pw::value(obj_labels) != pw::cst(arr_big[1]))).rw(), false);*/
// Histo
histo::array<mln_value(I)> bg_histo = histo::compute(input | pw::value(ima_bg) == true);
histo::array<mln_value(I)> obj_histo = histo::compute(input | pw::value(ima_obj) == true);
......@@ -236,23 +271,27 @@ find_threshold_value(const Image<I>& input, const Neighborhood<N>& nbh)
ima_bg_histo(point1d(0)) = 0;
unsigned bg_sum = level::compute(sum_accu, ima_bg_histo);
// We remove the 0 value because it is not part of the image.
std::ofstream fout_bg("bg_histo_norm.plot");
std::ofstream fout_bg("bg_histo.plot");
std::ofstream fout_p_bg("bg_histo_norm.plot");
for (unsigned int i = 0; i < ima_bg_histo.nelements(); ++i)
{
fout_bg << i << " " << ima_bg_histo(point1d(i)) << std::endl;
ima_bg_histo(point1d(i)) *= 10000;
ima_bg_histo(point1d(i)) /= bg_sum;
fout_bg << i << " " << ima_bg_histo(point1d(i)) << std::endl;
fout_p_bg << i << " " << ima_bg_histo(point1d(i)) << std::endl;
}
image1d<unsigned> ima_obj_histo;
convert::from_to(obj_histo, ima_obj_histo);
unsigned obj_sum = level::compute(sum_accu, ima_obj_histo);
std::ofstream fout_obj("obj_histo_norm.plot");
std::ofstream fout_obj("obj_histo.plot");
std::ofstream fout_p_obj("obj_histo_norm.plot");
for (unsigned int i = 0; i < ima_obj_histo.nelements(); ++i)
{
fout_obj << i << " " << ima_obj_histo(point1d(i)) << std::endl;
ima_obj_histo(point1d(i)) *= 10000;
ima_obj_histo(point1d(i)) /= obj_sum;
fout_obj << i << " " << ima_obj_histo(point1d(i)) << std::endl;
fout_p_obj << i << " " << ima_obj_histo(point1d(i)) << std::endl;
}
// Search for the index with the min distance between histogrammes.
......@@ -275,7 +314,7 @@ template <typename I, typename N>
unsigned
find_threshold_mean(const Image<I>& input, const Neighborhood<N>& nbh)
{
unsigned result = 0;
unsigned coef = 1;
accu::mean<unsigned> mean_accu;
unsigned mean = level::compute(mean_accu, (input | (pw::value(input) != 0)));
......@@ -283,13 +322,13 @@ find_threshold_mean(const Image<I>& input, const Neighborhood<N>& nbh)
accu::stat::deviation<unsigned, unsigned, float> dev_accu(mean);
float deviation = level::compute(dev_accu, (input | pw::value(input) != 0));
std::cout << "mean = " << mean << std::endl << "deviation = " << deviation << std::endl;
return floor(deviation);
std::cout << "[mean = " << mean << " | deviation = " << deviation << "]";
return floor(mean + coef * deviation);
}
template <typename I, typename N, typename L>
/*template <typename I, typename N, typename L>
mln_ch_value(I, bool)
igr_seg(const Image<I>& input_, const Neighborhood<N>& nbh_, L& nlabels)
{
......@@ -301,9 +340,10 @@ igr_seg(const Image<I>& input_, const Neighborhood<N>& nbh_, L& nlabels)
mln_ch_value(I, bool) ima_thres;
initialize(ima_thres, input);
data::fill(ima_thres, false);
//unsigned threshold_value = find_threshold_value(input, nbh);
std::cout << "double threshold value = " << find_threshold_value(input, nbh) << std::endl;
unsigned threshold_value = find_threshold_mean(input, nbh);
std::cout << " deviation threshold value = " << threshold_value << std::endl;
data::fill((ima_thres | pw::value(input) < pw::cst(threshold_value)).rw(), true);
return ima_thres;
}
}*/
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment