Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Minor refactoring of polarization values and defaults in the slicer function #1390

Merged
merged 4 commits into from
Nov 6, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
117 changes: 54 additions & 63 deletions src/Frame/Frame.cc
Original file line number Diff line number Diff line change
Expand Up @@ -85,6 +85,10 @@ Frame::Frame(uint32_t session_id, std::shared_ptr<FileLoader> loader, const std:
_depth = (_z_axis >= 0 ? _image_shape(_z_axis) : 1);
_num_stokes = (_stokes_axis >= 0 ? _image_shape(_stokes_axis) : 1);

_all_x = AxisRange(0, _width - 1);
_all_y = AxisRange(0, _height - 1);
_all_z = AxisRange(0, _depth - 1);

_use_tile_cache = _loader->UseTileCache();

// load full image cache for loaders that don't use the tile cache and mipmaps
Expand Down Expand Up @@ -207,7 +211,7 @@ bool Frame::GetBeams(std::vector<CARTA::Beam>& beams) {
}

StokesSlicer Frame::GetImageSlicer(const AxisRange& z_range, int stokes) {
return GetImageSlicer(AxisRange(ALL_X), AxisRange(ALL_Y), z_range, stokes);
return GetImageSlicer(_all_x, _all_y, z_range, stokes);
}

StokesSlicer Frame::GetImageSlicer(const AxisRange& x_range, const AxisRange& y_range, const AxisRange& z_range, int stokes) {
Expand All @@ -226,13 +230,6 @@ StokesSlicer Frame::GetImageSlicer(const AxisRange& x_range, const AxisRange& y_
int start_x(x_range.from), end_x(x_range.to);

// Normalize x constants
if (start_x == ALL_X) {
start_x = 0;
}
if (end_x == ALL_X) {
end_x = _width - 1;
}

if (stokes_source.IsOriginalImage()) {
start(_x_axis) = start_x;
end(_x_axis) = end_x;
Expand All @@ -247,13 +244,6 @@ StokesSlicer Frame::GetImageSlicer(const AxisRange& x_range, const AxisRange& y_
int start_y(y_range.from), end_y(y_range.to);

// Normalize y constants
if (start_y == ALL_Y) {
start_y = 0;
}
if (end_y == ALL_Y) {
end_y = _height - 1;
}

if (stokes_source.IsOriginalImage()) {
start(_y_axis) = start_y;
end(_y_axis) = end_y;
Expand Down Expand Up @@ -313,7 +303,7 @@ bool Frame::CheckZ(int z) {
}

bool Frame::CheckStokes(int stokes) {
return (((stokes >= 0) && (stokes < NumStokes())) || IsComputedStokes(stokes));
return (((stokes >= 0) && (stokes < NumStokes())) || Stokes::IsComputed(stokes));
}

bool Frame::ZStokesChanged(int z, int stokes) {
Expand Down Expand Up @@ -349,7 +339,7 @@ bool Frame::SetImageChannels(int new_z, int new_stokes, std::string& message) {
_z_index = new_z;
_stokes_index = new_stokes;

if (!(_use_tile_cache && _loader->HasMip(2)) || IsComputedStokes(_stokes_index)) {
if (!(_use_tile_cache && _loader->HasMip(2)) || Stokes::IsComputed(_stokes_index)) {
// Reload the full channel cache for loaders which use it
FillImageCache();
} else {
Expand Down Expand Up @@ -582,7 +572,7 @@ bool Frame::GetRasterTileData(std::shared_ptr<std::vector<float>>& tile_data_ptr
tile_data_ptr = _tile_pool->Pull();
bool loaded_data(0);

if (mip > 1 && !IsComputedStokes(_stokes_index)) {
if (mip > 1 && !Stokes::IsComputed(_stokes_index)) {
// Try to load downsampled data from the image file
loaded_data = _loader->GetDownsampledRasterData(*tile_data_ptr, _z_index, _stokes_index, bounds, mip, _image_mutex);
} else if (!_image_cache_valid && _use_tile_cache) {
Expand Down Expand Up @@ -1159,7 +1149,7 @@ bool Frame::FillSpatialProfileData(PointXy point, std::vector<CARTA::SetSpatialR
bool have_profile(false);
bool downsample(mip >= 2);

if (downsample && _loader->HasMip(2) && !IsComputedStokes(stokes)) { // Use a mipmap dataset to return downsampled data
if (downsample && _loader->HasMip(2) && !Stokes::IsComputed(stokes)) { // Use a mipmap dataset to return downsampled data
while (!_loader->HasMip(mip)) {
mip /= 2;
}
Expand Down Expand Up @@ -1447,8 +1437,8 @@ bool Frame::FillSpectralProfileData(std::function<void(CARTA::SpectralProfileDat

std::vector<float> spectral_data;
int xy_count(1);
if (!IsComputedStokes(stokes) && _loader->GetCursorSpectralData(spectral_data, stokes, (start_cursor.x + 0.5), xy_count,
(start_cursor.y + 0.5), xy_count, _image_mutex)) {
if (!Stokes::IsComputed(stokes) && _loader->GetCursorSpectralData(spectral_data, stokes, (start_cursor.x + 0.5), xy_count,
(start_cursor.y + 0.5), xy_count, _image_mutex)) {
// Use loader data
spectral_profile->set_raw_values_fp32(spectral_data.data(), spectral_data.size() * sizeof(float));
cb(profile_message);
Expand Down Expand Up @@ -2329,60 +2319,61 @@ casacore::Slicer Frame::GetExportRegionSlicer(const CARTA::SaveFile& save_file_m

bool Frame::GetStokesTypeIndex(const string& coordinate, int& stokes_index) {
// Coordinate could be profile (x, y, z), stokes string (I, Q, U), or combination (Ix, Qy)
bool is_stokes_string = StokesStringTypes.find(coordinate) != StokesStringTypes.end();
bool is_combination = (coordinate.size() > 1 && (coordinate.back() == 'x' || coordinate.back() == 'y' || coordinate.back() == 'z'));

if (is_combination || is_stokes_string) {
bool stokes_ok(false);
if (coordinate.empty() || coordinate == 'x' || coordinate == 'y' || coordinate == 'z') {
// Profile only or blank; use current Stokes
stokes_index = CurrentStokes();
return true;
}

std::string stokes_string;
if (coordinate.size() > 1 && (coordinate.back() == 'x' || coordinate.back() == 'y' || coordinate.back() == 'z')) {
// Combination Stokes and profile string
stokes_string = coordinate.substr(0, coordinate.size() - 1);
} else {
// Stokes string
stokes_string = coordinate;
}

bool stokes_ok(false);

std::string stokes_string;
if (is_stokes_string) {
stokes_string = coordinate;
auto stokes_type = Stokes::Get(stokes_string);
if (stokes_type) {
if (_loader->GetStokesTypeIndex(stokes_type, stokes_index)) {
stokes_ok = true;
} else if (Stokes::IsComputed(stokes_type)) {
stokes_index = stokes_type;
stokes_ok = true;
} else {
stokes_string = coordinate.substr(0, coordinate.size() - 1);
}

if (StokesStringTypes.count(stokes_string)) {
CARTA::PolarizationType stokes_type = StokesStringTypes[stokes_string];
if (_loader->GetStokesTypeIndex(stokes_type, stokes_index)) {
stokes_ok = true;
} else if (IsComputedStokes(stokes_string)) {
stokes_index = StokesStringTypes.at(stokes_string);
int assumed_stokes_index = (stokes_type - 1) % 4;
if (NumStokes() > assumed_stokes_index) {
stokes_index = assumed_stokes_index;
stokes_ok = true;
} else {
int assumed_stokes_index = (StokesValues[stokes_type] - 1) % 4;
if (NumStokes() > assumed_stokes_index) {
stokes_index = assumed_stokes_index;
stokes_ok = true;
spdlog::warn("Can not get stokes index from the header. Assuming stokes {} index is {}.", stokes_string, stokes_index);
}
spdlog::warn("Can not get stokes index from the header. Assuming stokes {} index is {}.", stokes_string, stokes_index);
}
}
if (!stokes_ok) {
spdlog::error("Spectral or spatial requirement {} failed: invalid stokes axis for image.", coordinate);
return false;
}
} else {
stokes_index = CurrentStokes(); // current stokes
}

if (!stokes_ok) {
spdlog::error("Spectral or spatial requirement {} failed: invalid stokes axis for image.", coordinate);
return false;
}

return true;
}

std::string Frame::GetStokesType(int stokes_index) {
for (auto stokes_type : StokesStringTypes) {
int tmp_stokes_index;
if (_loader->GetStokesTypeIndex(stokes_type.second, tmp_stokes_index) && (tmp_stokes_index == stokes_index)) {
std::string stokes = (stokes_type.first.length() == 1) ? fmt::format("Stokes {}", stokes_type.first) : stokes_type.first;
return stokes;
}
}
if (IsComputedStokes(stokes_index)) {
CARTA::PolarizationType stokes_type = StokesTypes[stokes_index];
if (ComputedStokesName.count(stokes_type)) {
return ComputedStokesName[stokes_type];
}
auto stokes_type = CARTA::PolarizationType::POLARIZATION_TYPE_NONE;

// Computed stokes: stokes index is equal to numeric value
if (Stokes::IsComputed(stokes_index)) {
stokes_type = Stokes::Get(stokes_index);
}
return "Unknown";

// Otherwise try to map index to type with loader
_loader->GetStokesType(stokes_index, stokes_type);

return Stokes::Description(stokes_type);
}

std::shared_mutex& Frame::GetActiveTaskMutex() {
Expand Down
7 changes: 1 addition & 6 deletions src/Frame/Frame.h
Original file line number Diff line number Diff line change
Expand Up @@ -75,12 +75,6 @@ static std::unordered_map<CARTA::FileType, string> FileTypeString{{CARTA::FileTy
{CARTA::FileType::DS9_REG, "DS9"}, {CARTA::FileType::FITS, "FITS"}, {CARTA::FileType::HDF5, "HDF5"},
{CARTA::FileType::MIRIAD, "MIRIAD"}, {CARTA::FileType::UNKNOWN, "Unknown"}};

static std::unordered_map<CARTA::PolarizationType, std::string> ComputedStokesName{
{CARTA::PolarizationType::Ptotal, "Total polarization intensity"}, {CARTA::PolarizationType::Plinear, "Linear polarization intensity"},
{CARTA::PolarizationType::PFtotal, "Fractional total polarization intensity"},
{CARTA::PolarizationType::PFlinear, "Fractional linear polarization intensity"},
{CARTA::PolarizationType::Pangle, "Polarization angle"}};

class Frame {
public:
// Load image cache for default_z, except for PV preview image which needs cube
Expand Down Expand Up @@ -286,6 +280,7 @@ class Frame {
int _spectral_axis, _stokes_axis;
int _z_index, _stokes_index; // current index
size_t _width, _height, _depth, _num_stokes;
AxisRange _all_x, _all_y, _all_z;

// Image settings
CARTA::AddRequiredTiles _required_animation_tiles;
Expand Down
14 changes: 0 additions & 14 deletions src/ImageData/FileInfo.h
Original file line number Diff line number Diff line change
Expand Up @@ -14,7 +14,6 @@
#include <casacore/images/Images/SubImage.h>

#include <carta-protobuf/defs.pb.h>
//#include <carta-protobuf/enums.pb.h>

namespace carta {
namespace FileInfo {
Expand Down Expand Up @@ -122,19 +121,6 @@ inline casacore::uInt GetFitsHdu(const std::string& hdu) {
return hdu_num;
}

// convert between CARTA::PolarizationType values and FITS standard stokes values
static bool ConvertFitsStokesValue(const int& in_stokes_value, int& out_stokes_value) {
if (in_stokes_value >= 1 && in_stokes_value <= 4) {
out_stokes_value = in_stokes_value;
return true;
} else if ((in_stokes_value >= 4 && in_stokes_value <= 12) || (in_stokes_value <= -1 && in_stokes_value >= -8)) {
// convert between [5, 6, ..., 12] and [-1, -2, ..., -8]
out_stokes_value = -in_stokes_value + 4;
return true;
}
return false;
}

} // namespace FileInfo
} // namespace carta

Expand Down
34 changes: 23 additions & 11 deletions src/ImageData/FileLoader.cc
Original file line number Diff line number Diff line change
Expand Up @@ -310,8 +310,10 @@ bool FileLoader::FindCoordinateAxes(casacore::IPosition& shape, std::vector<int>
for (int i = 0; i < _num_stokes; ++i) {
int stokes_fits_value = _stokes_crval + (i + 1 - _stokes_crpix) * _stokes_cdelt;
int stokes_value;
if (FileInfo::ConvertFitsStokesValue(stokes_fits_value, stokes_value)) {
_stokes_indices[GetStokesType(stokes_value)] = i;
if (Stokes::ConvertFits(stokes_fits_value, stokes_value)) {
auto stokes_type = static_cast<CARTA::PolarizationType>(stokes_value);
_stokes_indices[stokes_type] = i;
_stokes_types[i] = stokes_type;
}
}
}
Expand Down Expand Up @@ -846,7 +848,7 @@ void FileLoader::LoadImageStats(bool load_percentiles) {
}

FileInfo::ImageStats& FileLoader::GetImageStats(int current_stokes, int z) {
if (!IsComputedStokes(current_stokes)) { // Note: loader cache does not support the computed stokes
if (!Stokes::IsComputed(current_stokes)) { // Note: loader cache does not support the computed stokes
return (z >= 0 ? _z_stats[current_stokes][z] : _cube_stats[current_stokes]);
}
return _empty_stats;
Expand Down Expand Up @@ -912,11 +914,21 @@ double FileLoader::CalculateBeamArea() {
}

bool FileLoader::GetStokesTypeIndex(const CARTA::PolarizationType& stokes_type, int& stokes_index) {
if (_stokes_indices.count(stokes_type)) {
stokes_index = _stokes_indices[stokes_type];
try {
stokes_index = _stokes_indices.at(stokes_type);
return true;
} catch (const std::out_of_range& e) {
return false;
}
}

bool FileLoader::GetStokesType(const int& stokes_index, CARTA::PolarizationType& stokes_type) {
try {
stokes_type = _stokes_types.at(stokes_index);
return true;
} catch (const std::out_of_range& e) {
return false;
}
return false;
}

typename FileLoader::ImageRef FileLoader::GetStokesImage(const StokesSource& stokes_source) {
Expand All @@ -929,15 +941,15 @@ typename FileLoader::ImageRef FileLoader::GetStokesImage(const StokesSource& sto
carta::PolarizationCalculator polarization_calculator(
GetImage(), AxisRange(stokes_source.z_range), AxisRange(stokes_source.x_range), AxisRange(stokes_source.y_range));

if (stokes_source.stokes == COMPUTE_STOKES_PTOTAL) {
if (stokes_source.stokes == CARTA::PolarizationType::Ptotal) {
_computed_stokes_image = polarization_calculator.ComputeTotalPolarizedIntensity();
} else if (stokes_source.stokes == COMPUTE_STOKES_PFTOTAL) {
} else if (stokes_source.stokes == CARTA::PolarizationType::PFtotal) {
_computed_stokes_image = polarization_calculator.ComputeTotalFractionalPolarizedIntensity();
} else if (stokes_source.stokes == COMPUTE_STOKES_PLINEAR) {
} else if (stokes_source.stokes == CARTA::PolarizationType::Plinear) {
_computed_stokes_image = polarization_calculator.ComputePolarizedIntensity();
} else if (stokes_source.stokes == COMPUTE_STOKES_PFLINEAR) {
} else if (stokes_source.stokes == CARTA::PolarizationType::PFlinear) {
_computed_stokes_image = polarization_calculator.ComputeFractionalPolarizedIntensity();
} else if (stokes_source.stokes == COMPUTE_STOKES_PANGLE) {
} else if (stokes_source.stokes == CARTA::PolarizationType::Pangle) {
_computed_stokes_image = polarization_calculator.ComputePolarizedAngle();
} else {
spdlog::error("Unknown computed stokes index {}", stokes_source.stokes);
Expand Down
2 changes: 2 additions & 0 deletions src/ImageData/FileLoader.h
Original file line number Diff line number Diff line change
Expand Up @@ -121,6 +121,7 @@ class FileLoader {
virtual void SetStokesCrpix(float stokes_crpix);
virtual void SetStokesCdelt(int stokes_cdelt);
virtual bool GetStokesTypeIndex(const CARTA::PolarizationType& stokes_type, int& stokes_index);
virtual bool GetStokesType(const int& stokes_index, CARTA::PolarizationType& stokes_type);
std::unordered_map<CARTA::PolarizationType, int> GetStokesIndices() {
return _stokes_indices;
};
Expand Down Expand Up @@ -180,6 +181,7 @@ class FileLoader {

// Storage for the stokes type vs. stokes index
std::unordered_map<CARTA::PolarizationType, int> _stokes_indices;
std::unordered_map<int, CARTA::PolarizationType> _stokes_types;
float _stokes_crval;
float _stokes_crpix;
int _stokes_cdelt;
Expand Down
Loading
Loading