OTB  6.7.0
Orfeo Toolbox
otbLabelImageToOGRDataSourceFilter.hxx
Go to the documentation of this file.
1 /*
2  * Copyright (C) 2005-2019 Centre National d'Etudes Spatiales (CNES)
3  *
4  * This file is part of Orfeo Toolbox
5  *
6  * https://www.orfeo-toolbox.org/
7  *
8  * Licensed under the Apache License, Version 2.0 (the "License");
9  * you may not use this file except in compliance with the License.
10  * You may obtain a copy of the License at
11  *
12  * http://www.apache.org/licenses/LICENSE-2.0
13  *
14  * Unless required by applicable law or agreed to in writing, software
15  * distributed under the License is distributed on an "AS IS" BASIS,
16  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
17  * See the License for the specific language governing permissions and
18  * limitations under the License.
19  */
20 
21 #ifndef otbLabelImageToOGRDataSourceFilter_hxx
22 #define otbLabelImageToOGRDataSourceFilter_hxx
23 
25 #include "otbGdalDataTypeBridge.h"
26 
27 //gdal libraries
28 #include "gdal.h"
29 #include "gdal_priv.h"
30 #include "cpl_conv.h"
31 #include "gdal_alg.h"
32 
33 #include "stdint.h" //needed for uintptr_t
34 
35 namespace otb
36 {
37 template <class TInputImage>
39 ::LabelImageToOGRDataSourceFilter() : m_FieldName("DN"), m_Use8Connected(false)
40 {
41  this->SetNumberOfRequiredInputs(2);
42  this->SetNumberOfRequiredInputs(1);
43  this->SetNumberOfRequiredOutputs(1);
44 
45  GDALAllRegister();
46 
47  this->ProcessObject::SetNthOutput(0, this->MakeOutput(0) );
48 }
49 
50 
51 template <class TInputImage>
55 {
56  return static_cast< DataObjectPointer >(OGRDataSourceType::New().GetPointer());
57 }
58 
59 template <class TInputImage>
63 {
64  return static_cast< const OGRDataSourceType *>(
65  this->ProcessObject::GetOutput(0));
66 }
67 
68 template <class TInputImage>
69 void
72 {
73  this->Superclass::SetNthInput(0, const_cast<InputImageType *>(input));
74 }
75 
76 template <class TInputImage>
78 ::InputImageType *
81 {
82  if (this->GetNumberOfInputs() < 1)
83  {
84  return nullptr;
85  }
86 
87  return static_cast<const InputImageType *>(this->Superclass::GetInput(0));
88 }
89 
90 template <class TInputImage>
91 void
94 {
95  this->Superclass::SetNthInput(1, const_cast<InputImageType *>(input));
96 }
97 
98 template <class TInputImage>
100 ::InputImageType *
103 {
104  if (this->GetNumberOfInputs() < 2)
105  {
106  return nullptr;
107  }
108 
109  return static_cast<const InputImageType *>(this->Superclass::GetInput(1));
110 }
111 
112 template <class TInputImage>
113 void
116 {
117  // call the superclass' implementation of this method
118  Superclass::GenerateInputRequestedRegion();
119 
120  // get pointers to the inputs
121  typename InputImageType::Pointer input =
122  const_cast<InputImageType *> (this->GetInput());
123 
124  if ( !input )
125  {
126  return;
127  }
128  // The input is necessarily the largest possible region.
129  input->SetRequestedRegionToLargestPossibleRegion();
130 
131  typename InputImageType::Pointer mask =
132  const_cast<InputImageType *> (this->GetInputMask());
133  if(!mask)
134  {
135  return;
136  }
137  // The input is necessarily the largest possible region.
138  mask->SetRequestedRegionToLargestPossibleRegion();
139 }
140 
141 
142 template <class TInputImage>
143 void
146 {
147  if (this->GetInput()->GetRequestedRegion() != this->GetInput()->GetLargestPossibleRegion())
148  {
149  itkExceptionMacro(<< "Not streamed filter. ERROR : requested region is not the largest possible region.");
150  }
151 
152  SizeType size;
153  unsigned int nbBands = 0;
154  unsigned int bytePerPixel = 0;
155 
156  /* Convert the input image into a GDAL raster needed by GDALPolygonize */
157  size = this->GetInput()->GetLargestPossibleRegion().GetSize();
158  nbBands = this->GetInput()->GetNumberOfComponentsPerPixel();
159  bytePerPixel = sizeof(InputPixelType);
160 
161  // buffer casted in unsigned long cause under Win32 the address
162  // don't begin with 0x, the address in not interpreted as
163  // hexadecimal but alpha numeric value, then the conversion to
164  // integer make us pointing to an non allowed memory block => Crash.
165  std::ostringstream stream;
166  stream << "MEM:::"
167  << "DATAPOINTER=" << (uintptr_t)(this->GetInput()->GetBufferPointer()) << ","
168  << "PIXELS=" << size[0] << ","
169  << "LINES=" << size[1] << ","
170  << "BANDS=" << nbBands << ","
171  << "DATATYPE=" << GDALGetDataTypeName(GdalDataTypeBridge::GetGDALDataType<InputPixelType>()) << ","
172  << "PIXELOFFSET=" << bytePerPixel * nbBands << ","
173  << "LINEOFFSET=" << bytePerPixel * nbBands * size[0] << ","
174  << "BANDOFFSET=" << bytePerPixel;
175 
176  GDALDataset * dataset = static_cast<GDALDataset *> (GDALOpen(stream.str().c_str(), GA_ReadOnly));
177 
178  //Set input Projection ref and Geo transform to the dataset.
179  dataset->SetProjection(this->GetInput()->GetProjectionRef().c_str());
180 
181  unsigned int projSize = this->GetInput()->GetGeoTransform().size();
182  double geoTransform[6];
183 
184  //Set the geo transform of the input image (if any)
185  // Reporting origin and spacing of the buffered region
186  // the spacing is unchanged, the origin is relative to the buffered region
187  IndexType bufferIndexOrigin = this->GetInput()->GetBufferedRegion().GetIndex();
188  OriginType bufferOrigin;
189  this->GetInput()->TransformIndexToPhysicalPoint(bufferIndexOrigin, bufferOrigin);
190  geoTransform[0] = bufferOrigin[0] - 0.5 * this->GetInput()->GetSignedSpacing()[0];
191  geoTransform[3] = bufferOrigin[1] - 0.5 * this->GetInput()->GetSignedSpacing()[1];
192  geoTransform[1] = this->GetInput()->GetSignedSpacing()[0];
193  geoTransform[5] = this->GetInput()->GetSignedSpacing()[1];
194  // FIXME: Here component 1 and 4 should be replaced by the orientation parameters
195  if (projSize == 0)
196  {
197  geoTransform[2] = 0.;
198  geoTransform[4] = 0.;
199  }
200  else
201  {
202  geoTransform[2] = this->GetInput()->GetGeoTransform()[2];
203  geoTransform[4] = this->GetInput()->GetGeoTransform()[4];
204  }
205  dataset->SetGeoTransform(geoTransform);
206 
207  //Create the output layer for GDALPolygonize().
209 
210  OGRLayerType outputLayer = ogrDS->CreateLayer("layer",nullptr,wkbPolygon);
211 
212  OGRFieldDefn field(m_FieldName.c_str(),OFTInteger);
213  outputLayer.CreateField(field, true);
214 
215  //Call GDALPolygonize()
216  char ** options;
217  options = nullptr;
218  char * option[2]= { nullptr , nullptr } ;
219  if (m_Use8Connected == true)
220  {
221  std::string opt("8CONNECTED:8");
222  option[0] = const_cast<char *>(opt.c_str());
223  options=option;
224  }
225 
226  /* Convert the mask input into a GDAL raster needed by GDALPolygonize */
227  typename InputImageType::ConstPointer inputMask = this->GetInputMask();
228  if (!inputMask.IsNull())
229  {
230  size = this->GetInputMask()->GetLargestPossibleRegion().GetSize();
231  nbBands = this->GetInputMask()->GetNumberOfComponentsPerPixel();
232  bytePerPixel = sizeof(InputPixelType);
233  // buffer casted in unsigned long cause under Win32 the address
234  // don't begin with 0x, the address in not interpreted as
235  // hexadecimal but alpha numeric value, then the conversion to
236  // integer make us pointing to an non allowed memory block => Crash.
237  std::ostringstream maskstream;
238  maskstream << "MEM:::"
239  << "DATAPOINTER=" << (uintptr_t)(this->GetInputMask()->GetBufferPointer()) << ","
240  << "PIXELS=" << size[0] << ","
241  << "LINES=" << size[1] << ","
242  << "BANDS=" << nbBands << ","
243  << "DATATYPE=" << GDALGetDataTypeName(GdalDataTypeBridge::GetGDALDataType<InputPixelType>()) << ","
244  << "PIXELOFFSET=" << bytePerPixel * nbBands << ","
245  << "LINEOFFSET=" << bytePerPixel * nbBands * size[0] << ","
246  << "BANDOFFSET=" << bytePerPixel;
247 
248  GDALDataset * maskDataset = static_cast<GDALDataset *> (GDALOpen(maskstream.str().c_str(), GA_ReadOnly));
249 
250  //Set input Projection ref and Geo transform to the dataset.
251  maskDataset->SetProjection(this->GetInputMask()->GetProjectionRef().c_str());
252 
253  projSize = this->GetInputMask()->GetGeoTransform().size();
254 
255  //Set the geo transform of the input mask image (if any)
256  // Reporting origin and spacing of the buffered region
257  // the spacing is unchanged, the origin is relative to the buffered region
258  bufferIndexOrigin = this->GetInputMask()->GetBufferedRegion().GetIndex();
259  this->GetInputMask()->TransformIndexToPhysicalPoint(bufferIndexOrigin, bufferOrigin);
260  geoTransform[0] = bufferOrigin[0] - 0.5 * this->GetInputMask()->GetSignedSpacing()[0];
261  geoTransform[3] = bufferOrigin[1] - 0.5 * this->GetInputMask()->GetSignedSpacing()[1];
262  geoTransform[1] = this->GetInputMask()->GetSignedSpacing()[0];
263  geoTransform[5] = this->GetInputMask()->GetSignedSpacing()[1];
264  // FIXME: Here component 1 and 4 should be replaced by the orientation parameters
265  if (projSize == 0)
266  {
267  geoTransform[2] = 0.;
268  geoTransform[4] = 0.;
269  }
270  else
271  {
272  geoTransform[2] = this->GetInputMask()->GetGeoTransform()[2];
273  geoTransform[4] = this->GetInputMask()->GetGeoTransform()[4];
274  }
275  maskDataset->SetGeoTransform(geoTransform);
276 
277  GDALPolygonize(dataset->GetRasterBand(1), maskDataset->GetRasterBand(1), &outputLayer.ogr(), 0, options, nullptr, nullptr);
278  GDALClose(maskDataset);
279  }
280  else
281  {
282  GDALPolygonize(dataset->GetRasterBand(1), nullptr, &outputLayer.ogr(), 0, options, nullptr, nullptr);
283  }
284 
285  this->SetNthOutput(0,ogrDS);
286 
287  //Clear memory
288  GDALClose(dataset);
289 
290 }
291 
292 
293 } // end namespace otb
294 
295 #endif
Collection of geometric objects.
void CreateField(FieldDefn const &field, bool bApproxOK=true)
Layer of geometric objets.
virtual void SetInput(const InputImageType *input)
static Pointer New()
DataObjectPointerArray::size_type DataObjectPointerArraySizeType
OGRLayer & ogr()
virtual void SetInputMask(const InputImageType *input)
typedef::vcl_size_t uintptr_t
DataObjectPointer MakeOutput(DataObjectPointerArraySizeType idx) override