Apply transfer function to dots on map

This commit is contained in:
Roxeena
2025-02-20 17:49:51 +01:00
parent 24a73ba38a
commit c9e8cd135b
4 changed files with 105 additions and 83 deletions
+2 -2
View File
@@ -1,13 +1,13 @@
return {
--[[{
{
Type = "FindImpactsTask",
AsteroidName = "2023 CX1",
KernelDirectory = "C:/Users/malej60/Documents/develop/b612-oribt-uncertainty/data/2023 CX1/openspace_variants/",
OutputFilename = "C:/Users/malej60/Documents/develop/b612-oribt-uncertainty/data/2023 CX1/impacts/2023-CX1-impacts.txt",
TimeIntervalStart = "2023-02-13T02:39:00.000",
TimeIntervalEnd = "2023-02-13T03:39:00.000"
},]]--
},
{
Type = "FindImpactsTask",
AsteroidName = "2004 MN4",
+7 -23
View File
@@ -1,15 +1,5 @@
return {
-- 2023 CX1
--[[{
Type = "FindImpactsTask",
AsteroidName = "2023 CX1",
KernelDirectory = "C:/Users/malej60/Documents/develop/b612-oribt-uncertainty/data/2023 CX1/openspace_variants/",
OutputFilename = "C:/Users/malej60/Documents/develop/b612-oribt-uncertainty/data/2023 CX1/impacts/2023-CX1-impacts.txt",
TimeIntervalStart = "2023-02-13T02:39:00.000",
TimeIntervalEnd = "2023-02-13T03:39:00.000"
}]]--
{
Type = "ImpactCorridorTask",
AsteroidName = "2023 CX1",
@@ -17,23 +7,17 @@ return {
OutputFilename = "C:/Users/malej60/Documents/develop/b612-oribt-uncertainty/images/new/2023-CX1-test.png",
ImageWidth = 5400,
ImageHeight = 2700,
--ColorMap = "${SYNC}/http/default_colormaps_sequential/1/GnBu.cmap"
}
-- 2004 MN4
--[[{
Type = "FindImpactsTask",
AsteroidName = "2004 MN4",
KernelDirectory = "C:/Users/malej60/Documents/develop/b612-oribt-uncertainty/data/2004 MN4/openspace_variants_high_ip/",
OutputFilename = "C:/Users/malej60/Documents/develop/b612-oribt-uncertainty/data/2004 MN4/impacts/2004-MN4-impacts.txt",
TimeIntervalStart = "2029-02-01T00:00:00.000",
TimeIntervalEnd = "2029-08-17T00:00:00.000",
}
FilterStrength = 2.5,
ColorMap = "${SYNC}/http/default_colormaps_sequential/1/GnBu.cmap"
},
{
Type = "ImpactCorridorTask",
AsteroidName = "2004 MN4",
ImpactFile = "C:/Users/malej60/Documents/develop/b612-oribt-uncertainty/data/2004 MN4/impacts/2004-MN4-impacts.txt",
OutputFilename = "C:/Users/malej60/Documents/develop/b612-oribt-uncertainty/images/new/2004-MN4-test.png",
ImageWidth = 5400,
ImageHeight = 2700,
FilterStrength = 2.5,
ColorMap = "${SYNC}/http/default_colormaps_sequential/1/GnBu.cmap"
}]]--
}
}
+93 -57
View File
@@ -24,6 +24,7 @@
#include <modules/neoviz/tasks/impactcorridortask.h>
#include <openspace/data/dataloader.h>
#include <openspace/documentation/verifier.h>
#include <openspace/util/spicemanager.h>
#include <ghoul/filesystem/filesystem.h>
@@ -62,14 +63,22 @@ namespace {
for (int x = -halfSize; x <= halfSize; x++) {
for (int y = -halfSize; y <= halfSize; y++) {
distanceSqared = x * x + y * y;
kernel[x + halfSize][y + halfSize] =
(exp(-(distanceSqared) / sig)) / (std::numbers::pi * sig);
sum += kernel[x + halfSize][y + halfSize];
double value = (exp(-(distanceSqared) / sig)) / (std::numbers::pi * sig);
kernel[x + halfSize][y + halfSize] = value;
if (x == 0 && y == 0) {
sum = value;
}
//sum += kernel[x + halfSize][y + halfSize];
}
}
// No need to normalize the kernel, as that will be done before the image is
// saved to file anyway
// Normalize so the middle of the kernel is 1
for (int i = 0; i < kernelSize; ++i) {
for (int j = 0; j < kernelSize; ++j) {
kernel[i][j] /= sum;
}
}
return kernel;
}
@@ -161,42 +170,66 @@ void ImpactCorridorTask::perform(const Task::ProgressCallback& progressCallback)
progressCallback(0.0);
// Apply a transfer function to bring color to the image
std::vector<GLubyte> image = std::vector<GLubyte>(Size, 0);
if (_hasColorMapFile) {
LINFO("Processing image...");
int colorMapWidth = 0;
int colorMapHeight = 0;
int colorMapChannels = 0;
unsigned char* colorMapData = stbi_load(
_colorMap.string().c_str(),
&colorMapWidth,
&colorMapHeight,
&colorMapChannels,
0
);
ghoul_assert(colorMapHeight != 1, "Color map needs to be 1D");
progressCallback(0.0);
// Read the color map image
openspace::dataloader::ColorMap colorMap = dataloader::color::loadFile(_colorMap);
// Apply the transfer function for each pixel and store the resulting color in
// the image
for (int p = 0, i = 0; p < NumPixels && i < Size; ++p, i += NChannels) {
glm::vec4 color = glm::ivec4(0);
double normalizedValue = rawData[p] / maxRawSaturation;
if (std::abs(normalizedValue) < std::numeric_limits<double>::epsilon()) {
progressCallback((p + 1) / static_cast<float>(NumPixels));
continue;
}
// Find the color in the color map that cooresponds to the value
int colorMapIndex = static_cast<int>(std::round(
normalizedValue * (colorMap.entries.size() - 1)
));
color = colorMap.entries[colorMapIndex];
color.a = normalizedValue;
// Save the color after the transferfunction to the image
for (int c = 0; c < NChannels; ++c) {
// Watch out for overflow
double clampedSaturation = std::min(std::round(color[c] * 255.0), 255.0);
image[i + c] = static_cast<GLubyte>(clampedSaturation);
}
progressCallback((p + 1)/ static_cast<float>(NumPixels));
}
}
else {
// Convert the raw image to a normalized image that can be written to file
LINFO("Finalizing image...");
for (int p = 0, i = 0; p < NumPixels && i < Size; ++p, i += NChannels) {
// Note that this normalizes all channels equally
double saturation = rawData[p] / maxRawSaturation;
if (std::abs(saturation) < std::numeric_limits<double>::epsilon()) {
progressCallback((p + 1) / static_cast<float>(NumPixels));
continue;
}
// Watch out for overflow
double clampedSaturation = std::min(std::round(saturation * 255.0), 255.0);
// Save the normalized saturation to all color channels
for (int c = 0; c < NChannels; ++c) {
image[i + c] = static_cast<GLubyte>(clampedSaturation);
}
progressCallback((p + 1) / static_cast<float>(NumPixels));
}
}
// @TODO: Use the night layer image to estimate impact risk due to population density
// Convert the raw image to a normalized image that can be written to file
LINFO("Finalizing image...");
std::vector<GLubyte> image = std::vector<GLubyte>(Size, 0);
for (unsigned int i = 0, p = 0; p < NumPixels && i < Size; ++p, i += NChannels) {
// Note that this normalizes all channels equally
double saturation = rawData[p] / maxRawSaturation;
// Watch out for overflow
double clampedSaturation = std::min(std::round(saturation * 255.0), 255.0);
// Save the normalized saturation to all color channels
for (int c = 0; c < NChannels; ++c) {
image[i + c] = static_cast<GLubyte>(clampedSaturation);
}
progressCallback(p / static_cast<float>(NumPixels));
}
// Create and save the resulting impact map image
stbi_write_png(
_outputFilename.string().c_str(),
@@ -258,33 +291,35 @@ ImpactCorridorTask::ImpactCoordinate ImpactCorridorTask::readImpactCoordinate(
}
int ImpactCorridorTask::pixelIndex(int pixelW, int pixelH, int numChannels, int imageWidth,
int imageHeight)
int imageHeight, bool allowWrap)
{
if (pixelH < 0) {
// Mirror in the noth pole line
pixelH = std::abs(pixelH);
if (allowWrap) {
if (pixelH < 0) {
// Mirror in the noth pole line
pixelH = std::abs(pixelH);
// Mirror in the 0 longitude line (Greenwich)
int greenwich = static_cast<int>(std::round(imageWidth / 2.0));
pixelW = greenwich + (greenwich - pixelW);
// Mirror in the 0 longitude line (Greenwich)
int greenwich = static_cast<int>(std::round(imageWidth / 2.0));
pixelW = greenwich + (greenwich - pixelW);
}
else if (pixelH >= imageHeight) {
// Mirror in the south pole line
pixelH = imageHeight - (pixelH - imageHeight);
}
else if (pixelH >= imageHeight) {
// Mirror in the south pole line
pixelH = imageHeight - (pixelH - imageHeight);
// Mirror in the 0 longitude line (Greenwich)
int greenwich = static_cast<int>(std::round(imageWidth / 2.0));
pixelW = greenwich + (greenwich - pixelW);
}
// Mirror in the 0 longitude line (Greenwich)
int greenwich = static_cast<int>(std::round(imageWidth / 2.0));
pixelW = greenwich + (greenwich - pixelW);
}
if (pixelW < 0) {
// Go over to the left side of the international daytime line
pixelW = imageWidth - std::abs(pixelW);
}
else if (pixelW >= imageWidth) {
// Go over to the right side of the international daytime line
pixelW = pixelW - imageWidth;
if (pixelW < 0) {
// Go over to the left side of the international daytime line
pixelW = imageWidth - std::abs(pixelW);
}
else if (pixelW >= imageWidth) {
// Go over to the right side of the international daytime line
pixelW = pixelW - imageWidth;
}
}
int pixelIndex = numChannels * pixelW + numChannels * imageWidth * pixelH;
@@ -325,7 +360,8 @@ std::vector<double> ImpactCorridorTask::rawImpactData(
pixelH + h,
1,
_imageWidth,
_imageHeight
_imageHeight,
true
);
// Get the eucledian distance from the center of the brush and calculate
+3 -1
View File
@@ -77,10 +77,12 @@ private:
* \param numChannels The number of color channels in the image
* \param imageWidth The total width of the image in pixels
* \param imageHeight The total height of the image in pixels
* \param allowWrap Whether to allow the pixel to wrap around the image if it is
* located outside
* \return The index of the first channel of the given pixel in the flat data list
*/
int pixelIndex(int pixelW, int pixelH, int numChannels, int imageWidth,
int imageHeight);
int imageHeight, bool allowWrap);
/**
* Create a raw impact map image data list that only have information with impact