otbSARConcatenateBursts.cxx 7.96 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64
/*
 * 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.
 */

#include "otbWrapperApplication.h"
#include "otbWrapperApplicationFactory.h"

#include "otbMultiChannelExtractROI.h"
#include "otbImageList.h"

#include "otbSarConcatenateBurstsImageFilter.h"

namespace otb
{
namespace Wrapper
{
class SARConcatenateBursts : public Application
{
public:
  /** Standard class typedefs. */
  typedef SARConcatenateBursts     Self;
  typedef Application                   Superclass;
  typedef itk::SmartPointer<Self>       Pointer;
  typedef itk::SmartPointer<const Self> ConstPointer;

  /** Standard macro */
  itkNewMacro(Self);

  itkTypeMacro(SARConcatenateBursts, otb::Application);

  typedef otb::ImageList<FloatVectorImageType>  ImageListType;
  typedef MultiChannelExtractROI<FloatVectorImageType::InternalPixelType,
				 FloatVectorImageType::InternalPixelType>       ExtractROIFilterType;
  typedef ObjectList<ExtractROIFilterType>                                      ExtractROIFilterListType;
 
  typedef otb::SarConcatenateBurstsImageFilter<FloatVectorImageType> BurstFusionFilterType;

private:
  void DoInit() override
  {
    SetName("SARConcatenateBursts");
    SetDescription("Concatenate a list of Bursts to provide a whole SAR Deburst Image.");

    // Documentation
    SetDocName("SAR Concatenate Bursts");
    SetDocLongDescription("This application performs a burst concatenation and provides a SAR Deburst Image. "
			  "It reads the input image list (single bursts) "
			  "and generates a whole SAR image with deburst operations.");
65 66 67
    SetDocLimitations("Only Sentinel1 IW SLC products are supported for now. In order to concatenate several" 
		      " bursts, all valid lines of each burst are required as inputs." 
		      "ie : Careful with ROI extraction inside a Burst.");
68 69 70
    SetDocAuthors("OTB-Team");

    AddDocTag(Tags::SAR);
71 72
    AddDocTag(Tags::Calibration);
    
73 74 75 76 77 78 79 80 81
    AddParameter(ParameterType_InputImageList,  "il",   "Input bursts list");
    SetParameterDescription("il", "The list of bursts to concatenate.");

    AddParameter(ParameterType_InputImage,  "insar",   "Input Sentinel1 IW SLC Image (only metadata used)");
    SetParameterDescription("insar", "Raw Sentinel1 IW SLC image, or any extract of such made by OTB (geom file needed).");

    AddParameter(ParameterType_OutputImage, "out",  "Output Image");
    SetParameterDescription("out", "The concatenated and debursted output image.");

82
    AddParameter(ParameterType_Int,  "burstindex", "Index of the first Burst");
Gaëlle USSEGLIO's avatar
Gaëlle USSEGLIO committed
83
    SetParameterDescription("burstindex", "Index for the first required Burst");
84 85 86
    MandatoryOff("burstindex");
    SetDefaultParameterInt("burstindex", 0);

87 88 89
    AddRAMParameter();

    // Doc example parameter settings
90
    SetDocExampleParameterValue("insar","s1_iw_slc.tif");
91 92 93 94 95 96 97 98 99 100 101
    SetDocExampleParameterValue("il", "Burst0.png Burst1.png");
    SetDocExampleParameterValue("out", "otbConcatenateBursts.tif");

    SetOfficialDocLink();
  }

  void DoUpdateParameters() override
  {}

