Orfeo Toolbox  4.0
otbStreamingImageToOGRLayerSegmentationFilter.txx
Go to the documentation of this file.
1 /*=========================================================================
2 
3  Program: ORFEO Toolbox
4  Language: C++
5  Date: $Date$
6  Version: $Revision$
7 
8 
9  Copyright (c) Centre National d'Etudes Spatiales. All rights reserved.
10  See OTBCopyright.txt for details.
11 
12  Some parts of this code are derived from ITK. See ITKCopyright.txt
13  for details.
14 
15 
16  This software is distributed WITHOUT ANY WARRANTY; without even
17  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
18  PURPOSE. See the above copyright notices for more information.
19 
20 =========================================================================*/
21 #ifndef __otbStreamingImageToOGRLayerSegmentationFilter_txx
22 #define __otbStreamingImageToOGRLayerSegmentationFilter_txx
23 
25 
27 #include "itkAffineTransform.h"
28 
29 #include "itkTimeProbe.h"
30 #include "otbMacro.h"
31 #include <cassert>
32 
33 namespace otb
34 {
35 
36 template <class TImageType, class TSegmentationFilter>
38 ::PersistentImageToOGRLayerSegmentationFilter() : m_TileMaxLabel(0), m_StartLabel(0), m_SegmentationFilter(), m_FieldName("DN"), m_Use8Connected(false),
39  m_FilterSmallObject(false), m_MinimumObjectSize(1),m_Simplify(false), m_SimplificationTolerance(0.3)
40 {
41  this->SetNumberOfRequiredInputs(2);
42  this->SetNumberOfRequiredInputs(1);
43  m_SegmentationFilter = SegmentationFilterType::New();
44  m_TileNumber = 1;
45 }
46 
47 template <class TImageType, class TSegmentationFilter>
50 {
51 }
52 
53 template <class TImageType, class TSegmentationFilter>
54 void
57 {
58  this->itk::ProcessObject::SetNthInput(1, const_cast<LabelImageType *>(mask));
59 }
60 
61 template <class TImageType, class TSegmentationFilter>
63 ::LabelImageType *
66 {
67  return static_cast<const LabelImageType *>(this->itk::ProcessObject::GetInput(1));
68 }
69 
70 
71 template <class TImageType, class TSegmentationFilter>
75 {
76  otbMsgDebugMacro(<< "tile number : " << m_TileNumber);
77  ++m_TileNumber;
78  itk::TimeProbe tileChrono;
79  tileChrono.Start();
80 
81  // Apply an ExtractImageFilter to avoid problems with filters asking for the LargestPossibleRegion
82  typedef itk::ExtractImageFilter<InputImageType, InputImageType> ExtractImageFilterType;
83  typename ExtractImageFilterType::Pointer extract = ExtractImageFilterType::New();
84  extract->SetInput( this->GetInput() );
85  extract->SetExtractionRegion( this->GetInput()->GetRequestedRegion() );
86  extract->Update();
87 
88  // WARNING: itk::ExtractImageFilter does not copy the MetadataDictionnary
89  extract->GetOutput()->SetMetaDataDictionary(this->GetInput()->GetMetaDataDictionary());
90 
91  const unsigned int labelImageIndex = LabeledOutputAccessor<SegmentationFilterType>::LabeledOutputIndex;
92 
93  typename LabelImageToOGRDataSourceFilterType::Pointer labelImageToOGRDataFilter =
94  LabelImageToOGRDataSourceFilterType::New();
95 
96  itk::TimeProbe chrono1;
97  chrono1.Start();
98  m_SegmentationFilter->SetInput(extract->GetOutput());
99  m_SegmentationFilter->UpdateLargestPossibleRegion();
100 
101  chrono1.Stop();
102  otbMsgDebugMacro(<<"segmentation took " << chrono1.GetTotal() << " sec");
103 
104  itk::TimeProbe chrono2;
105  chrono2.Start();
106  typename LabelImageType::ConstPointer inputMask = this->GetInputMask();
107  if (!inputMask.IsNull())
108  {
109  // Apply an ExtractImageFilter to avoid problems with filters asking for the LargestPossibleRegion
110  typedef itk::ExtractImageFilter<LabelImageType, LabelImageType> ExtractLabelImageFilterType;
111  typename ExtractLabelImageFilterType::Pointer maskExtract = ExtractLabelImageFilterType::New();
112  maskExtract->SetInput( this->GetInputMask() );
113  maskExtract->SetExtractionRegion( this->GetInput()->GetRequestedRegion() );
114  maskExtract->Update();
115  // WARNING: itk::ExtractImageFilter does not copy the MetadataDictionnary
116  maskExtract->GetOutput()->SetMetaDataDictionary(this->GetInputMask()->GetMetaDataDictionary());
117 
118  labelImageToOGRDataFilter->SetInputMask(maskExtract->GetOutput());
119  }
120 
121  labelImageToOGRDataFilter->SetInput(dynamic_cast<LabelImageType *>(m_SegmentationFilter->GetOutputs().at(labelImageIndex).GetPointer()));
122  labelImageToOGRDataFilter->SetFieldName(m_FieldName);
123  labelImageToOGRDataFilter->SetUse8Connected(m_Use8Connected);
124  labelImageToOGRDataFilter->Update();
125 
126  chrono2.Stop();
127  otbMsgDebugMacro(<< "vectorization took " << chrono2.GetTotal() << " sec");
128 
129  //Relabeling & simplication of geometries & filtering small objects
130  itk::TimeProbe chrono3;
131  chrono3.Start();
132  OGRDataSourcePointerType tmpDS = const_cast<OGRDataSourceType *>(labelImageToOGRDataFilter->GetOutput());
133  OGRLayerType tmpLayer = tmpDS->GetLayer(0);
134 
135  const typename InputImageType::SpacingType inSpacing = this->GetInput()->GetSpacing();
136  const double tol = m_SimplificationTolerance * std::max(vcl_abs(inSpacing[0]),vcl_abs(inSpacing[1]));
137 
138  typename OGRLayerType::iterator featIt = tmpLayer.begin();
139  for(featIt = tmpLayer.begin(); featIt!=tmpLayer.end(); ++featIt)
140  {
141  ogr::Field field = (*featIt)[0];
142  // field.Unset();
143  field.SetValue(m_TileMaxLabel);
144  m_TileMaxLabel++;
145 
146  //Simplify the geometry
147  if (m_Simplify)
148  {
149  const OGRGeometry * geom = (*featIt).GetGeometry();
150  assert(geom && "geometry is NULL ! Can't simplify it.");
151 
152  (*featIt).SetGeometryDirectly(ogr::Simplify(*geom,tol));
153  }
154  //Need to rewrite the feature otherwise changes are not considered.
155  tmpLayer.SetFeature(*featIt);
156 
157  //Filter small objects.
158  if(m_FilterSmallObject)
159  {
160  double area = static_cast<const OGRPolygon *>((*featIt).GetGeometry())->get_Area();
161  //convert into pixel coordinates
162  typename InputImageType::SpacingType spacing = this->GetInput()->GetSpacing();
163  double pixelsArea = area / (vcl_abs(spacing[0]*spacing[1]));
164  otbMsgDebugMacro(<<"DN = "<<field.GetValue<int>()<<", area = "<<pixelsArea);
165  if( pixelsArea < m_MinimumObjectSize )
166  {
167  tmpLayer.DeleteFeature((*featIt).GetFID());
168  }
169  }
170  }
171  chrono3.Stop();
172  otbMsgDebugMacro(<< "relabeling, filtering small objects and simplifying geometries took " << chrono3.GetTotal() << " sec");
173 
174  return tmpDS;
175 }
176 
177 
178 } // end namespace otb
179 #endif

Generated at Sat Mar 8 2014 16:19:51 for Orfeo Toolbox with doxygen 1.8.3.1