diff --git a/app/otbSARAltAmbig.cxx b/app/otbSARAltAmbig.cxx
index 49e3a90a1da2aacaa6fa578643fc4fb5c6dc005e..9a01b22701840ad06d7932f43b01d26b91b960df 100644
--- a/app/otbSARAltAmbig.cxx
+++ b/app/otbSARAltAmbig.cxx
@@ -28,6 +28,7 @@
#include <iomanip>
#include <string>
#include <fstream>
+#include <vector>
namespace otb
{
@@ -48,12 +49,12 @@ private:
void DoInit() override
{
SetName("SARAltAmbig");
- SetDescription("Evaluation of altitudes of ambiguity of orbits");
+ SetDescription("Evaluation of ambiguity height");
SetDocLongDescription("Calculate the angle alpha between "
"ground-master and ground-slave vectors. "
- "Calculates the local incidence.\nThe ambiguity "
- "altitude is given by: "
+ "Calculates the local incidence.\n The ambiguity "
+ "heigth is given by: "
"lambda * sin (incidence) / 2 * tg (alpha)");
//Optional descriptors
@@ -63,7 +64,7 @@ private:
AddDocTag("DiapOTB");
//Parameter declarations
- AddParameter(ParameterType_InputImage, "inmaster", "Input Master geosar image");
+ AddParameter(ParameterType_InputImage, "inmaster", "Input Master image");
SetParameterDescription("inmaster", "SAR Image master");
AddParameter(ParameterType_Float, "lat", "Tartget lat");
@@ -81,7 +82,7 @@ private:
SetMinimumParameterFloatValue("height", -100.);
SetMaximumParameterFloatValue("height", 9000.);
- AddParameter(ParameterType_InputImage, "inslave", "Input Slave geosar image");
+ AddParameter(ParameterType_InputImage, "inslave", "Input Slave image");
SetParameterDescription("inslave", "SAR image slave");
AddParameter(ParameterType_OutputFilename, "outfile", "Output file to store the results");
@@ -115,27 +116,24 @@ private:
// Start the target pipeline (Read SAR image metadata)
ComplexFloatImageType::Pointer SARPtr_slave = GetParameterComplexFloatImage("inslave");
- // Open the output file
- std::ofstream output(GetParameterString("outfile"), std::ios::out | std::ios::trunc);
-
+
// Get lat, lon and height
float lat = GetParameterFloat("lat");
otbAppLogINFO(<<"Lat : "<<lat);
- output << "Latitude of the target = " << lat << std::endl;
-
+
float lon = GetParameterFloat("lon");
otbAppLogINFO(<<"Lon : "<<lon);
- output << "Longitude of the target = " << lon << std::endl;
+
if (lon < 0)
lon = 360.0 + lon;
float height = GetParameterFloat("height");
otbAppLogINFO(<<"Height : "<<height);
- output << "Altitude of the target = " << height << std::endl;
-
+
bool bistatic = IsParameterEnabled("bistatic");
otbAppLogINFO(<<"Bistatic : " << bistatic);
+ // Get image keywordlists
otb::ImageKeywordlist sarKWL = SARPtr->GetImageKeywordlist();
otb::ImageKeywordlist sarKWL_slave = SARPtr_slave->GetImageKeywordlist();
@@ -147,56 +145,94 @@ private:
SarSensorModelAdapter::Point3DType xyzCart;
SarSensorModelAdapter::Point3DType satellitePosition[2];
SarSensorModelAdapter::Point3DType satelliteVelocity[2];
- double incidence[2];
demGeoPoint[0] = lon;
demGeoPoint[1] = lat;
demGeoPoint[2] = height;
+ // Results allocations
+ double * resultsLocalVector = new double[3];
+ double * resultsAltAmbig = new double[4];
+
+
+ // Cartesian conversion
SarSensorModelAdapter::WorldToCartesian(demGeoPoint, xyzCart);
- otbAppLogINFO(<<"Cartesian coords : " << xyzCart[0] << " " << xyzCart[1] << " " << xyzCart[2]);
- output << "GROUND COORDINATES = " << xyzCart[0] << " " << xyzCart[1] << " " << xyzCart[2] << std::endl;
+ otbAppLogDEBUG(<<"Cartesian coords : " << xyzCart[0] << " " << xyzCart[1] << " " << xyzCart[2]);
+
+ // Image names
+ std::string master = getFilename(GetParameterString("inmaster"));
+ std::string slave = getFilename(GetParameterString("inslave"));
+ std::vector<std::string> imgNames;
+ imgNames.push_back(master);
+ imgNames.push_back(slave);
+
+ // Satellite position and velocity
// Master
- otbAppLogINFO(<<"Inmaster Image : ");
+ otbAppLogDEBUG(<<"Inmaster Image : ");
m_SarSensorModelAdapter->WorldToSatPositionAndVelocity(demGeoPoint, satellitePosition[0], satelliteVelocity[0]);
- otbAppLogINFO(<<"Satellite Position : " << satellitePosition[0][0] << " " << satellitePosition[0][1] << " " << satellitePosition[0][2]);
- otbAppLogINFO(<<"Satellite Velocity : " << satelliteVelocity[0][0] << " " << satelliteVelocity[0][1] << " " << satelliteVelocity[0][2]);
+ otbAppLogDEBUG(<<"Satellite Position : " << satellitePosition[0][0] << " " << satellitePosition[0][1] << " " << satellitePosition[0][2]);
+ otbAppLogDEBUG(<<"Satellite Velocity : " << satelliteVelocity[0][0] << " " << satelliteVelocity[0][1] << " " << satelliteVelocity[0][2]);
+ // Slave
+ loadOk = m_SarSensorModelAdapter->LoadState(sarKWL_slave);
+ otbAppLogDEBUG(<<"Target Image : ");
+ m_SarSensorModelAdapter->WorldToSatPositionAndVelocity(demGeoPoint, satellitePosition[1], satelliteVelocity[1]);
+ otbAppLogDEBUG(<<"Satellite Position : " << satellitePosition[1][0] << " " << satellitePosition[1][1] << " " << satellitePosition[1][2]);
+ otbAppLogDEBUG(<<"Satellite Velocity : " << satelliteVelocity[1][0] << " " << satelliteVelocity[1][1] << " " << satelliteVelocity[1][2]);
+
// Calculate lambda
double radarFreq = atof(sarKWL.GetMetadataByKey("support_data.radar_frequency").c_str());
const double C = 299792458.;
-
float lambda = C/radarFreq;
- otbAppLogINFO(<<"Lambda : " << lambda);
+ otbAppLogDEBUG(<<"Lambda : " << lambda);
+
+ // Estimations : mean incidence, ambig_heigth, radial coord, lateral coord and local projections for master
+ otbAppLogDEBUG(<< "Compute Alt Ambig");
+ computeAltAmbig(xyzCart, satellitePosition, satelliteVelocity, lambda, imgNames, resultsLocalVector, resultsAltAmbig);
+
+
+ // Write into output file
+ // Open the output file
+ std::ofstream output(GetParameterString("outfile"), std::ios::out | std::ios::trunc);
+
+ output << "Latitude of the target = " << lat << std::endl;
+ output << "Longitude of the target = " << lon << std::endl;
+ output << "Altitude of the target = " << height << std::endl;
+
+ output << "GROUND COORDINATES = " << xyzCart[0] << " " << xyzCart[1] << " " << xyzCart[2] << std::endl;
- // Target
- loadOk = m_SarSensorModelAdapter->LoadState(sarKWL_slave);
- otbAppLogINFO(<<"Target Image : ");
- m_SarSensorModelAdapter->WorldToSatPositionAndVelocity(demGeoPoint, satellitePosition[1], satelliteVelocity[1]);
- otbAppLogINFO(<<"Satellite Position : " << satellitePosition[1][0] << " " << satellitePosition[1][1] << " " << satellitePosition[1][2]);
- otbAppLogINFO(<<"Satellite Velocity : " << satelliteVelocity[1][0] << " " << satelliteVelocity[1][1] << " " << satelliteVelocity[1][2]);
-
- otbAppLogINFO(<< "Compute Alt Ambig");
- otbAppLogINFO(<< std::left
- << std::setw(60) << "MASTER IMAGE"
- << std::setw(60) << "SLAVE IMAGE"
- << std::setw(12) << "Incidence"
- << std::setw(12) << "ALT AMBIG"
- << std::setw(12) << "Radial"
- << std::setw(12) << "Lateral");
output << std::endl << std::left
- << std::setw(60) << "MASTER IMAGE"
- << std::setw(60) << "SLAVE IMAGE"
+ << std::setw(70) << "MASTER IMAGE"
+ << std::setw(70) << "SLAVE IMAGE"
<< std::setw(12) << "Incidence"
<< std::setw(12) << "ALT AMBIG"
<< std::setw(12) << "Radial"
<< std::setw(12) << "Lateral" << std::endl;
- computeAltAmbig(xyzCart, satellitePosition, satelliteVelocity, lambda, output);
+
+
+ output << std::left
+ << std::setw(70) << master
+ << std::setw(70) << slave
+ << std::setw(12) << resultsAltAmbig[0]
+ << std::setw(12) << resultsAltAmbig[1]
+ << std::setw(12) <<resultsAltAmbig[2]
+ << std::setw(12) << resultsAltAmbig[3] << std::endl;
+
+ output << "Ground target --> Satellite Vector" << std::endl;
+ output << "North Projection : " << resultsLocalVector[0] << std::endl;
+ output << "East Projection : " << resultsLocalVector[1] << std::endl;
+ output << "Vertical Projection : " << resultsLocalVector[2] << std::endl;
// Close the output file
output.close();
+
+ // Free Memory
+ delete resultsLocalVector;
+ delete resultsAltAmbig;
+ resultsLocalVector = 0;
+ resultsAltAmbig = 0;
}
std::string getFilename(std::string filePath){
@@ -205,7 +241,13 @@ private:
return filePath.substr(sepPos + 1, dotPos - sepPos - 1);
}
- double computeIncidence(SarSensorModelAdapter::Point3DType satellitePosition, SarSensorModelAdapter::Point3DType xyzCart, SarSensorModelAdapter::Point3DType rm){
+
+ ///////// compute functions /////////
+ // Incidence
+ double computeIncidence(SarSensorModelAdapter::Point3DType satellitePosition,
+ SarSensorModelAdapter::Point3DType xyzCart,
+ SarSensorModelAdapter::Point3DType rm)
+ {
// Incidence
double normeS = sqrt(satellitePosition[0]*satellitePosition[0] +
satellitePosition[1]*satellitePosition[1] +
@@ -218,13 +260,20 @@ private:
return incidence;
}
- void searchPassageCloserToGround(SarSensorModelAdapter::Point3DType satellitePosition, SarSensorModelAdapter::Point3DType xyzCart, std::ofstream& output, SarSensorModelAdapter::Point3DType rm){
+ // Local Vector
+ void computeLocalVector(SarSensorModelAdapter::Point3DType satellitePosition,
+ SarSensorModelAdapter::Point3DType xyzCart,
+ SarSensorModelAdapter::Point3DType rm,
+ double * resultsLocalVector)
+ {
double normeRm = sqrt(rm[0] * rm[0] + rm[1] * rm[1] + rm[2] * rm[2]);
// Calculate the local vectors
SarSensorModelAdapter::Point3DType vecti;
SarSensorModelAdapter::Point3DType vectj;
SarSensorModelAdapter::Point3DType vectk;
+
+ // lon an lat into radians
double lon = GetParameterFloat("lon")*M_PI/180.;
double lat = GetParameterFloat("lat")*M_PI/180.;
@@ -241,19 +290,27 @@ private:
vectk[2] = sin(lat);
// Rm projection
- output << "Ground target --> Satellite Vector" << std::endl;
+ // Use -Rm because the given rm, here, is egal to : -(xyzCart - satellitePosition)
double northProjection = -(vecti[0]*(-rm[0])+vecti[1]*(-rm[1])+vecti[2]*(-rm[2]))/normeRm;
double eastProjection = -(vectj[0]*(-rm[0])+vectj[1]*(-rm[1])+vectj[2]*(-rm[2]))/normeRm;
double verticalProjection = -(vectk[0]*(-rm[0])+vectk[1]*(-rm[1])+vectk[2]*(-rm[2]))/normeRm;
- otbAppLogINFO(<< "North Projection : " << northProjection);
- otbAppLogINFO(<< "East Projection : " << eastProjection);
- otbAppLogINFO(<< "Vertical Projection : " << verticalProjection);
- output << "North Projection : " << northProjection << std::endl;
- output << "East Projection : " << eastProjection << std::endl;
- output << "Vertical Projection : " << verticalProjection << std::endl;
+ otbAppLogDEBUG(<< "North Projection : " << northProjection);
+ otbAppLogDEBUG(<< "East Projection : " << eastProjection);
+ otbAppLogDEBUG(<< "Vertical Projection : " << verticalProjection);
+
+ // Assign output
+ resultsLocalVector[0] = northProjection;
+ resultsLocalVector[1] = eastProjection;
+ resultsLocalVector[2] = verticalProjection;
}
- void computeAltAmbig(SarSensorModelAdapter::Point3DType xyzCart, SarSensorModelAdapter::Point3DType* satellitePosition, SarSensorModelAdapter::Point3DType* satelliteVelocity, double lambda, std::ofstream& output){
+ // Alt Ambig
+ void computeAltAmbig(SarSensorModelAdapter::Point3DType xyzCart,
+ SarSensorModelAdapter::Point3DType* satellitePosition,
+ SarSensorModelAdapter::Point3DType* satelliteVelocity,
+ double lambda, std::vector<std::string> imgNames,
+ double * resultsLocalVector, double * results)
+ {
// Adjust the formula in case of non-monostatic mode
double factorBperp = (!IsParameterEnabled("bistatic"))? 1.0:2.0;
double factorHamb = (!IsParameterEnabled("bistatic"))? 2.0:1.0;
@@ -313,24 +370,24 @@ private:
double inc_moy = ((incidenceRe + incidenceRm) / 2) * M_PI / 180.;
double alt_amb = lambda * sin(inc_moy) / (factorHamb * tan(alpha));
- std::string master = getFilename(GetParameterString("inmaster"));
- std::string slave = getFilename(GetParameterString("inslave"));
- otbAppLogINFO(<< std::left
- << std::setw(60) << master
- << std::setw(60) << slave
+ std::string master = imgNames[0];
+ std::string slave = imgNames[1];
+ otbAppLogDEBUG(<< std::left
+ << std::setw(70) << master
+ << std::setw(70) << slave
<< std::setw(12) << inc_moy * 180. / M_PI
<< std::setw(12) << alt_amb
<< std::setw(12) << xe
<< std::setw(12) << ye / factorBperp);
- output << std::left
- << std::setw(60) << master
- << std::setw(60) << slave
- << std::setw(12) << inc_moy * 180. / M_PI
- << std::setw(12) << alt_amb
- << std::setw(12) << xe
- << std::setw(12) << ye / factorBperp << std::endl;
-
- searchPassageCloserToGround(satellitePosition[0], xyzCart, output, rm);
+
+ // Assign output
+ results[0] = inc_moy * 180. / M_PI; // Mean incidence (in degrees)
+ results[1] = alt_amb; // Ambiguity height
+ results[2] = xe; // Radial coordinate
+ results[3] = ye / factorBperp; // Lateral coordinate
+
+ // Local vector estimation (for master only)
+ computeLocalVector(satellitePosition[0], xyzCart, rm, resultsLocalVector);
}
};