Orfeo Toolbox  4.0
otbVectorDataToLabelImageFilter.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 
13  This software is distributed WITHOUT ANY WARRANTY; without even
14  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
15  PURPOSE. See the above copyright notices for more information.
16 
17 =========================================================================*/
18 #ifndef __otbVectorDataToLabelImageFilter_txx
19 #define __otbVectorDataToLabelImageFilter_txx
20 
22 #include "otbOGRIOHelper.h"
23 #include "otbGdalDataTypeBridge.h"
24 #include "otbMacro.h"
25 
26 #include "gdal_alg.h"
27 #include "ogr_srs_api.h"
28 
29 namespace otb
30 {
31 template<class TVectorData, class TOutputImage>
34  : m_OGRDataSourcePointer(0),
35  m_BurnAttribute("FID")
36 {
37  this->SetNumberOfRequiredInputs(1);
38 
39  // Output parameters initialization
40  m_OutputSpacing.Fill(1.0);
41  m_OutputSize.Fill(0);
42  m_OutputStartIndex.Fill(0);
43 
44  // This filter is intended to work with otb::Image
45  m_BandsToBurn.push_back(1);
46 
47  // Default burn value if no burnAttribute availabke
48  m_DefaultBurnValue = 1.;
49 }
50 
51 template<class TVectorData, class TOutputImage>
52 void
55 {
57 }
58 
59 template <class TVectorData, class TOutputImage>
62 ::GetInput(unsigned int idx)
63 {
64  return static_cast<const TVectorData *>
65  (this->itk::ProcessObject::GetInput(idx));
66 }
67 
68 template <class TVectorData, class TOutputImage>
69 void
72 {
73  if (this->m_OutputSpacing != spacing)
74  {
75  this->m_OutputSpacing = spacing;
76  this->Modified();
77  }
78 }
79 
80 template <class TVectorData, class TOutputImage>
81 void
83 ::SetOutputSpacing(const double spacing[2])
84 {
85  OutputSpacingType s(spacing);
86  this->SetOutputSpacing(s);
87 }
88 
89 template <class TVectorData, class TOutputImage>
90 void
92 ::SetOutputSpacing(const float spacing[2])
93 {
94  itk::Vector<float, 2> sf(spacing);
96  s.CastFrom(sf);
97  this->SetOutputSpacing(s);
98 }
99 
100 template <class TVectorData, class TOutputImage>
101 void
103 ::SetOutputOrigin(const double origin[2])
104 {
105  OutputOriginType p(origin);
106  this->SetOutputOrigin(p);
107 }
108 
109 template <class TVectorData, class TOutputImage>
110 void
112 ::SetOutputOrigin(const float origin[2])
113 {
114  itk::Point<float, 2> of(origin);
116  p.CastFrom(of);
117  this->SetOutputOrigin(p);
118 }
119 
120 template <class TVectorData, class TOutputImage>
121 void
124 {
125  this->SetOutputOrigin ( src->GetOrigin() );
126  this->SetOutputSpacing ( src->GetSpacing() );
127  this->SetOutputSize ( src->GetLargestPossibleRegion().GetSize() );
129  this->SetOutputProjectionRef(imi->GetProjectionRef());
130 }
131 
132 template<class TVectorData, class TOutputImage>
133 void
136 {
137  // get pointer to the output
138  OutputImagePointer outputPtr = this->GetOutput();
139  if (!outputPtr)
140  {
141  return;
142  }
143 
144  // Set the size of the output region
145  typename TOutputImage::RegionType outputLargestPossibleRegion;
146  outputLargestPossibleRegion.SetSize(m_OutputSize);
147  //outputLargestPossibleRegion.SetIndex(m_OutputStartIndex);
148  outputPtr->SetLargestPossibleRegion(outputLargestPossibleRegion);
149 
150  // Set spacing and origin
151  outputPtr->SetSpacing(m_OutputSpacing);
152  outputPtr->SetOrigin(m_OutputOrigin);
153 
154  itk::MetaDataDictionary& dict = outputPtr->GetMetaDataDictionary();
155  itk::EncapsulateMetaData<std::string> (dict, MetaDataKey::ProjectionRefKey,
156  static_cast<std::string>(this->GetOutputProjectionRef()));
157 
158  // Generate the OGRLayers from the input VectorDatas
159  // iteration begin from 1 cause the 0th input is a image
160  for (unsigned int idx = 0; idx < this->GetNumberOfInputs(); ++idx)
161  {
162  const VectorDataType* vd = dynamic_cast< const VectorDataType*>(this->itk::ProcessObject::GetInput(idx));
163 
164  // Get the projection ref of the current VectorData
165  std::string projectionRefWkt = vd->GetProjectionRef();
166  bool projectionInformationAvailable = !projectionRefWkt.empty();
167  OGRSpatialReference * oSRS = NULL;
168 
169  if (projectionInformationAvailable)
170  {
171  oSRS = static_cast<OGRSpatialReference *>(OSRNewSpatialReference(projectionRefWkt.c_str()));
172  }
173  else
174  {
175  otbMsgDevMacro(<< "Projection information unavailable");
176  }
177 
178  // Retrieving root node
179  DataTreeConstPointerType tree = vd->GetDataTree();
180 
181  // Get the input tree root
182  InternalTreeNodeType * inputRoot = const_cast<InternalTreeNodeType *>(tree->GetRoot());
183 
184  // Iterative method to build the layers from a VectorData
185  OGRLayer * ogrCurrentLayer = NULL;
186  std::vector<OGRLayer *> ogrLayerVector;
188 
189  // The method ConvertDataTreeNodeToOGRLayers create the
190  // OGRDataSource but don t release it. Destruction is done in the
191  // desctructor
192  m_OGRDataSourcePointer = NULL;
193  ogrLayerVector = IOConversion->ConvertDataTreeNodeToOGRLayers(inputRoot,
194  m_OGRDataSourcePointer,
195  ogrCurrentLayer,
196  oSRS);
197 
198  // From OGRLayer* to OGRGeometryH vector
199  for (unsigned int idx = 0; idx < ogrLayerVector.size(); ++idx)
200  {
201  // test if the layers contain a field m_BurnField;
202  int burnField = -1;
203 
204  if( !m_BurnAttribute.empty() )
205  {
206  burnField = OGR_FD_GetFieldIndex( OGR_L_GetLayerDefn( (OGRLayerH)(ogrLayerVector[idx]) ),
207  m_BurnAttribute.c_str() );
208 
209  // Get the geometries of the layer
210  OGRFeatureH hFeat;
211  OGR_L_ResetReading( (OGRLayerH)(ogrLayerVector[idx]) );
212  while( ( hFeat = OGR_L_GetNextFeature( (OGRLayerH)(ogrLayerVector[idx]) )) != NULL )
213  {
214  OGRGeometryH hGeom;
215  if( OGR_F_GetGeometryRef( hFeat ) == NULL )
216  {
217  OGR_F_Destroy( hFeat );
218  continue;
219  }
220 
221  hGeom = OGR_G_Clone( OGR_F_GetGeometryRef( hFeat ) );
222  m_SrcDataSetGeometries.push_back( hGeom );
223 
224  if (burnField == -1 )
225  {
226  // TODO : if no burnAttribute available, warning or raise an exception??
227  m_FullBurnValues.push_back(m_DefaultBurnValue++);
228  itkWarningMacro(<<"Failed to find attribute "<<m_BurnAttribute << " in layer "
229  << OGR_FD_GetName( OGR_L_GetLayerDefn( (OGRLayerH)(ogrLayerVector[idx]) ))
230  <<" .Setting burn value to default = "
231  << m_DefaultBurnValue);
232  }
233  else
234  {
235  m_FullBurnValues.push_back( OGR_F_GetFieldAsDouble( hFeat, burnField ) );
236  }
237 
238  OGR_F_Destroy( hFeat );
239  }
240  }
241 
242  // Destroy the oSRS
243  if (oSRS != NULL)
244  {
245  OSRRelease(oSRS);
246  }
247  }
248  }
249 }
250 
251 template<class TVectorData, class TOutputImage>
252 void
254 {
255  // Call Superclass GenerateData
256  this->AllocateOutputs();
257 
258  // Get the buffered region
259  OutputImageRegionType bufferedRegion = this->GetOutput()->GetBufferedRegion();
260 
261  // nb bands
262  unsigned int nbBands = this->GetOutput()->GetNumberOfComponentsPerPixel();
263 
264  // register drivers
265  GDALAllRegister();
266 
267  std::ostringstream stream;
268  stream << "MEM:::"
269  << "DATAPOINTER=" << (unsigned long)(this->GetOutput()->GetBufferPointer()) << ","
270  << "PIXELS=" << bufferedRegion.GetSize()[0] << ","
271  << "LINES=" << bufferedRegion.GetSize()[1]<< ","
272  << "BANDS=" << nbBands << ","
273  << "DATATYPE=" << GDALGetDataTypeName(GdalDataTypeBridge::GetGDALDataType<OutputImageInternalPixelType>()) << ","
274  << "PIXELOFFSET=" << sizeof(OutputImageInternalPixelType) * nbBands << ","
275  << "LINEOFFSET=" << sizeof(OutputImageInternalPixelType)*nbBands*bufferedRegion.GetSize()[0] << ","
276  << "BANDOFFSET=" << sizeof(OutputImageInternalPixelType);
277 
278  GDALDatasetH dataset = GDALOpen(stream.str().c_str(), GA_Update);
279 
280  // Add the projection ref to the dataset
281  GDALSetProjection (dataset, this->GetOutput()->GetProjectionRef().c_str());
282 
283  // add the geoTransform to the dataset
284  itk::VariableLengthVector<double> geoTransform(6);
285 
286  // Reporting origin and spacing of the buffered region
287  // the spacing is unchanged, the origin is relative to the buffered region
288  OutputIndexType bufferIndexOrigin = bufferedRegion.GetIndex();
289  OutputOriginType bufferOrigin;
290  this->GetOutput()->TransformIndexToPhysicalPoint(bufferIndexOrigin, bufferOrigin);
291  geoTransform[0] = bufferOrigin[0];
292  geoTransform[3] = bufferOrigin[1];
293  geoTransform[1] = this->GetOutput()->GetSpacing()[0];
294  geoTransform[5] = this->GetOutput()->GetSpacing()[1];
295 
296  // FIXME: Here component 1 and 4 should be replaced by the orientation parameters
297  geoTransform[2] = 0.;
298  geoTransform[4] = 0.;
299  GDALSetGeoTransform(dataset,const_cast<double*>(geoTransform.GetDataPointer()));
300 
301  // Burn the geometries into the dataset
302  if (dataset != NULL)
303  {
304  GDALRasterizeGeometries( dataset, m_BandsToBurn.size(),
305  &(m_BandsToBurn[0]),
306  m_SrcDataSetGeometries.size(),
307  &(m_SrcDataSetGeometries[0]),
308  NULL, NULL, &(m_FullBurnValues[0]),
309  NULL,
310  GDALDummyProgress, NULL );
311 
312  // release the dataset
313  GDALClose( dataset );
314  }
315 }
316 
317 template<class TVectorData, class TOutputImage>
318 void
320 ::PrintSelf(std::ostream& os, itk::Indent indent) const
321 {
322  Superclass::PrintSelf(os, indent);
323 }
324 
325 } // end namespace otb
326 
327 #endif
328 

Generated at Sat Mar 8 2014 16:24:39 for Orfeo Toolbox with doxygen 1.8.3.1