Skip to content

Commit

Permalink
Merge pull request pioneerspacesim#5758 from Gliese852/holes-in-space
Browse files Browse the repository at this point in the history 
Fix star background generation
  • Loading branch information
Webster Sheets authored Feb 15, 2024
2 parents a3339aa + eada00c commit 60f643c
Show file tree
Hide file tree
Showing 3 changed files with 124 additions and 54 deletions.
160 changes: 108 additions & 52 deletions src/Background.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -30,6 +30,7 @@
#include <SDL_stdinc.h>
#include <iostream>
#include <sstream>
#include <numeric>

using namespace Graphics;

Expand Down Expand Up @@ -250,7 +251,6 @@ namespace Background {

struct StarQueryInfo {
const SystemPath *systemPath;
uint32_t numStars;
int32_t sectorMin;
int32_t sectorMax;
int32_t visibleRadiusSqr;
Expand All @@ -259,20 +259,22 @@ namespace Background {
float brightnessFactor;
};

class SampleStarTask : public Task {
class SampleStarsTask : public Task {
public:
SampleStarTask(RefCountedPtr<Galaxy> galaxy, StarQueryInfo info, StarInfo *outStars, TaskRange range) :
SampleStarsTask(RefCountedPtr<Galaxy> galaxy, StarQueryInfo info, int32_t starsLimit, StarInfo &stars, double &medianBrightness, TaskRange range) :
Task(range),
galaxy(galaxy),
info(info),
outStars(outStars)
starsLimit(starsLimit),
stars(stars),
medianBrightness(medianBrightness)
{
stars.pos.reserve(info.numStars);
stars.color.reserve(info.numStars);
stars.brightness.reserve(info.numStars);
stars.pos.reserve(starsLimit);
stars.color.reserve(starsLimit);
stars.brightness.reserve(starsLimit);
}

void SampleStars(TaskRange range)
void OnExecute(TaskRange range) override
{
PROFILE_SCOPED()
const SystemPath *systemPath = info.systemPath;
Expand All @@ -286,7 +288,7 @@ namespace Background {
for (Sint32 z = minZ; z <= maxZ; ++z) {
SystemPath sys(systemPath->sectorX + x, systemPath->sectorY + y, systemPath->sectorZ + z);

if (SystemPath::SectorDistanceSqr(sys, *systemPath) * Sector::SIZE >= info.visibleRadiusSqr)
if (SystemPath::SectorDistanceSqr(sys, *systemPath) * Sector::SIZE * Sector::SIZE >= info.visibleRadiusSqr)
continue; // early out

// TODO: we're generating these sectors manually and not caching for two reasons:
Expand All @@ -299,7 +301,7 @@ namespace Background {
RefCountedPtr<const Sector> sec = galaxy->GetGenerator()->Generate<Sector, SectorCache>(galaxy, sys, nullptr);

// add as many systems as we can
const size_t numSystems = std::min(info.numStars - stars.pos.size(), sec->m_systems.size());
const size_t numSystems = std::min(starsLimit - stars.pos.size(), sec->m_systems.size());
for (size_t systemIndex = 0; systemIndex < numSystems; systemIndex++) {
const Sector::System *ss = &(sec->m_systems[systemIndex]);

Expand Down Expand Up @@ -332,21 +334,14 @@ namespace Background {
}

// Don't process any more sectors if we've generated our quota of stars.
if (stars.pos.size() >= info.numStars)
if (stars.pos.size() >= starsLimit)
break;
}
}
}
}

// Do the brightness sort on worker threads using the worker's subset
// of stars rather than on the main thread with all stars.
void SortStars()
{
PROFILE_SCOPED()
const size_t numStars = stars.pos.size();

// find the median brightness of all visible stars
const size_t numStars = stars.pos.size();
std::vector<uint32_t> sortedBrightnessIndex;
sortedBrightnessIndex.reserve(numStars);

Expand All @@ -358,11 +353,33 @@ namespace Background {
return stars.brightness[a] > stars.brightness[b];
});

double medianBrightness = 0.0;
constexpr float medianPosition = 0.7;
if (numStars > 0) {
medianBrightness = stars.brightness[sortedBrightnessIndex[Clamp<uint32_t>(medianPosition * numStars, 0, numStars - 1)]];
}
}

RefCountedPtr<Galaxy> galaxy;
const StarQueryInfo info;
const int32_t starsLimit;
StarInfo &stars;
double &medianBrightness;
};

class SortStarsTask : public Task {
public:
SortStarsTask(StarQueryInfo info, StarInfo &stars, double medianBrightness) :
info(info),
stars(stars),
medianBrightness(medianBrightness)
{}

// Do the brightness sort on worker threads using the worker's subset
// of stars rather than on the main thread with all stars.
virtual void OnExecute(TaskRange) override
{
PROFILE_SCOPED()
const size_t numStars = stars.pos.size();

for (size_t i = 0; i < numStars; ++i) {
// dividing through the median helps bringing the logarithmic brightnesses to a scala that is easier to work with
Expand Down Expand Up @@ -390,27 +407,31 @@ namespace Background {
}
}

virtual void OnExecute(TaskRange range) override
{
SampleStars(range);
SortStars();
}

virtual void OnComplete() override
{
PROFILE_SCOPED();
const StarQueryInfo info;
StarInfo &stars;
const double medianBrightness;
};

outStars->pos.insert(outStars->pos.end(), stars.pos.begin(), stars.pos.end());
outStars->color.insert(outStars->color.end(), stars.color.begin(), stars.color.end());
outStars->brightness.insert(outStars->brightness.end(), stars.brightness.begin(), stars.brightness.end());
}
// https://en.wikipedia.org/wiki/Spherical_segment
static double spherical_segment_volume(double h, double r1_sq, double r2_sq)
{
return M_PI * h / 6.0 * (3 * r1_sq + 3 * r2_sq + h * h);
}

RefCountedPtr<Galaxy> galaxy;
const StarQueryInfo info;
static double spherical_circle_radius_sq(double h, double R)
{
// distance to sphere center
double l = std::abs(R - h);
return R * R - l * l;
}

StarInfo stars;
StarInfo *outStars;
};
static double task_spherical_segment_volume(uint32_t begin, uint32_t end, double r)
{
double h = end - begin;
double r1_sq = spherical_circle_radius_sq(begin, r);
double r2_sq = spherical_circle_radius_sq(end, r);
return spherical_segment_volume(h, r1_sq, r2_sq);
}

void Starfield::Fill(Random &rand, const SystemPath *const systemPath, RefCountedPtr<Galaxy> galaxy)
{
Expand Down Expand Up @@ -450,40 +471,69 @@ namespace Background {
// but don't split the number of stars too much that we have visible brightness "patches"
const uint32_t numTasks = std::min(graph->GetNumWorkerThreads() + 1, 8U);

/* the number of visible systems is in a cubic relationship with the visible radius,
i.e. visibleRadius = x * numberSystems^(1/3)
I experimentally determined that x is approximately 3.89
and that stays probably the same as long as the galaxy has the same system density */
const Sint32 visibleRadius = std::min<Sint32>(BG_STAR_RADIUS_MAX, 3.89 * pow((float)NUM_BG_STARS, 1.0 / 3.0));
// judging by the current sector generator, maximum average number
// of stars in a sector is 6
// It’s easy to express what the radius of a ball should be so that
// at such a density it would contain approximately NUM_BG_STARS stars:
const double density = 6.0;
const double maxBall = pow(3.0 / 4.0 / M_PI * (double)NUM_BG_STARS / density, 1.0 / 3.0);
const int32_t visibleRadius = std::min<int32_t>(BG_STAR_RADIUS_MAX, maxBall * Sector::SIZE);

StarQueryInfo info;
info.systemPath = systemPath;
info.numStars = NUM_BG_STARS / numTasks;
info.sectorMin = -(visibleRadius / Sector::SIZE); // lyrs_radius / sector_size_in_lyrs
info.sectorMax = visibleRadius / Sector::SIZE; // lyrs_radius / sector_size_in_lyrs
info.visibleRadiusSqr = (visibleRadius * visibleRadius);
info.colorMin = Color((Uint8)(m_rMin * 255), (Uint8)(m_gMin * 255), (Uint8)(m_rMin * 255));
info.colorMax = Color((Uint8)(m_rMax * 255), (Uint8)(m_gMax * 255), (Uint8)(m_rMax * 255));
info.brightnessFactor = brightnessApparentSizeFactor;

TaskSet *pickStarTaskSet = new TaskSet();
TaskSet *sampleStarsTaskSet = new TaskSet();

int32_t starsLeft = NUM_BG_STARS;
const double realRadius = info.sectorMax + 0.5;
const double realDensity = NUM_BG_STARS / (M_PI / 0.75 * realRadius * realRadius * realRadius);

std::vector<StarInfo> taskStars(numTasks);
std::vector<double> taskMedians(numTasks);

// Split the visible area of the galaxy up into separate tasks
uint32_t current = 0;
int32_t range_step = (info.sectorMax - info.sectorMin) / numTasks;
for (size_t i = 0; i < numTasks; i++) {
int32_t starsLimit;
uint32_t end = current + range_step;
if (i + 1 == numTasks)
if (i + 1 == numTasks) {
end = (info.sectorMax - info.sectorMin);
starsLimit = starsLeft;
} else {
starsLimit = realDensity * task_spherical_segment_volume(current, end + 1, realRadius);
starsLeft -= starsLimit;
}

pickStarTaskSet->AddTask(new SampleStarTask(galaxy, info, &stars, { current, end }));
current = end;
sampleStarsTaskSet->AddTask(new SampleStarsTask(galaxy, info, starsLimit, taskStars[i], taskMedians[i], { current, end }));
current = end + 1;
}

// We can't make progress until all stars are gathered, so run the
// star collection on the 'main' thread as well.
auto handle = graph->QueueTaskSet(pickStarTaskSet);
graph->WaitForTaskSet(handle);
auto sampleHandle = graph->QueueTaskSet(sampleStarsTaskSet);
graph->WaitForTaskSet(sampleHandle);

double medianBrightness = std::reduce(taskMedians.begin(), taskMedians.end()) / taskMedians.size();

TaskSet *sortStarsTaskSet = new TaskSet();
for (size_t i = 0; i < numTasks; i++) {
sortStarsTaskSet->AddTask(new SortStarsTask(info, taskStars[i], medianBrightness));
}
auto sortHandle = graph->QueueTaskSet(sortStarsTaskSet);
graph->WaitForTaskSet(sortHandle);

for(auto &item : taskStars) {
stars.pos.insert(stars.pos.end(), item.pos.begin(), item.pos.end());
stars.color.insert(stars.color.end(), item.color.begin(), item.color.end());
stars.brightness.insert(stars.brightness.end(), item.brightness.begin(), item.brightness.end());
}
}
num = stars.pos.size();
Output("Stars picked from galaxy: %d\n", stars.pos.size());
Expand All @@ -510,7 +560,7 @@ namespace Background {
const float u = float(rand.Double(-1.0, 1.0));

// squeeze the starfield a bit to get more density near horizon using matrix3x3f::Scale
const auto star = matrix3x3f::Scale(1.0, 0.4, 1.0) * (vector3f(sqrt(1.0f - u * u) * cos(theta), u, sqrt(1.0f - u * u) * sin(theta)).Normalized() * 1000.0f);
const auto star = matrix3x3f::Scale(1.0, 1.0, 0.4) * (vector3f(sqrt(1.0f - u * u) * cos(theta), u, sqrt(1.0f - u * u) * sin(theta)).Normalized() * 1000.0f);

stars.pos[i] = star;
stars.color[i] = col;
Expand Down Expand Up @@ -679,7 +729,13 @@ namespace Background {
m_milkyWay.Draw();
}
if (DRAW_STARS & m_drawFlags) {
m_renderer->SetTransform(matrix4x4f(transform));
auto Zup_to_Yup = matrix4x4f::FromRowMajor({
0, 1, 0, 0,
0, 0, 1, 0,
1, 0, 0, 0,
0, 0, 0, 1
});
m_renderer->SetTransform(matrix4x4f(transform) * Zup_to_Yup);
m_starField.Draw();
}
}
Expand Down
4 changes: 2 additions & 2 deletions src/galaxy/SystemPath.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -85,15 +85,15 @@ class SystemPath : public LuaWrappable {
{
const Sint32 x = b.sectorX - a.sectorX;
const Sint32 y = b.sectorY - a.sectorY;
const Sint32 z = b.sectorZ - b.sectorZ;
const Sint32 z = b.sectorZ - a.sectorZ;
return sqrt(x * x + y * y + z * z); // sqrt is slow
}

static inline double SectorDistanceSqr(const SystemPath &a, const SystemPath &b)
{
const Sint32 x = b.sectorX - a.sectorX;
const Sint32 y = b.sectorY - a.sectorY;
const Sint32 z = b.sectorZ - b.sectorZ;
const Sint32 z = b.sectorZ - a.sectorZ;
return (x * x + y * y + z * z); // return the square of the distance
}

Expand Down
14 changes: 14 additions & 0 deletions src/matrix4x4.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -105,6 +105,20 @@ class matrix4x4 {
r[7] = cell[6];
r[8] = cell[10];
}
// matrix4x4 is column-major
// but it seems more natural to write the matrix by row
template <std::size_t N>
static matrix4x4 FromRowMajor(const T (&a)[N])
{
static_assert(N == 16, "Incorrect number of elements specified");
matrix4x4 m;
auto &c = m.cell;
c[0] = a[0]; c[4] = a[1]; c[8] = a[2]; c[12] = a[3];
c[1] = a[4]; c[5] = a[5]; c[9] = a[6]; c[13] = a[7];
c[2] = a[8]; c[6] = a[9]; c[10] = a[10]; c[14] = a[11];
c[3] = a[12]; c[7] = a[13]; c[11] = a[14]; c[15] = a[15];
return m;
}
static matrix4x4 Identity()
{
matrix4x4 m = matrix4x4(0.0);
Expand Down

0 comments on commit 60f643c

Please sign in to comment.