diff --git a/Code/Common/otbRemoteSensingRegion.h b/Code/Common/otbRemoteSensingRegion.h index c6612fdca4f0b852e2ccfd9eb0623c1bf60cd159..9ead9c344c3224fa56a035fbef679251f43d221c 100644 --- a/Code/Common/otbRemoteSensingRegion.h +++ b/Code/Common/otbRemoteSensingRegion.h @@ -381,6 +381,7 @@ TransformPhysicalRegionToIndexRegion(const RemoteSensingRegionType& region, cons point[0] = region.GetIndex()[0]; point[1] = region.GetIndex()[1]; image->TransformPhysicalPointToIndex(point, index); + std::cout <<" physicalPoint : "<< point << " --> index "<< index << std::endl; point[0] = region.GetIndex()[0] + region.GetSize()[0]; point[1] = region.GetIndex()[1] + region.GetSize()[1]; @@ -398,42 +399,6 @@ TransformPhysicalRegionToIndexRegion(const RemoteSensingRegionType& region, cons return outputRegion; } -template<class ImageType, class RemoteSensingRegionType> -typename ImageType::RegionType -TransformContinousRegionToIndexRegion(const RemoteSensingRegionType& region, const ImageType* image) -{ - std::cout << "TransformContinousRegionToIndexRegion" <<std::endl; - typename ImageType::RegionType outputRegion; - typename ImageType::RegionType::IndexType index; - typename ImageType::RegionType::IndexType index2; - - typedef typename ImageType::RegionType::IndexType ImageIndexType; - typedef typename ImageIndexType::IndexValueType ImageIndexValueType; - - typename ImageType::PointType point; - point[0] = region.GetIndex()[0]; - point[1] = region.GetIndex()[1]; - - index[0] = static_cast<ImageIndexValueType>(point[0]); - index[1] = static_cast<ImageIndexValueType>(point[1]); - - point[0] = region.GetIndex()[0] + region.GetSize()[0]; - point[1] = region.GetIndex()[1] + region.GetSize()[1]; - - index2[0] = static_cast<ImageIndexValueType>(point[0]); - index2[1] = static_cast<ImageIndexValueType>(point[1]); - - typename ImageType::RegionType::SizeType size; - size[0] = std::abs(index2[0] - index[0]); - size[1] = std::abs(index2[1] - index[1]); - - index[0] = std::min(index2[0], index[0]); - index[1] = std::min(index2[1], index[1]); - - outputRegion.SetIndex(index); - outputRegion.SetSize(size); - return outputRegion; -} } // end namespace otb diff --git a/Code/Learning/otbListSampleGenerator.h b/Code/Learning/otbListSampleGenerator.h index 444ccae5227ad4139b96531e3d39aa77b7276cda..19d2bf4d969c7bf342283a361fcefd4627d59ee2 100644 --- a/Code/Learning/otbListSampleGenerator.h +++ b/Code/Learning/otbListSampleGenerator.h @@ -24,10 +24,6 @@ #include "itkPreOrderTreeIterator.h" #include "itkMersenneTwisterRandomVariateGenerator.h" -#include "otbVectorDataProjectionFilter.h" -#include "otbVectorDataExtractROI.h" - - namespace otb { /** \class ListSampleGenerator @@ -84,12 +80,6 @@ public: typedef itk::Statistics::ListSample<LabelType> ListLabelType; typedef typename ListLabelType::Pointer ListLabelPointerType; - // Reprojection filter - typedef VectorDataProjectionFilter<VectorDataType,VectorDataType> - VectorDataProjectionFilterType; - typedef VectorDataExtractROI<VectorDataType> VectorDataExtractROIType; - typedef typename VectorDataExtractROIType::RegionType RemoteSensingRegionType; - /** Connects the input image for which the sample list is going to be extracted */ void SetInput(const ImageType *); const ImageType * GetInput() const; @@ -144,8 +134,6 @@ public: void GenerateClassStatistics(); - void ReprojectVectorData(const ImagePointerType image, VectorDataPointerType vectorDataInput); - protected: ListSampleGenerator(); virtual ~ListSampleGenerator() {} diff --git a/Code/Learning/otbListSampleGenerator.txx b/Code/Learning/otbListSampleGenerator.txx index e3f665490a7cec165768f74607b05e19befeac99..3dd93396159725e26694f67d798fac39dd998e5c 100644 --- a/Code/Learning/otbListSampleGenerator.txx +++ b/Code/Learning/otbListSampleGenerator.txx @@ -22,9 +22,40 @@ #include "itkImageRegionConstIteratorWithIndex.h" #include "otbRemoteSensingRegion.h" +#include "otbVectorDataProjectionFilter.h" namespace otb { + +template <class TVectorData> +void printVectorData(TVectorData * vectorData, string msg = "") +{ + typedef TVectorData VectorDataType; + typedef itk::PreOrderTreeIterator<typename VectorDataType::DataTreeType> TreeIteratorType; + + TreeIteratorType itVector(vectorData->GetDataTree()); + itVector.GoToBegin(); + + if (!msg.empty()) + { + std::cout<< msg << std::endl; + } + + while (!itVector.IsAtEnd()) + { + if (itVector.Get()->IsPolygonFeature()) + { + std::cout << itVector.Get()->GetNodeTypeAsString() << std::endl; + for (unsigned int itPoints = 0; itPoints < itVector.Get()->GetPolygonExteriorRing()->GetVertexList()->Size(); itPoints++) + { + std::cout << "vertex[" << itPoints << "]: " << itVector.Get()->GetPolygonExteriorRing()->GetVertexList()->GetElement(itPoints) <<std::endl; + } + std::cout << "Polygon bounding region:\n" << itVector.Get()->GetPolygonExteriorRing()->GetBoundingRegion() << std::endl; + } + ++itVector; + } +} + template<class TImage, class TVectorData> ListSampleGenerator<TImage, TVectorData> ::ListSampleGenerator() : @@ -77,25 +108,8 @@ ListSampleGenerator<TImage, TVectorData> this->ProcessObject::SetNthInput(1, const_cast<VectorDataType *>(vectorData)); - TreeIteratorType itVector(vectorData->GetDataTree()); - itVector.GoToBegin(); + printVectorData(vectorData); - int idPolygon=1; - while (!itVector.IsAtEnd()) - { - if (itVector.Get()->IsPolygonFeature()) - { - std::cout << "POLYGON " << idPolygon << ": " <<std::endl; - for (unsigned int itPoints = 0; itPoints < itVector.Get()->GetPolygonExteriorRing()->GetVertexList()->Size(); itPoints++) - { - std::cout << "vertex[" << itPoints << "]: " << itVector.Get()->GetPolygonExteriorRing()->GetVertexList()->GetElement(itPoints) <<std::endl; - } - std::cout << "Polygon region: " << itVector.Get()->GetPolygonExteriorRing()->GetBoundingRegion() << std::endl; - ++idPolygon; - } - ++itVector; - - } std::cout << "ListSampleGenerator::SetInputVectorData ... END" << std::endl; } @@ -145,9 +159,6 @@ ListSampleGenerator<TImage, TVectorData> ImagePointerType image = const_cast<ImageType*>(this->GetInput()); - // Reproject VectorData - this->ReprojectVectorData(image, const_cast<VectorDataType*>(this->GetInputVectorData())); - VectorDataPointerType vectorData = const_cast<VectorDataType*>(this->GetInputVectorData()); //Gather some information about the relative size of the classes @@ -173,21 +184,16 @@ ListSampleGenerator<TImage, TVectorData> { if (itVector.Get()->IsPolygonFeature()) { - /*typename ImageType::RegionType polygonRegion = - otb::TransformPhysicalRegionToIndexRegion(itVector.Get()->GetPolygonExteriorRing()->GetBoundingRegion(), - image.GetPointer());*/ - typename ImageType::RegionType polygonRegion = - otb::TransformContinousRegionToIndexRegion(itVector.Get()->GetPolygonExteriorRing()->GetBoundingRegion(), - image.GetPointer()); - + otb::TransformPhysicalRegionToIndexRegion(itVector.Get()->GetPolygonExteriorRing()->GetBoundingRegion(), + image.GetPointer()); - std::cout << "Image region from polygon: " << polygonRegion << std::endl; - std::cout << "Image region : " << image->GetLargestPossibleRegion() << std::endl; + std::cout << "Image region from polygon:\n" << polygonRegion << std::endl; + std::cout << "Image largest possible region:\n" << image->GetLargestPossibleRegion() << std::endl; image->SetRequestedRegion(polygonRegion); image->PropagateRequestedRegion(); image->UpdateOutputData(); - std::cout << "Image region requested: " << image->GetRequestedRegion() << std::endl; + std::cout << "Image region requested:\n" << image->GetRequestedRegion() << std::endl; typedef itk::ImageRegionConstIteratorWithIndex<ImageType> IteratorType; IteratorType it(image, polygonRegion); @@ -197,8 +203,8 @@ ListSampleGenerator<TImage, TVectorData> itk::ContinuousIndex<double, 2> point; image->TransformIndexToPhysicalPoint(it.GetIndex(), point); //std::cout << it.GetIndex() << " -> " << point << std::endl; - if (itVector.Get()->GetPolygonExteriorRing()->IsInside(it.GetIndex())) - //if (itVector.Get()->GetPolygonExteriorRing()->IsInside(point)) + //if (itVector.Get()->GetPolygonExteriorRing()->IsInside(it.GetIndex())) + if (itVector.Get()->GetPolygonExteriorRing()->IsInside(point)) { double randomValue = m_RandomGenerator->GetUniformVariate(0.0, 1.0); if (randomValue < m_ClassesProbTraining[itVector.Get()->GetFieldAsInt(m_ClassKey)]) @@ -230,81 +236,6 @@ ListSampleGenerator<TImage, TVectorData> std::cout << "ListSampleGenerator::GenerateData() : ... END" <<std::endl; } -template <class TImage, class TVectorData> -void -ListSampleGenerator<TImage, TVectorData> -::ReprojectVectorData(const ImagePointerType image, VectorDataPointerType vectorROIs) -{ - typedef typename ImageType::IndexType IndexType; - typedef typename ImageType::PointType PointType; - - if(image.IsNull()) - { - itkExceptionMacro("Invalid input image."); - } - - // Vector data reprojection - typename VectorDataProjectionFilterType::Pointer vproj; - typename VectorDataExtractROIType::Pointer vdextract; - - // Extract The part of the VectorData that actually overlaps with - // the image extent - vdextract = VectorDataExtractROIType::New(); - vdextract->SetInput(vectorROIs); - - // Find the geographic region of interest - - // Ge the index of the corner of the image - IndexType ul, ur, ll, lr; - PointType pul, pur, pll, plr; - ul = image->GetLargestPossibleRegion().GetIndex(); - ur = ul; - ll = ul; - lr = ul; - ur[0] += image->GetLargestPossibleRegion().GetSize()[0]; - lr[0] += image->GetLargestPossibleRegion().GetSize()[0]; - lr[1] += image->GetLargestPossibleRegion().GetSize()[1]; - ll[1] += image->GetLargestPossibleRegion().GetSize()[1]; - - // Transform to physical point - image->TransformIndexToPhysicalPoint(ul, pul); - image->TransformIndexToPhysicalPoint(ur, pur); - image->TransformIndexToPhysicalPoint(ll, pll); - image->TransformIndexToPhysicalPoint(lr, plr); - - // Build the cartographic region - RemoteSensingRegionType rsRegion; - typename RemoteSensingRegionType::IndexType rsOrigin; - typename RemoteSensingRegionType::SizeType rsSize; - rsOrigin[0] = min(pul[0], plr[0]); - rsOrigin[1] = min(pul[1], plr[1]); - rsSize[0] = vcl_abs(pul[0] - plr[0]); - rsSize[1] = vcl_abs(pul[1] - plr[1]); - - rsRegion.SetOrigin(rsOrigin); - rsRegion.SetSize(rsSize); - rsRegion.SetRegionProjection(image->GetProjectionRef()); - rsRegion.SetKeywordList(image->GetImageKeywordlist()); - - // Set the cartographic region to the extract roi filter - vdextract->SetRegion(rsRegion); - - // Reproject VectorData in image projection - vproj = VectorDataProjectionFilterType::New(); - vproj->SetInput(vdextract->GetOutput()); - vproj->SetInputProjectionRef(vectorROIs->GetProjectionRef()); - vproj->SetOutputKeywordList(image->GetImageKeywordlist()); - vproj->SetOutputProjectionRef(image->GetProjectionRef()); - vproj->SetOutputOrigin(image->GetOrigin()); - vproj->SetOutputSpacing(image->GetSpacing()); - - vproj->Update(); - - // Update InputData - std::cout<< "-----REPROJECTED VECTOR DATA-----" << std::endl; - this->SetInputVectorData(vproj->GetOutput()); -} - template <class TImage, class TVectorData> void ListSampleGenerator<TImage, TVectorData> @@ -325,7 +256,6 @@ ListSampleGenerator<TImage, TVectorData> { if (itVector.Get()->IsPolygonFeature()) { - std::cout << itVector.Get()->GetNodeTypeAsString() << std::endl; m_ClassesSize[itVector.Get()->GetFieldAsInt(m_ClassKey)] += itVector.Get()->GetPolygonExteriorRing()->GetArea() / pixelArea; // in pixel std::cout << "Area = "<<itVector.Get()->GetPolygonExteriorRing()->GetArea() << std::endl;