From 125ff464801b71f9b1dc83bbb39a458408dc5c82 Mon Sep 17 00:00:00 2001
From: Robert Butora <robert.butora@inaf.it>
Date: Thu, 12 Jun 2025 18:03:29 +0300
Subject: [PATCH] cutout: fixes double->long conversion for pixel-grid in
overlap calculation on domain=GRID
---
.../engine/src/common/include/ast4vl.hpp | 13 +---
data-access/engine/src/common/src/ast4vl.cpp | 64 +------------------
.../engine/src/common/src/ast_frameset.cpp | 27 ++++----
.../engine/src/common/src/ast_frameset.hpp | 4 +-
data-access/engine/src/common/src/cutout.cpp | 5 +-
data-access/engine/src/vlkb/src/ast.cpp | 4 +-
6 files changed, 27 insertions(+), 90 deletions(-)
diff --git a/data-access/engine/src/common/include/ast4vl.hpp b/data-access/engine/src/common/include/ast4vl.hpp
index 7c40f8a..d297a2a 100644
--- a/data-access/engine/src/common/include/ast4vl.hpp
+++ b/data-access/engine/src/common/include/ast4vl.hpp
@@ -5,7 +5,7 @@
#include <string>
#include <vector>
-#include "cutout.hpp" // coordinate needed
+#include "cutout.hpp" // uint_bounds, coordinate needed
struct point2d {double lon; double lat;};
@@ -21,17 +21,10 @@ struct Bounds
int naxis;
};
-struct double_xy
-{
- double x;
- double y;
- /*unsigned*/ char type;
-};
-
struct overlap_ranges
{
int ov_code;
- std::vector<double_xy> pixel_ranges;
+ std::vector<uint_bounds> pixel_ranges;
};
@@ -47,7 +40,7 @@ std::vector<point2d> calc_skyvertices(std::string header, std::string skysys);
std::vector<Bounds> calc_bounds(std::string header, std::string skysys_str, std::string specsys_str);
-std::vector<uint_bounds> calc_overlap(const std::string header, const coordinates coord, int& ov_code);
+overlap_ranges calc_overlap(const std::string header, const coordinates coord);
int header_modif_coordsys(const std::string& skysys_str, const std::string& specsys_str,
const std::string& pathname);
diff --git a/data-access/engine/src/common/src/ast4vl.cpp b/data-access/engine/src/common/src/ast4vl.cpp
index eb8a0ba..484c367 100644
--- a/data-access/engine/src/common/src/ast4vl.cpp
+++ b/data-access/engine/src/common/src/ast4vl.cpp
@@ -38,7 +38,7 @@ std::ostream& operator<<( std::ostream &out, struct uint_bounds const& p)
std::ostream& operator<<( std::ostream &out, overlap_ranges const& p)
{
out << p.ov_code;
- for(double_xy r : p.pixel_ranges) out <<" ("<< r.x << ", " << r.y << ")";
+ for(uint_bounds r : p.pixel_ranges) out <<" ("<< r.pix1 << ", " << r.pix2 << ")";
return out;
}
@@ -106,7 +106,7 @@ std::vector<Bounds> calc_bounds(std::string header, std::string skysys, std::str
*/
-std::vector<uint_bounds> calc_overlap(const std::string header, const coordinates coord, int& ov_code)
+overlap_ranges calc_overlap(const std::string header, const coordinates coord)
{
LOG_trace(__func__);
@@ -124,67 +124,9 @@ std::vector<uint_bounds> calc_overlap(const std::string header, const coordinate
overlap_ranges pix_ranges = frm_set.overlap(coord);
- ov_code = pix_ranges.ov_code;
-
LOG_STREAM << "ov-code & pix ranges[double]: " << pix_ranges << endl;
- /* convert to uint */
-
- vector<uint_bounds> uint_bounds_vec;
-
- LOG_STREAM << "uint_bounds: " << uint_bounds_vec.size()
- << " pix_ranges:: " << pix_ranges.pixel_ranges.size()
- << endl;
-
-
- LOG_STREAM << "pix ranges[uint]:";
- int ix = 0;
- for(double_xy dbl_range : pix_ranges.pixel_ranges)
- {
- if(dbl_range.x < 0)
- throw out_of_range(string{__FILE__} + ":" + to_string(__LINE__)
- + " pixel axis from overlap x is negative " + to_string(dbl_range.x));
-
- if(dbl_range.y < 0)
- throw out_of_range(string{__FILE__} + ":" + to_string(__LINE__)
- + " pixel axis from overlap y is negative " + to_string(dbl_range.y));
-
- // conversion double -> uint: result must be: 1 <= result <= NAXISn
- // because NAXISn start with 1 (FITS standard) which corresponds to ASTlib's GRID-domain
- // FitsChan uses GRID Domain for FITS-pixel coords
- // in GRID first element represents a pixel-range: [0.5..1.5) (0.5 in pixel_1 ; 1.5 in pixel_2)
- //
- // dbl-range is checked to be within 0.5 .. NAXIS+0.5 in ast_frameset::overlap()
- if(dbl_range.x <= dbl_range.y)
- {
- long lowb = lround(dbl_range.x);
- long upb = lround(dbl_range.y);
- long NAXISi = frm_set.get_naxis(ix++);
-
- if(lowb == (NAXISi + 1)) lowb--;
- if(upb == (NAXISi + 1)) upb--;
-
- uint_bounds ui_range{lowb, upb, dbl_range.type};
- uint_bounds_vec.push_back(ui_range);
- LOG_STREAM << " " << ui_range;
- }
- else
- {
- long lowb = lround(dbl_range.y);
- long upb = lround(dbl_range.x);
- long NAXISi = frm_set.get_naxis(ix++);
-
- if(lowb == (NAXISi + 1)) lowb--;
- if(upb == (NAXISi + 1)) upb--;
-
- uint_bounds ui_range{lowb, upb, dbl_range.type};
- uint_bounds_vec.push_back(ui_range);
- LOG_STREAM << "s" << ui_range;
- }
- }
- LOG_STREAM << endl;
-
- return uint_bounds_vec;
+ return pix_ranges;
}
diff --git a/data-access/engine/src/common/src/ast_frameset.cpp b/data-access/engine/src/common/src/ast_frameset.cpp
index 0c4064d..f887a3c 100644
--- a/data-access/engine/src/common/src/ast_frameset.cpp
+++ b/data-access/engine/src/common/src/ast_frameset.cpp
@@ -147,11 +147,6 @@ ast::frameset::frameset(string header)
}
-long ast::frameset::get_naxis(int ix)
-{
- return m_NAXISn[ix];
-}
-
bool is_specdomain(AstFrame * frm)
{
@@ -882,8 +877,12 @@ void ast::frameset::set_bounds_to_rect(coordinates& coord)
}
}
-
-
+long lround_grid(double num, long NAXISi)
+{
+ long lnum = lround(num);
+ if(lnum == (NAXISi+1)) lnum--;
+ return lnum;
+}
/* calcs overlap in pixel coords */
overlap_ranges ast::frameset::overlap(coordinates coord)
@@ -994,7 +993,6 @@ overlap_ranges ast::frameset::overlap(coordinates coord)
up[ix_lat] = D2R * (coord.lat_deg + coord.dlat_deg/2.0);
}
- vector<double_xy> pix_ranges;
if(m_has_specframe && (coord.specsys != specsystem::NONE))
{
@@ -1020,7 +1018,9 @@ overlap_ranges ast::frameset::overlap(coordinates coord)
up [m_spec_axis-1] = coord.vu_kmps;
log_bounds("BOUNDS client WCS: ", naxes, low, up);
LOG_STREAM << "BOUNDS no overlap in spectrum axis, returning ov_code=1" << endl;
- return overlap_ranges{1, pix_ranges };
+
+ vector<uint_bounds> empty_pix_ranges;
+ return overlap_ranges{1, empty_pix_ranges };
}
else // at least partial overlap -> cut coord to min/max of header values
{
@@ -1079,7 +1079,7 @@ overlap_ranges ast::frameset::overlap(coordinates coord)
/* if overlap, calc pixel ranges */
-
+ vector<uint_bounds> pix_ranges;
bool no_overlap = ((ov_code == 1) || (ov_code == 6));
if(!no_overlap)
{
@@ -1118,17 +1118,18 @@ overlap_ranges ast::frameset::overlap(coordinates coord)
+ to_string(ix+1) + "] + 0.5 = " + to_string(p1[ix])
+ ": " + to_string(ubnd[ix]));
- pix_ranges.push_back({lbnd[ix],ubnd[ix],m_axis_type[ix]});
+ pix_ranges.push_back( { lround_grid(lbnd[ix],m_NAXISn[ix]),
+ lround_grid(ubnd[ix],m_NAXISn[ix]),
+ m_axis_type[ix]} );
}
}
- return overlap_ranges{ov_code, pix_ranges };
+ return overlap_ranges{ ov_code, pix_ranges };
}
#endif
-
bool almost_equal(double x, double y, int ulp)
{
// the machine epsilon has to be scaled to the magnitude of the values used
diff --git a/data-access/engine/src/common/src/ast_frameset.hpp b/data-access/engine/src/common/src/ast_frameset.hpp
index aa8c111..c4b86fb 100644
--- a/data-access/engine/src/common/src/ast_frameset.hpp
+++ b/data-access/engine/src/common/src/ast_frameset.hpp
@@ -7,7 +7,7 @@ extern "C" {
}
#include "cutout.hpp" // coordinates needed
-#include "ast4vl.hpp" // struct decls needed
+#include "ast4vl.hpp" // overlap_ranges, struct decls needed
#include <array>
#include <vector>
@@ -24,8 +24,6 @@ class frameset
frameset(std::string header);
~frameset();
- long get_naxis(int ix);
-
bool has_specframe(void);
bool has_timeaxis(void);
diff --git a/data-access/engine/src/common/src/cutout.cpp b/data-access/engine/src/common/src/cutout.cpp
index f32256d..77f1c7c 100644
--- a/data-access/engine/src/common/src/cutout.cpp
+++ b/data-access/engine/src/common/src/cutout.cpp
@@ -304,8 +304,9 @@ std::uintmax_t cutout_file(
}
int ov_code;
- vector<uint_bounds> bounds = calc_overlap(header, coord, ov_code);
- //vector<uint_bounds> bounds = legacy::call_AST4VL_overlap(coord, header, ov_code);
+ overlap_ranges ov_ranges = calc_overlap(header, coord);
+ ov_code = ov_ranges.ov_code;
+ vector<uint_bounds> bounds = ov_ranges.pixel_ranges;
if((ov_code==1)||(ov_code==6))
throw invalid_argument("given coordinates do not overlap with given fits-hdu area (ov_code = " + to_string(ov_code) +")");
diff --git a/data-access/engine/src/vlkb/src/ast.cpp b/data-access/engine/src/vlkb/src/ast.cpp
index deb0a4d..e176da7 100644
--- a/data-access/engine/src/vlkb/src/ast.cpp
+++ b/data-access/engine/src/vlkb/src/ast.cpp
@@ -191,7 +191,9 @@ int vlkb_overlap(const string& pathname, const string& region, vector<uint_bound
for(unsigned int i=0; i<allHdus.size(); i++)
{
fitsfiles::Hdu hd = allHdus.at(i);
- bnds = calc_overlap(hd.m_header, coord, ov_code);
+ overlap_ranges ov_bnds = calc_overlap(hd.m_header, coord);
+ ov_code = ov_bnds.ov_code;
+ bnds = ov_bnds.pixel_ranges;
}
update_with_pixel_i(bnds, coord.pixel_i);
return ov_code;
--
GitLab