diff --git a/src/vw/Cartography/GeoTransform.cc b/src/vw/Cartography/GeoTransform.cc index dc1f2518d..d1e8cf4c9 100644 --- a/src/vw/Cartography/GeoTransform.cc +++ b/src/vw/Cartography/GeoTransform.cc @@ -36,43 +36,7 @@ namespace cartography { // Constructor GeoTransform::GeoTransform(GeoReference const& src_georef, GeoReference const& dst_georef, BBox2i const& src_bbox, BBox2i const& dst_bbox) : - m_src_georef(src_georef), m_dst_georef(dst_georef), m_offset(Vector2(0, 0)) { - - // Deal with the fact that longitudes could differ by 360 degrees - - // between src and dst. - if (src_bbox.empty() || dst_bbox.empty()) { - // If we don't know the image areas, simply use the 0, 0 corner - try { - Vector2 src_origin = src_georef.pixel_to_lonlat(Vector2(0, 0)); - Vector2 dst_origin = dst_georef.pixel_to_lonlat(Vector2(0, 0)); - m_offset = Vector2( 360.0*round( (dst_origin[0] - src_origin[0])/360.0 ), 0.0 ); - }catch ( const std::exception & e ) { - m_offset = Vector2(0, 0); - } - } else{ - - // Try to offset by 360 degrees until the lon-lat boxes are most compatible, - // measured by how close the centers of the boxes are. - - BBox2 src_lonlat_box = src_georef.pixel_to_lonlat_bbox(src_bbox); - BBox2 dst_lonlat_box = dst_georef.pixel_to_lonlat_bbox(dst_bbox); - std::vector shift, dist_vec; - for (int val = -360; val <= 360; val+= 360) { - shift.push_back(val); - BBox2 shifted_src = src_lonlat_box + Vector2(val, 0); - - // Distance between centers of boxes. The smaller it is, the better. - double dist = norm_2((shifted_src.min() + shifted_src.max())/2.0 - - (dst_lonlat_box.min() + dst_lonlat_box.max())/2.0); - dist_vec.push_back(dist); - } - int min_index = std::distance(dist_vec.begin(), min_element(dist_vec.begin(), dist_vec.end())); - m_offset = Vector2(shift[min_index], 0); - } - - //std::cout << "Geotransform: m_offset = " << m_offset << std::endl; - m_offset = Vector2(0, 0); // Disable offset code! + m_src_georef(src_georef), m_dst_georef(dst_georef) { const std::string src_datum = m_src_georef.datum().proj4_str(); const std::string dst_datum = m_dst_georef.datum().proj4_str(); @@ -85,14 +49,6 @@ namespace cartography { else m_skip_map_projection = false; - // Bugfix: If one of the projections is lon-lat, and the offset is 360 degrees, - // don't skip map-projection, as have to correct for the offset. - bool has_lonat = (m_src_georef.proj4_str().find("+proj=longlat") != std::string::npos || - m_dst_georef.proj4_str().find("+proj=longlat") != std::string::npos ); - if (has_lonat && m_offset != Vector2()) { - m_skip_map_projection = false; - } - // This optimizes the case where the two datums are the same, // and thus we don't need to call proj to convert between them // as we transform. @@ -119,14 +75,11 @@ namespace cartography { set_tolerance( 0.1 ); } - void GeoTransform::set_offset(Vector2 const& offset){ - m_offset = offset; - } Vector2 GeoTransform::reverse(Vector2 const& v) const { if (m_skip_map_projection) return m_src_georef.point_to_pixel(m_dst_georef.pixel_to_point(v)); - Vector2 src_lonlat = m_dst_georef.pixel_to_lonlat(v) - m_offset; + Vector2 src_lonlat = m_dst_georef.pixel_to_lonlat(v); if(m_skip_datum_conversion) return m_src_georef.lonlat_to_pixel(src_lonlat); src_lonlat = lonlat_to_lonlat(src_lonlat, false); @@ -137,7 +90,7 @@ namespace cartography { Vector2 GeoTransform::forward(Vector2 const& v) const { if (m_skip_map_projection) return m_dst_georef.point_to_pixel(m_src_georef.pixel_to_point(v)); - Vector2 dst_lonlat = m_src_georef.pixel_to_lonlat(v) + m_offset; + Vector2 dst_lonlat = m_src_georef.pixel_to_lonlat(v); if(m_skip_datum_conversion) return m_dst_georef.lonlat_to_pixel(dst_lonlat); dst_lonlat = lonlat_to_lonlat(dst_lonlat, true); @@ -192,7 +145,7 @@ namespace cartography { Vector2 GeoTransform::pixel_to_point( Vector2 const& v ) const { if (m_skip_map_projection) return m_src_georef.pixel_to_point(v); - Vector2 src_lonlat = m_src_georef.pixel_to_lonlat(v) + m_offset; + Vector2 src_lonlat = m_src_georef.pixel_to_lonlat(v); if(m_skip_datum_conversion) return m_dst_georef.lonlat_to_point(src_lonlat); Vector2 dst_lonlat = lonlat_to_lonlat(src_lonlat, true); @@ -202,7 +155,7 @@ namespace cartography { Vector2 GeoTransform::point_to_pixel( Vector2 const& v ) const { if (m_skip_map_projection) return m_src_georef.point_to_pixel(v); - Vector2 src_lonlat = m_src_georef.point_to_lonlat(v) + m_offset; + Vector2 src_lonlat = m_src_georef.point_to_lonlat(v); if(m_skip_datum_conversion) return m_dst_georef.lonlat_to_pixel(src_lonlat); Vector2 dst_lonlat = lonlat_to_lonlat(src_lonlat, true); diff --git a/src/vw/Cartography/GeoTransform.h b/src/vw/Cartography/GeoTransform.h index afcf5553f..e667629dd 100644 --- a/src/vw/Cartography/GeoTransform.h +++ b/src/vw/Cartography/GeoTransform.h @@ -44,7 +44,6 @@ namespace cartography { m_dst_datum_proj; bool m_skip_map_projection; bool m_skip_datum_conversion; - Vector2 m_offset; // TODO delete public: /// Normal constructor @@ -102,10 +101,6 @@ namespace cartography { /// Convert a point bounding box in the source to a pixel bounding box in the destination. BBox2 point_to_pixel_bbox( BBox2 const& point_bbox ) const; - - /// Sometimes the offset is not computed correctly, for example, when one - /// of the two georeferenced images encompasses the whole globe. - void set_offset(Vector2 const& offset); }; // End class GeoTransform diff --git a/src/vw/Cartography/tests/TestGeoReferenceUtils.cxx b/src/vw/Cartography/tests/TestGeoReferenceUtils.cxx index c89ce07f7..827127693 100644 --- a/src/vw/Cartography/tests/TestGeoReferenceUtils.cxx +++ b/src/vw/Cartography/tests/TestGeoReferenceUtils.cxx @@ -97,7 +97,7 @@ TEST( Common, gdal_write_checks ) { double nodata = -1000; bool has_nodata = true, has_georef = true; TerminalProgressCallback tpc("vw", ""); - BaseOptions opt; + GdalWriteOptions opt; block_write_gdal_image("img1.tif", dem, has_georef, georef, has_nodata, nodata, opt, tpc); block_write_gdal_image("img2.tif", dem, has_georef, georef, has_nodata, nodata, opt, tpc); diff --git a/src/vw/Cartography/tests/TestGeoTransform.cxx b/src/vw/Cartography/tests/TestGeoTransform.cxx index 021ce830e..a841d3e9d 100644 --- a/src/vw/Cartography/tests/TestGeoTransform.cxx +++ b/src/vw/Cartography/tests/TestGeoTransform.cxx @@ -147,7 +147,7 @@ TEST(GeoTransform, RefToRef) { //std::cout << "pixel1 = " << pixel1 << std::endl; //std::cout << "point1 = " << point1 << std::endl; - Vector2 lonlat = georef1.point_to_lonlat(point1); + //Vector2 lonlat = georef1.point_to_lonlat(point1); //Vector2 point2 = georef2.lonlat_to_point(lonlat); //Vector2 pixel2 = georef2.point_to_pixel(point2); //std::cout << "lonlat = " << lonlat << std::endl; diff --git a/src/vw/tools/image2qtree.h b/src/vw/tools/image2qtree.h index 70eab5121..28d4551b2 100644 --- a/src/vw/tools/image2qtree.h +++ b/src/vw/tools/image2qtree.h @@ -274,10 +274,6 @@ void do_mosaic(const Options& opt, const vw::ProgressCallback *progress) { boost::shared_ptr file( DiskImageResource::open(filename) ); GeoTransform geotx( input_georef, output_georef ); - // Even though the output georeference starts at -180, and input georef may start - // close to 180, here we don't want to correct for this discrepancy. - geotx.set_offset(Vector2(0, 0)); - ImageViewRef source = DiskImageView( file ); // Handle nodata values/mask