From 1f11a9d648925335666303378cd2c486fcb26fc1 Mon Sep 17 00:00:00 2001
From: aregimbeau <antoine.regimbeau@c-s.fr>
Date: Wed, 5 Jul 2017 18:14:16 +0200
Subject: [PATCH] First draw of the histogram equalisation working on
black&white pictures and color ones
---
Modules/Feature/Contrast/CMakeLists.txt | 15 ++
.../Contrast/otbContrastEnhancementFilter.cxx | 90 +++++++++++
.../Contrast/otbContrastEnhancementFilter.h | 0
Modules/Filtering/Contrast/CMakeLists.txt | 22 +++
.../include/otbContrastEnhancementFilter.h | 0
Modules/Filtering/Contrast/otb-module.cmake | 38 +++++
.../Filtering/Contrast/test/CMakeLists.txt | 29 ++++
.../test/otbContrastEnhancementFilter.cxx | 141 ++++++++++++++++++
8 files changed, 335 insertions(+)
create mode 100644 Modules/Feature/Contrast/CMakeLists.txt
create mode 100644 Modules/Feature/Contrast/otbContrastEnhancementFilter.cxx
create mode 100644 Modules/Feature/Contrast/otbContrastEnhancementFilter.h
create mode 100644 Modules/Filtering/Contrast/CMakeLists.txt
create mode 100644 Modules/Filtering/Contrast/include/otbContrastEnhancementFilter.h
create mode 100644 Modules/Filtering/Contrast/otb-module.cmake
create mode 100644 Modules/Filtering/Contrast/test/CMakeLists.txt
create mode 100644 Modules/Filtering/Contrast/test/otbContrastEnhancementFilter.cxx
diff --git a/Modules/Feature/Contrast/CMakeLists.txt b/Modules/Feature/Contrast/CMakeLists.txt
new file mode 100644
index 0000000000..87d47ca879
--- /dev/null
+++ b/Modules/Feature/Contrast/CMakeLists.txt
@@ -0,0 +1,15 @@
+cmake_minimum_required (VERSION 2.6)
+project (projchambolle)
+
+find_package(OTB REQUIRED)
+include(${OTB_USE_FILE})
+
+
+set(SOURCES otbContrastEnhancementFilter.cxx)
+
+add_executable(Contrast otbContrastEnhancementFilter.cxx)
+
+target_link_libraries(Contrast ${OTB_LIBRARIES})
+
+
+# message(WARNING "Include dirs : ${OTB_INCLUDE_DIRS}")
\ No newline at end of file
diff --git a/Modules/Feature/Contrast/otbContrastEnhancementFilter.cxx b/Modules/Feature/Contrast/otbContrastEnhancementFilter.cxx
new file mode 100644
index 0000000000..37efcc9775
--- /dev/null
+++ b/Modules/Feature/Contrast/otbContrastEnhancementFilter.cxx
@@ -0,0 +1,90 @@
+#include <iostream>
+#include <math.h>
+#include "otbImage.h"
+#include "otbVectorImage.h"
+#include "otbImageFileReader.h"
+#include "otbImageFileWriter.h"
+#include "itkAdditiveGaussianNoiseImageFilter.h"
+#include "otbStreamingStatisticsImageFilter.h"
+#include "itkImageLinearIteratorWithIndex.h"
+#include "itkMultiplyImageFilter.h"
+#include "itkDivideImageFilter.h"
+#include "itkAbsImageFilter.h"
+#include "itkSubtractImageFilter.h"
+#include "itkAddImageFilter.h"
+#include "itkSquareImageFilter.h"
+#include "itkSqrtImageFilter.h"
+#include "itkMaximumImageFilter.h"
+
+typedef float PixelType;
+typedef itk::Image< PixelType , 2> ImageType;
+typedef itk::AddImageFilter<ImageType, ImageType, ImageType> AddFilterType;
+typedef otb::StreamingStatisticsImageFilter<ImageType> StatFilter;
+
+ImageType::Pointer
+createImagefrom( ImageType::Pointer input )
+{
+ ImageType::Pointer cloneImage( ImageType::New() );
+ cloneImage->SetRegions( input->GetRequestedRegion() );
+ cloneImage->Allocate();
+ cloneImage->SetOrigin( input->GetOrigin() );
+ return cloneImage ;
+
+}
+void
+copyAtoB( const ImageType::Pointer input ,
+ ImageType::Pointer output )
+{
+ itk::ImageRegionConstIterator< ImageType > itinp( input , input->GetRequestedRegion() );
+ itk::ImageRegionIterator< ImageType > itout( output , output->GetRequestedRegion() );
+ itinp.GoToBegin();
+ itout.GoToBegin();
+ while (!itinp.IsAtEnd())
+ {
+ itout.Set( itinp.Get() );
+ ++itout;
+ ++itinp;
+ }
+}
+
+void
+normL( ImageType::Pointer inputX ,
+ ImageType::Pointer inputY ,
+ ImageType::Pointer output )
+{
+ typedef itk::SquareImageFilter< ImageType , ImageType > SquareFilterType;
+ typedef itk::SqrtImageFilter< ImageType , ImageType > SqrtFitlerType;
+ SquareFilterType::Pointer square( SquareFilterType::New() );
+ SqrtFitlerType::Pointer sqrt( SqrtFitlerType::New() );
+ AddFilterType::Pointer add( AddFilterType::New() );
+ square->SetInput( inputX );
+ square->Update();
+ copyAtoB( square->GetOutput() , output );
+ add->SetInput1( output );
+ square->SetInput( inputY );
+ square->Update();
+ add->SetInput2( square->GetOutput() );
+ sqrt->SetInput( add->GetOutput() );
+ sqrt->Update();
+ copyAtoB( sqrt->GetOutput() , output );
+}
+
+int
+main( int argc,
+ char const *argv[] )
+{
+ const char * inputFilename = argv[ 1 ];
+ const char * outputFilename = argv[ 2 ];
+
+
+ typedef otb::ImageFileReader< ImageType > ReaderType;
+
+ ReaderType::Pointer reader( ReaderType::New() );
+ reader->SetFileName( argv[ 1 ] );
+
+ ImageType::Pointer noisyImage( ImageType::New() );
+
+ typedef otb::ImageFileWriter<ImageType> WriterType;
+ WriterType::Pointer writer( WriterType::New());
+ return 0;
+}
\ No newline at end of file
diff --git a/Modules/Feature/Contrast/otbContrastEnhancementFilter.h b/Modules/Feature/Contrast/otbContrastEnhancementFilter.h
new file mode 100644
index 0000000000..e69de29bb2
diff --git a/Modules/Filtering/Contrast/CMakeLists.txt b/Modules/Filtering/Contrast/CMakeLists.txt
new file mode 100644
index 0000000000..45a89f243f
--- /dev/null
+++ b/Modules/Filtering/Contrast/CMakeLists.txt
@@ -0,0 +1,22 @@
+#
+# Copyright (C) 2005-2017 Centre National d'Etudes Spatiales (CNES)
+#
+# This file is part of Orfeo Toolbox
+#
+# https://www.orfeo-toolbox.org/
+#
+# Licensed under the Apache License, Version 2.0 (the "License");
+# you may not use this file except in compliance with the License.
+# You may obtain a copy of the License at
+#
+# http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+#
+
+project(OTBContrast)
+otb_module_impl()
\ No newline at end of file
diff --git a/Modules/Filtering/Contrast/include/otbContrastEnhancementFilter.h b/Modules/Filtering/Contrast/include/otbContrastEnhancementFilter.h
new file mode 100644
index 0000000000..e69de29bb2
diff --git a/Modules/Filtering/Contrast/otb-module.cmake b/Modules/Filtering/Contrast/otb-module.cmake
new file mode 100644
index 0000000000..84da787da6
--- /dev/null
+++ b/Modules/Filtering/Contrast/otb-module.cmake
@@ -0,0 +1,38 @@
+#
+# Copyright (C) 2005-2017 Centre National d'Etudes Spatiales (CNES)
+#
+# This file is part of Orfeo Toolbox
+#
+# https://www.orfeo-toolbox.org/
+#
+# Licensed under the Apache License, Version 2.0 (the "License");
+# you may not use this file except in compliance with the License.
+# You may obtain a copy of the License at
+#
+# http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+#
+
+set(DOCUMENTATION "work in progress")
+
+otb_module(OTBContrast
+ DEPENDS
+ OTBITK
+ OTBCommon
+ OTBImageBase
+
+
+ TEST_DEPENDS
+ OTBTestKernel
+ OTBImageIO
+ OTBImageBase
+
+
+ DESCRIPTION
+ "${DOCUMENTATION}"
+)
diff --git a/Modules/Filtering/Contrast/test/CMakeLists.txt b/Modules/Filtering/Contrast/test/CMakeLists.txt
new file mode 100644
index 0000000000..e6d3bf6920
--- /dev/null
+++ b/Modules/Filtering/Contrast/test/CMakeLists.txt
@@ -0,0 +1,29 @@
+#
+# Copyright (C) 2005-2017 Centre National d'Etudes Spatiales (CNES)
+#
+# This file is part of Orfeo Toolbox
+#
+# https://www.orfeo-toolbox.org/
+#
+# Licensed under the Apache License, Version 2.0 (the "License");
+# you may not use this file except in compliance with the License.
+# You may obtain a copy of the License at
+#
+# http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+#
+
+otb_module_test()
+
+set(OTBContrastTests
+otbContrastEnhancementFilter.cxx
+)
+
+add_executable(otbContrastTestDriver ${OTBContrastTests})
+target_link_libraries(otbContrastTestDriver ${OTBContrast-Test_LIBRARIES})
+otb_module_target_label(otbContrastTestDriver)
diff --git a/Modules/Filtering/Contrast/test/otbContrastEnhancementFilter.cxx b/Modules/Filtering/Contrast/test/otbContrastEnhancementFilter.cxx
new file mode 100644
index 0000000000..2ac67c20b9
--- /dev/null
+++ b/Modules/Filtering/Contrast/test/otbContrastEnhancementFilter.cxx
@@ -0,0 +1,141 @@
+#include <iostream>
+#include <math.h>
+#include <vector>
+#include "otbImage.h"
+#include "otbVectorImage.h"
+#include "otbImageFileReader.h"
+#include "otbImageFileWriter.h"
+#include "itkImageLinearIteratorWithIndex.h"
+#include "itkImageToHistogramFilter.h"
+#include "itkVectorIndexSelectionCastImageFilter.h"
+#include "itkMultiplyImageFilter.h"
+
+
+typedef int PixelType;
+typedef itk::Image< PixelType , 2 > ImageType;
+typedef itk::VectorImage< PixelType , 2 > VectorImageType;
+
+void
+computeLuminance( ImageType::Pointer luminance ,
+ VectorImageType::Pointer input )
+{
+ itk::ImageRegionIterator< VectorImageType > itinp( input , input->GetRequestedRegion() );
+ itk::ImageRegionIterator< ImageType > itlum( luminance , luminance->GetRequestedRegion() );
+ itinp.GoToBegin();
+ itlum.GoToBegin();
+ while ( !itinp.IsAtEnd() )
+ {
+ // 0.3R 0.6G 0.1B
+ itlum.Set( round(0.3 * itinp.Get()[0] + 0.6 * itinp.Get()[1] + 0.1 * itinp.Get()[2]) );
+ ++itinp;
+ ++itlum;
+ }
+}
+
+int
+main( int argc,
+ char const *argv[] )
+{
+ // const char * inputfilename = argv[ 1 ];
+ // const char * outputfilename = argv[ 2 ];
+
+
+ typedef otb::ImageFileReader< VectorImageType > ReaderType;
+ typedef itk::Statistics::ImageToHistogramFilter< ImageType > HistogramFilterType;
+ typedef HistogramFilterType::HistogramType HistogramType;
+ typedef itk::VectorIndexSelectionCastImageFilter< VectorImageType , ImageType > VectorToImageFilterType;
+
+
+ ReaderType::Pointer reader( ReaderType::New() );
+ reader->SetFileName( argv[ 1 ] );
+ reader->Update();
+ VectorImageType::Pointer input = reader->GetOutput();
+ ImageType::Pointer lum ( ImageType::New() );
+ lum->SetRegions( input->GetRequestedRegion() );
+ lum->Allocate();
+ lum->SetOrigin( input->GetOrigin() );
+ VectorToImageFilterType::Pointer vectorToImageFilter ( VectorToImageFilterType::New() );
+ if ( input->GetVectorLength() == 3 )
+ {
+ computeLuminance(lum, input);
+ }
+ else
+ {
+ // default is a black and white image, if length>3 and is not equal to one
+ // other case : how to apply the equalization :
+ // - channel per channel?
+ // - with a general luminance?
+ vectorToImageFilter->SetInput( input );
+ vectorToImageFilter->SetIndex( 0 );
+ vectorToImageFilter->Update();
+ lum = vectorToImageFilter->GetOutput();
+ }
+
+
+ HistogramFilterType::Pointer histofilter( HistogramFilterType::New() );
+ histofilter->SetInput( lum );
+ histofilter->Update();
+ // Compute the histogram of the bw image
+ HistogramType * histo = histofilter->GetOutput();
+
+ // Define propertise of the target histogram : nomber of bin, height per bin
+ ImageType::SizeType lumSize = lum->GetLargestPossibleRegion().GetSize();
+ int npx = lumSize[0] * lumSize[1];
+ int nbin = histo->GetSize()[0];
+ float height = npx / nbin;
+
+ // The transfer function will be a look up table
+ std::vector< int > lut;
+ float nbOri, nbNew ;
+ int countnew, countold;
+ countnew = 0;
+ countold = 0;
+ nbNew = height;
+ HistogramType::Iterator it = histo->Begin();
+ nbOri = it.GetFrequency();
+ while ( it != histo->End() )
+ {
+ if (nbOri> nbNew)
+ {
+ nbNew += height;
+ // index[0] += 1;
+ ++countnew;
+ }
+ else
+ {
+ ++it;
+ nbOri += it.GetFrequency();
+ ++countold;
+ lut.push_back(countnew);
+ }
+ }
+
+ // Here we compute a gain image corresponding to the associated gain of the equalization
+ typedef itk::Image< float , 2 > ImageGainType;
+ ImageGainType::Pointer gainImage ( ImageGainType::New() );
+ gainImage->SetRegions( lum->GetRequestedRegion() );
+ gainImage->Allocate();
+ gainImage->SetOrigin( lum->GetOrigin() );
+ itk::ImageRegionIterator< ImageType > itlum( lum , lum->GetRequestedRegion() );
+ itk::ImageRegionIterator< ImageGainType > itgain( gainImage , gainImage->GetRequestedRegion() );
+ itlum.GoToBegin();
+ itgain.GoToBegin();
+ while( !itlum.IsAtEnd() )
+ {
+ itgain.Set( static_cast<float>(lut[itlum.Get()]) / static_cast<float>( itlum.Get() ) ) ;
+ ++itlum;
+ ++itgain;
+ }
+
+ typedef itk::MultiplyImageFilter< VectorImageType, ImageGainType, VectorImageType > MultiplyImageFilterType;
+ MultiplyImageFilterType::Pointer gainMultiplyer ( MultiplyImageFilterType::New() );
+ gainMultiplyer->SetInput1( input );
+ gainMultiplyer->SetInput2( gainImage );
+ gainMultiplyer->Update();
+ typedef otb::ImageFileWriter<VectorImageType> writertype;
+ writertype::Pointer writer( writertype::New());
+ writer->SetFileName( argv[2] );
+ writer->SetInput( gainMultiplyer->GetOutput() );
+ writer->Update();
+ return 0;
+}
--
GitLab