Files
OpenSpace/modules/solarbrowsing/util/spacecraftimagerymanager.cpp
T
Andreas Engberg febd298556 Update Solarbrowsing branch to compile on latest master
Also includes some additional minor changes from the original feature/solarbrowsing branch
2025-12-03 10:48:30 +01:00

607 lines
22 KiB
C++

/*****************************************************************************************
* *
* OpenSpace *
* *
* Copyright (c) 2014-2018 *
* *
* Permission is hereby granted, free of charge, to any person obtaining a copy of this *
* software and associated documentation files (the "Software"), to deal in the Software *
* without restriction, including without limitation the rights to use, copy, modify, *
* merge, publish, distribute, sublicense, and/or sell copies of the Software, and to *
* permit persons to whom the Software is furnished to do so, subject to the following *
* conditions: *
* *
* The above copyright notice and this permission notice shall be included in all copies *
* or substantial portions of the Software. *
* *
* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, *
* INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A *
* PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT *
* HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF *
* CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE *
* OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. *
****************************************************************************************/
#include <modules/solarbrowsing/util/spacecraftimagerymanager.h>
#include <openspace/engine/globals.h>
#include <openspace/engine/openspaceengine.h>
#include <openspace/rendering/transferfunction.h>
#include <openspace/util/spicemanager.h>
#include <openspace/util/timemanager.h>
#include <ghoul/filesystem/filesystem.h>
#include <ghoul/misc/threadpool.h>
#include <ghoul/opengl/texture.h>
#include <ghoul/logging/logmanager.h>
#include <ghoul/filesystem/file.h>
#include <nlohmann/json.hpp>
#include <chrono>
#include <format>
#include <fstream>
#include <future>
#include <string>
#include <iostream>
#include <optional>
namespace {
constexpr const char* _loggerCat = "SpacecraftImageryManager";
constexpr const double SUN_RADIUS = 1391600000.0 * 0.5;
} // namespace
namespace openspace {
// Conversion needed before passing dates into the spice manager
std::string SpacecraftImageryManager::ISO8601(std::string& datetime) {
std::string month = datetime.substr(5, 3);
std::string MM = "";
if (month == "JAN") MM = "01";
else if (month == "FEB") MM = "02";
else if (month == "MAR") MM = "03";
else if (month == "APR") MM = "04";
else if (month == "MAY") MM = "05";
else if (month == "JUN") MM = "06";
else if (month == "JUL") MM = "07";
else if (month == "AUG") MM = "08";
else if (month == "SEP") MM = "09";
else if (month == "OCT") MM = "10";
else if (month == "NOV") MM = "11";
else if (month == "DEC") MM = "12";
else ghoul_assert(false, "Bad month");
datetime.replace(4, 5, "-" + MM + "-");
return datetime;
}
void SpacecraftImageryManager::loadTransferFunctions(const std::string& path,
std::unordered_map<std::string, std::shared_ptr<TransferFunction>>& _tfMap)
{
std::filesystem::path sequenceDir(path);
if (!std::filesystem::is_directory(sequenceDir)) {
LERROR(std::format("Could not load directory '{}'", sequenceDir.string()));
}
std::vector<std::filesystem::path> sequencePaths = ghoul::filesystem::walkDirectory(
sequenceDir,
ghoul::filesystem::Recursive::Yes,
ghoul::filesystem::Sorted::Yes
);
for (const std::filesystem::path& seqPath : sequencePaths) {
if (seqPath.extension() == ".txt") {
std::string key = seqPath.stem().string();
_tfMap[key] = std::make_shared<TransferFunction>(seqPath);
}
}
}
bool SpacecraftImageryManager::loadMetadataFromDisk(const std::string& rootPath,
std::unordered_map<std::string, Timeline<ImageMetadata>>& _imageMetadataMap)
{
std::filesystem::path sequenceDir(rootPath);
if (!std::filesystem::is_directory(sequenceDir)) {
LERROR(std::format("Could not load directory '{}'", sequenceDir.string()));
}
bool metadataLoaded = false;
std::vector<std::filesystem::path> sequencePaths = ghoul::filesystem::walkDirectory(
sequenceDir,
ghoul::filesystem::Recursive::No,
ghoul::filesystem::Sorted::No
);
for (const std::filesystem::path& seqPath : sequencePaths) {
const std::string extension = seqPath.extension().string();
const std::string base = seqPath.filename().string();
const std::size_t foundCachedImageData = base.find("_cached");
if (extension != ".txt" || foundCachedImageData == std::string::npos) {
continue;
}
const std::string::size_type separator = base.rfind("_");
const std::string instrument = base.substr(0, separator);
LDEBUG(std::format("Loading instrument: {}", instrument));
metadataLoaded = true;
std::ifstream myfile(seqPath);
if (!myfile.is_open()) {
LERROR("Failed to open metadata file");
return false;
}
int numStates;
myfile >> numStates;
for (int i = 0; i < numStates; i++) {
ImageMetadata im;
myfile >> std::ws; // Skip the rest of the line
std::string date;
std::getline(myfile, date);
if (date.empty()) {
LERROR("Failed to open state metadata: date");
return false;
}
double timeObserved =
SpiceManager::ref().ephemerisTimeFromDate(ISO8601(date));
std::string relPath;
myfile >> relPath;
if (myfile.bad()) {
LERROR("Failed to open state metadata: relPath");
return false;
}
im.filename = rootPath + "/" + relPath;
myfile >> im.fullResolution;
if (myfile.bad()) {
LERROR("Failed to open state metadata: fullResolution");
return false;
}
myfile >> im.scale;
if (myfile.bad()) {
LERROR("Failed to open state metadata: scale");
return false;
}
float x, y;
myfile >> x >> y;
im.centerPixel = glm::vec2(x,y);
myfile >> im.isCoronaGraph;
if (myfile.bad()) {
LERROR("Failed to open state metadata: isCoronaGraph");
return false;
}
_imageMetadataMap[instrument].addKeyframe(timeObserved, std::move(im));
}
myfile.close();
}
return metadataLoaded;
}
void SpacecraftImageryManager::saveMetadataToDisk(const std::string& rootPath,
std::unordered_map<std::string, Timeline<ImageMetadata>>& _imageMetadataMap)
{
using K = std::string;
using V = Timeline<ImageMetadata>;
for (const std::pair<K, V>& instrument : _imageMetadataMap) {
std::ofstream ofs(rootPath + "/" + instrument.first + "_cached" + ".txt");
if (!ofs.is_open()) {
LERROR("Failed to open file");
}
const Timeline<ImageMetadata>& sequence = instrument.second;
ofs << sequence.nKeyframes() << "\n";
for (const Keyframe<ImageMetadata>& metadata : sequence.keyframes()) {
ofs << SpiceManager::ref().dateFromEphemerisTime(
metadata.timestamp
) << "\n";
const ImageMetadata& im = metadata.data;
std::filesystem::path filePath(im.filename);
std::filesystem::path fileName = filePath.filename();
ofs << fileName.string() << "\n";
ofs << im.fullResolution << "\n";
ofs << im.scale << "\n";
ofs << im.centerPixel.x << "\n";
ofs << im.centerPixel.y << "\n";
ofs << im.isCoronaGraph << "\n";
}
ofs.close();
}
}
// @TODO emiax: If openjpeg ever starts supporting reading XML metadata,
// this implementation should be improved in order not to search the entire buffer for
// XML data. There is an issue here:
// (https://github.com/uclouvain/openjpeg/issues/929)
ImageMetadata SpacecraftImageryManager::parseJ2kMetadata(
const ghoul::filesystem::File& file)
{
ImageMetadata im;
im.filename = file.path().string();
std::ifstream stream(file.path(), std::ios::binary | std::ios::ate);
std::streamsize size = stream.tellg();
stream.seekg(0, std::ios::beg);
std::vector<char> buffer(size);
if (!stream.read(buffer.data(), size)) {
LERROR(std::format("Failed to read data from '{}' ", file.path()));
return im;
}
std::string_view bufferView(buffer.data(), size);
auto extractInnerXml = [](std::string_view view, const std::string& elementName) ->
std::optional<std::string_view>
{
const std::string startTag =
std::string("<") + elementName + std::string(">");
const std::string endTag = std::string("</") + elementName + std::string(">");
const auto begin =
std::search(view.begin(), view.end(), startTag.begin(), startTag.end());
if (begin == view.end()) {
return std::nullopt;
; }
const auto afterBeginTag = begin + startTag.size();
const auto end =
std::search(afterBeginTag, view.end(), endTag.begin(), endTag.end());
if (end == view.end()) {
return std::nullopt;
}
return std::string_view(&*afterBeginTag, end - afterBeginTag);
};
std::optional<std::string_view> metaData = extractInnerXml(bufferView, "meta");
if (!metaData.has_value()) {
LERROR(std::format("Could not find metadata in {}", file.path()));
return im;
}
std::optional<std::string_view> telescop =
extractInnerXml(metaData.value(), "TELESCOP");
if (!telescop.has_value()) {
LERROR(std::format("Could not find TELESCOP tag {}", file.path()));
return im;
}
if (std::optional<std::string_view> naxis =
extractInnerXml(metaData.value(), "NAXIS1"))
{
im.fullResolution = std::stoi(std::string(naxis.value()));
}
else {
LERROR(std::format("Could not find NAXIS1 tag {}", file.path()));
return im;
}
const float halfRes = im.fullResolution / 2.f;
glm::vec2 centerPixel;
if (std::optional<std::string_view> centerPixelX =
extractInnerXml(bufferView, "CRPIX1"))
{
centerPixel.x = std::stof(std::string(centerPixelX.value()));
}
else {
LERROR(std::format("Could not find CRPIX1 tag {}", file.path()));
return im;
}
if (std::optional<std::string_view> centerPixelY =
extractInnerXml(bufferView, "CRPIX2"))
{
centerPixel.y = std::stof(std::string(centerPixelY.value()));
}
else {
LERROR(std::format("Could not find CRPIX2 tag {}", file.path()));
return im;
}
const glm::vec2 offset = ((halfRes - centerPixel) / halfRes) * glm::vec2(SUN_RADIUS);
im.centerPixel = offset;
if (telescop.value() == "SOHO") {
if (std::optional<std::string_view> plateScl =
extractInnerXml(metaData.value(), "PLATESCL"))
{
const float plateScale = std::stof(std::string(plateScl.value()));
im.scale = 1.f / (plateScale / 2.f);
im.isCoronaGraph = true;
}
else {
LERROR(std::format("Could not find NAXIS1 tag {}", file.path()));
return im;
}
}
else if (telescop.value() == "SDO") {
std::optional<std::string_view> rsunObs = extractInnerXml(bufferView, "RSUN_OBS");
std::optional<std::string_view> cDelt1 = extractInnerXml(bufferView, "CDELT1");
if (!rsunObs.has_value()) {
LERROR(std::format("Could not find RSUN_OBS tag {}", file.path()));
return im;
}
if (!cDelt1.has_value()) {
LERROR(std::format("Could not find CDELT1 tag {}", file.path()));
return im;
}
im.scale = (std::stof(std::string(rsunObs.value())) /
std::stof(std::string(cDelt1.value()))) /
(im.fullResolution / 2.f);
im.isCoronaGraph = false;
}
else { // Telescope is assumed to be STEREO
std::optional<std::string_view> rsun = extractInnerXml(bufferView, "RSUN");
std::optional<std::string_view> cDelt1 = extractInnerXml(bufferView, "CDELT1");
if (!rsun.has_value()) {
LERROR(std::format("Could not find RSUN_OBS tag {}", file.path()));
return im;
}
if (!cDelt1.has_value()) {
LERROR(std::format("Could not find CDELT1 tag {}", file.path()));
return im;
}
im.scale = (std::stof(std::string(rsun.value())) /
std::stof(std::string(cDelt1.value()))) /
(im.fullResolution / 2.f);
im.isCoronaGraph = false;
if (std::optional<std::string_view> detector =
extractInnerXml(bufferView, "DETECTOR"))
{
im.isCoronaGraph =
detector.value() == "COR1" || detector.value() == "COR2";
}
else {
LWARNING(std::format(
"Could not find DETECTOR tag {}", file.path()
));
}
}
return im;
}
// This is currently not used. Instead, the parseJ2kMetadata is used,
// extracting the data directoy from the JPEG2000 file by naively searching the entire
// buffer for metadata, avoiding pre-processing steps.
// If you want to use this, you need to extract metadata to json first,
// for example using: https://github.com/novalain/j2kcodec
ImageMetadata SpacecraftImageryManager::parseJsonMetadata(
const ghoul::filesystem::File& file)
{
ImageMetadata im;
im.filename = file.path().string();
const std::string filename = std::string(file.path().string() + ".json");
if (!std::filesystem::exists(filename)) {
LERROR(std::format("'{}' has no specified json metadata", file.path()));
// TODO: Hardcoded values, when there are no json files.
return im;
}
// Parse JSON metadata
using json = nlohmann::json;
std::ifstream i(filename);
if (i.fail()) {
LERROR(std::format("Error opening file '{}'", filename));
}
json j;
i >> j;
const json data = j["meta"]["fits"];
if (data.is_null()) {
LERROR(std::format("Error in metadata '{}'", filename));
}
json value = data["TELESCOP"];
if (value.is_null() && !value.is_string()) {
LERROR("Metadata did contain information about type of spacecraft");
return im;
}
// TODO: value might not exist
std::string spacecraftType = value;
value = data["NAXIS1"];
if (value.is_null()) {
LERROR("Metadata did not contain information about resolution");
}
const std::string sFullResolution = value;
im.fullResolution = std::stoi(sFullResolution);
// Special case of sdo - RSUN is given in pixels
// For SOHO the radius of the sun is not specified - we instead use platescl
if (spacecraftType == "SOHO") {
const std::string sScale = data["PLATESCL"];
const float plateScale = stof(sScale);
im.scale = 1.f / (plateScale / 2.f);
im.isCoronaGraph = true;
} else {
float sunRadiusPixels = 0.f;
// SDO has RSUN specified in pixels
if (spacecraftType == "SDO") {
value = data["RSUN_OBS"];
if (value.is_null()) {
LERROR("SDO Metadata: RSUN_OBS missing!");
}
std::string sSunRadiusArcsec = value;
value = data["CDELT1"];
if (value.is_null()) {
LERROR("SDO Metadata: CDELT1 missing!");
}
std::string sCdelt1 = value;
const float cdelt1 = stof(sCdelt1);
const float sunRadiusArcsec = stof(sSunRadiusArcsec);
sunRadiusPixels = sunRadiusArcsec / cdelt1;
im.isCoronaGraph = false;
} else {
// STEREO has RSUN specified in arcsecs - need to divide by factor
std::string sCdelt1 = data["CDELT1"];
const float cdelt1 = stof(sCdelt1);
std::string sSunRadiusArcsec = data["RSUN"];
const float sunRadiusArcsec = stof(sSunRadiusArcsec);
sunRadiusPixels = sunRadiusArcsec / cdelt1;
value = data["DETECTOR"];
if (value.is_null()) {
LERROR("No observer specified in Stereo");
}
std::string sObserver = value;
if (sObserver == "COR1" || sObserver == "COR2") {
im.isCoronaGraph = true;
} else {
im.isCoronaGraph = false;
}
}
float scale = sunRadiusPixels / (im.fullResolution / 2.f);
im.scale = scale;
}
const std::string sCenterPixelX = data["CRPIX1"];
const std::string sCenterPixelY = data["CRPIX2"];
const float centerPixelX = stof(sCenterPixelX);
const float centerPixelY = stof(sCenterPixelY);
const float halfRes = im.fullResolution / 2.f;
const float offsetX = ((halfRes - centerPixelX) / halfRes) * SUN_RADIUS;
const float offsetY = ((halfRes - centerPixelY) / halfRes) * SUN_RADIUS;
im.centerPixel = glm::vec2(offsetX, offsetY);
return im;
}
void SpacecraftImageryManager::loadImageMetadata(const std::string& path,
std::unordered_map<std::string, Timeline<ImageMetadata>>& _imageMetadataMap)
{
std::unordered_map<std::string, Timeline<ImageMetadata>> result;
LDEBUG("Begin loading spacecraft imagery metadata");
// Pre-processed data
if (loadMetadataFromDisk(path, result)) {
_imageMetadataMap.insert(result.begin(), result.end());
return;
}
std::filesystem::path sequenceDir(path);
if (!std::filesystem::is_directory(sequenceDir)) {
LERROR(std::format("Could not load directory '{}'", sequenceDir.string()));
}
LDEBUG("Loading sequence directory");
std::vector<std::filesystem::path> sequencePaths = ghoul::filesystem::walkDirectory(
sequenceDir,
ghoul::filesystem::Recursive::Yes,
ghoul::filesystem::Sorted::Yes
);
LDEBUG("Filtering data values");
sequencePaths.erase(
std::remove_if(
sequencePaths.begin(),
sequencePaths.end(),
[](const std::filesystem::path& path) {
const std::string& ext = path.extension().string();
return (ext != ".jp2") && (ext != ".j2k");
}
),
sequencePaths.end()
);
LDEBUG("Reading meta data");
size_t count = 0;
std::mutex spiceAndPushMutex;
std::vector<std::future<void>> futures;
futures.reserve(sequencePaths.size());
std::cout << '\n';
auto exec = [&](const std::filesystem::path& seqPath) {
ghoul::filesystem::File currentFile(seqPath);
std::string fileName = seqPath.filename().string();
size_t posSatelliteInfoStart = fileName.rfind("__") + 2;
std::string satelliteInfo = fileName.substr(posSatelliteInfoStart);
// Name
size_t posSatelliteNameEnd = satelliteInfo.find_first_of("_");
std::string satelliteName = satelliteInfo.substr(0, posSatelliteNameEnd);
// Instrument
size_t posInstrumentNameStart = posSatelliteNameEnd + 1;
std::string instrumentName = satelliteInfo.substr(posInstrumentNameStart);
size_t dot = instrumentName.rfind(".");
instrumentName = instrumentName.substr(0, dot);
// Time @TODO (anden88 2025-12-02): Parse the time string from file in a more
// robust way?
std::vector<std::string> tokens;
std::stringstream ss;
ss.str(seqPath.filename().string());
std::string item;
while (std::getline(ss, item, '_')) {
tokens.push_back(item);
}
if (tokens.size() >= 8) {
std::string time = tokens[0] + "-" + tokens[1] + "-" +
tokens[2] + "T" + tokens[4] + ":" +
tokens[5] + ":" + tokens[6] + "." + tokens[7];
ImageMetadata im = parseJ2kMetadata(currentFile);
spiceAndPushMutex.lock();
result[instrumentName].addKeyframe(
global::timeManager->time().convertTime(time),
std::move(im)
);
spiceAndPushMutex.unlock();
}
++count;
if (count % 10000 == 0) {
LINFO(std::format(
"Processing image {} out of {} ", count, sequencePaths.size()
));
}
};
for (const std::filesystem::path& seqPath : sequencePaths) {
exec(seqPath);
}
saveMetadataToDisk(path, result);
_imageMetadataMap.insert(result.begin(), result.end());
LDEBUG("Finish loading imagery metadata");
LDEBUG("Saving imagery metadata");
LDEBUG(std::format("{} images loaded", count));
LDEBUG(std::format("{} values in metamap", _imageMetadataMap.size()));
}
} //namespace openspace