  void DoExecute() override
  {
102 103 104
    // Get the first Burst index 
    int burst_index = GetParameterInt("burstindex");

105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127
    // Instanciate filters
    ExtractROIFilterListType::Pointer m_ExtractorList = ExtractROIFilterListType::New();
    ImageListType::Pointer m_ImageList = ImageListType::New();
    BurstFusionFilterType::Pointer fusionFilter = BurstFusionFilterType::New();
    
    // Get the input complex image list
    FloatVectorImageListType*  inList = GetParameterImageList("il");
    // Get the input complex image
    FloatVectorImageType*  in = GetParameterImage("insar");
    in->UpdateOutputInformation();

    std::vector<std::pair<unsigned long,unsigned long> > lines;
    std::vector<std::pair<unsigned long,unsigned long> > samples;

    // Check number of inputs with number of bursts into SLC image (must be the same)
    unsigned int nbBursts = 1;
    try
      {
	nbBursts = std::stoi(in->GetImageKeywordlist().GetMetadataByKey("support_data.geom.bursts.number"));
      }
    catch( ... )
      {
	// Throw an execption
Gaëlle USSEGLIO's avatar
Gaëlle USSEGLIO committed
128
	otbAppLogFATAL(<< "Failed to retrieve bursts.number value from .geom file.");
129 130
      }

131 132 133
    nbBursts = inList->Size();
    
    /*if (inList->Size() != nbBursts)
134 135
      {
	throw std::runtime_error("Failed to concatenate bursts. Some bursts are missing.");
136
	}*/ 
137 138 139

    // Coniguration for fusion filter
    fusionFilter->SetSLCImageKeyWorList(in->GetImageKeywordlist());
140 141 142
    // Get Invalid pixel Key (from Image 0) 
    FloatVectorImageType::Pointer Im0 = inList->GetNthElement(0);
    Im0->UpdateOutputInformation();
Gaëlle USSEGLIO's avatar
Gaëlle USSEGLIO committed
143 144 145 146 147 148 149
    
    auto const& kwl = Im0->GetImageKeywordlist();

    const bool inputWithInvalidPixels = kwl.HasKey("support_data.invalid_pixels")
      && kwl.GetMetadataByKey("support_data.invalid_pixels") == "yes";
 

150
    fusionFilter->getDeburstLinesAndSamples(lines, samples, burst_index, inputWithInvalidPixels);
151 152 153 154 155 156 157
       
    // Split each input burst to keep only interested region
    for( unsigned int i=0; i<inList->Size(); i++ )
      {
	FloatVectorImageType::Pointer vectIm = inList->GetNthElement(i);
	vectIm->UpdateOutputInformation();

158
	// Check invalid Pixel Key
Gaëlle USSEGLIO's avatar
Gaëlle USSEGLIO committed
159 160 161
	const bool inputWithInvalidPixels_loop = vectIm->GetImageKeywordlist().HasKey("support_data.invalid_pixels")
	  && vectIm->GetImageKeywordlist().GetMetadataByKey("support_data.invalid_pixels") == "yes";

162 163 164
	if (inputWithInvalidPixels_loop != inputWithInvalidPixels)
	  {
	    // Throw an execption
Gaëlle USSEGLIO's avatar
Gaëlle USSEGLIO committed
165
	    otbAppLogFATAL(<< "Incoherency between input images (for support_data.invalid_pixels key).");
166 167
	  }

168 169 170
	unsigned long originOffset_samples = static_cast<long>(vectIm->GetOrigin()[0]-0.5);
	unsigned long originOffset_lines = static_cast<long>(vectIm->GetOrigin()[1]-0.5);

171
	// Retrieve start and size for each burst
Gaëlle USSEGLIO's avatar
Gaëlle USSEGLIO committed
172 173
	std::pair<unsigned long,unsigned long> line = lines[i];
	std::pair<unsigned long,unsigned long> sample = samples[i];
174

Gaëlle USSEGLIO's avatar
Gaëlle USSEGLIO committed
175 176
	unsigned long minSamples = std::min(sample.second, static_cast<unsigned long>(vectIm->GetLargestPossibleRegion().GetSize()[0]-1));
	unsigned long minLines = std::min(line.second, static_cast<unsigned long>(vectIm->GetLargestPossibleRegion().GetSize()[1]-1));
177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192

	unsigned long startL = line.first - originOffset_lines;
	unsigned long sizeL = minLines - line.first + 1;
	unsigned long startS = sample.first - originOffset_samples;
	unsigned long sizeS = minSamples - sample.first + 1;
	// Readjust if origin is superior to the first selected line/sample for the current burst
	if (line.first < originOffset_lines) 
	  {
	    startL = 0;
	    sizeL =  minLines - line.first - originOffset_lines + 1;
	  }
	if (sample.first < originOffset_samples) 
	  {
	    startS = 0;
	    sizeS =  minSamples - sample.first - originOffset_samples + 1;
	  }
193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227

	ExtractROIFilterType::Pointer extractor = ExtractROIFilterType::New();
	extractor->SetInput(vectIm);
	extractor->SetStartX(startS);
	extractor->SetStartY(startL);
	extractor->SetSizeX(sizeS);
	extractor->SetSizeY(sizeL);
	extractor->UpdateOutputInformation();

	m_ExtractorList->PushBack( extractor );
	m_ImageList->PushBack( extractor->GetOutput() );
	
      }

    BurstFusionFilterType::SizeType layout;
    layout[0] = 1;
    layout[1] = nbBursts;
    fusionFilter->SetLayout(layout);

    for (unsigned int i=0; i<(layout[0]*layout[1]); i++)
      {
	fusionFilter->SetInput(i,m_ImageList->GetNthElement(i));
      }

    SetParameterOutputImage("out", fusionFilter->GetOutput());

    RegisterPipeline();
  }


};
}
}

OTB_APPLICATION_EXPORT(otb::Wrapper::SARConcatenateBursts)