diff --git a/src/asp/Sessions/RPC/RPCModel.cc b/src/asp/Sessions/RPC/RPCModel.cc index 0a0797b57..223d7600b 100644 --- a/src/asp/Sessions/RPC/RPCModel.cc +++ b/src/asp/Sessions/RPC/RPCModel.cc @@ -14,6 +14,12 @@ using namespace asp; using namespace vw; void asp::RPCModel::initialize( DiskImageResourceGDAL* resource ) { + // Extract Datum (by means of GeoReference) + GeoReference georef; + read_georeference( georef, *resource ); + m_datum = georef.datum(); + + // Extract RPC Info boost::shared_ptr dataset = resource->get_dataset_ptr(); if ( !dataset ) vw_throw( LogicErr() << "RPCModel: Could not read data. No file has been opened." ); @@ -24,19 +30,19 @@ void asp::RPCModel::initialize( DiskImageResourceGDAL* resource ) { vw_throw( LogicErr() << "RPCModel: GDAL resource appears not to have RPC metadata." ); // Copy information over to our data structures. - lonlatheight_offset = Vector3(gdal_rpc.dfLONG_OFF, - gdal_rpc.dfLAT_OFF, - gdal_rpc.dfHEIGHT_OFF); - lonlatheight_scale = Vector3(gdal_rpc.dfLONG_SCALE, - gdal_rpc.dfLAT_SCALE, - gdal_rpc.dfHEIGHT_SCALE); - line_os = Vector2(gdal_rpc.dfLINE_OFF, gdal_rpc.dfLINE_SCALE); - sample_os = Vector2(gdal_rpc.dfSAMP_OFF, gdal_rpc.dfSAMP_SCALE); - - m_line_num_coeff = Vector(gdal_rpc.adfLINE_NUM_COEFF); - m_line_den_coeff = Vector(gdal_rpc.adfLINE_DEN_COEFF); - m_sample_num_coeff = Vector(gdal_rpc.adfSAMP_NUM_COEFF); - m_sample_den_coeff = Vector(gdal_rpc.adfSAMP_DEN_COEFF); + m_lonlatheight_offset = Vector3(gdal_rpc.dfLONG_OFF, + gdal_rpc.dfLAT_OFF, + gdal_rpc.dfHEIGHT_OFF); + m_lonlatheight_scale = Vector3(gdal_rpc.dfLONG_SCALE, + gdal_rpc.dfLAT_SCALE, + gdal_rpc.dfHEIGHT_SCALE); + m_xy_offset = Vector2(gdal_rpc.dfSAMP_OFF, gdal_rpc.dfLINE_OFF); + m_xy_scale = Vector2(gdal_rpc.dfSAMP_SCALE, gdal_rpc.dfLINE_SCALE); + + m_line_num_coeff = CoeffVec(gdal_rpc.adfLINE_NUM_COEFF); + m_line_den_coeff = CoeffVec(gdal_rpc.adfLINE_DEN_COEFF); + m_sample_num_coeff = CoeffVec(gdal_rpc.adfSAMP_NUM_COEFF); + m_sample_den_coeff = CoeffVec(gdal_rpc.adfSAMP_DEN_COEFF); m_lonlat_bbox = BBox2(Vector2(gdal_rpc.dfMIN_LONG, gdal_rpc.dfMIN_LAT), Vector2(gdal_rpc.dfMAX_LONG, gdal_rpc.dfMAX_LAT)); } @@ -58,7 +64,7 @@ Vector2 asp::RPCModel::point_to_pixel( Vector3 const& point ) const { Vector3 lonlatrad = cartography::xyz_to_lon_lat_radius( point ); - Vector term = + CoeffVec term = calculate_terms( elem_quot(lonlatrad - lonlatheight_offset, lonlatheight_scale) ); diff --git a/src/asp/Sessions/RPC/RPCModel.h b/src/asp/Sessions/RPC/RPCModel.h index df7e0fdb4..7bf3f9060 100644 --- a/src/asp/Sessions/RPC/RPCModel.h +++ b/src/asp/Sessions/RPC/RPCModel.h @@ -10,6 +10,12 @@ #include #include +namespace vw { + namespace cartography { + class Datum; + } +} + namespace asp { // TODO: Need to work out a different way to triangulate. Our @@ -17,6 +23,15 @@ namespace asp { // * Pull Datum from file. LLH measurement are geodetic. class RPCModel : public vw::camera::CameraModel { + Datum m_datum; + + // Scaling parameters + vw::Vector m_line_num_coeff, m_line_den_coeff, + m_sample_num_coeff, m_sample_den_coeff; + vw::Vector2 m_xy_offset; + vw::Vector2 m_xy_scale; + vw::Vector3 m_lonlatheight_offset; + vw::Vector3 m_lonlatheight_scale; vw::BBox2 m_lonlat_bbox; void initialize( vw::DiskImageResourceGDAL* resource ); @@ -37,9 +52,21 @@ namespace asp { void inverse_transform( vw::Vector2 const& pix, vw::Vector3& point, vw::Vector3& direction ) const; - static vw::Vector - calculate_terms( vw::Vector3 const& v ) { - vw::Vector result; + // Access to constants + typedef vw::Vector CoeffVec; + Datum const& datum() const { return m_datum; } + CoeffVec const& line_num_coeff() const { return m_line_num_coeff; } + CoeffVec const& line_den_coeff() const { return m_line_den_coeff; } + CoeffVec const& sample_num_coeff() const { return m_sample_num_coeff; } + CoeffVec const& sample_den_coeff() const { return m_sample_den_coeff; } + vw::Vector2 const& xy_offset() const { return m_xy_offset; } + vw::Vector2 const& xy_scale() const { return m_xy_scale; } + vw::Vector3 const& lonlatheight_offset() const { return m_lonlatheight_offset; } + vw::Vector3 const& lonlatheight_scale() const { return m_lonlatheight_scale; } + vw::BBox2 const& lonlat_bbox() const { return m_lonlat_bbox; } + + static CoeffVec calculate_terms( vw::Vector3 const& v ) { + CoeffVec result; result[0] = 1.0; result[1] = v.x(); result[2] = v.y(); @@ -64,13 +91,6 @@ namespace asp { return result; } - // Scaling parameters - vw::Vector m_line_num_coeff, m_line_den_coeff, - m_sample_num_coeff, m_sample_den_coeff; - vw::Vector2 line_os; // [OFFSET, SCALE] - vw::Vector2 sample_os; - vw::Vector3 lonlatheight_offset; - vw::Vector3 lonlatheight_scale; }; } diff --git a/src/asp/Sessions/RPC/StereoSessionRPC.h b/src/asp/Sessions/RPC/StereoSessionRPC.h index 45d076ae0..4b2b5df04 100644 --- a/src/asp/Sessions/RPC/StereoSessionRPC.h +++ b/src/asp/Sessions/RPC/StereoSessionRPC.h @@ -10,11 +10,11 @@ #ifndef __STEREO_SESSION_RPC_H__ #define __STEREO_SESSION_RPC_H__ -#include +#include namespace asp { - class StereoSessionRPC : public StereoSession { + class StereoSessionRPC : public StereoSessionPinhole { public: virtual ~StereoSessionRPC() {} diff --git a/src/asp/Tools/stereo.cc b/src/asp/Tools/stereo.cc index 62b2bf316..bfcc40976 100644 --- a/src/asp/Tools/stereo.cc +++ b/src/asp/Tools/stereo.cc @@ -97,11 +97,18 @@ void vw::handle_arguments( int argc, char *argv[], Options& opt ) { // Some specialization here so that the user doesn't need to list // camera models on the command line for certain stereo session // types. (e.g. isis). + // + // TODO: This modification of arguments should probably happen in + // initialization and not be dependent on Stereo knowing what + // session it is in. bool check_for_camera_models = true; - if ( opt.stereo_session_string == "isis" ) { + if ( opt.stereo_session_string == "isis" || + opt.stereo_session_string == "rpc" ) { // Fix the ordering of the arguments if the user only supplies 3 - if (opt.out_prefix.empty()) + if (opt.out_prefix.empty()) { opt.out_prefix = opt.cam_file1; + opt.cam_file1.clear(); + } check_for_camera_models = false; } diff --git a/src/asp/Tools/stereo.in b/src/asp/Tools/stereo.in index 10037e97e..346ddd6fb 100644 --- a/src/asp/Tools/stereo.in +++ b/src/asp/Tools/stereo.in @@ -96,7 +96,10 @@ if __name__ == '__main__': if ( opt.entry_point <= 3 ): run('stereo_fltr', args, msg='3: Filtering') if ( opt.entry_point <= 4 ): - run('stereo_tri', args, msg='4: Triangulation') + tri = 'stereo_tri' + if ( opt.session.lowercase == "rpc" ): + tri += '_rpc' + run(tri, args, msg='4: Triangulation') except Exception, e: if not opt.debug: die(e) diff --git a/src/asp/Tools/stereo_tri_rpc.cc b/src/asp/Tools/stereo_tri_rpc.cc index 1203de0dc..0afd96e0b 100644 --- a/src/asp/Tools/stereo_tri_rpc.cc +++ b/src/asp/Tools/stereo_tri_rpc.cc @@ -45,6 +45,7 @@ int main( int argc, char* argv[] ) { // Triangulation with RPC is a little tricky. Here's the plan. // 1.) Solve for the mean XYZ for both cameras by looking at the // LLH boxes. + // 2.) Solve for one ray's intersection at the mean height. This // will seed the rigoruos algorithm. I believe this is equation 24 // in Grodecki 2004. We only need the first iteration's result.