mirror of
https://bitbucket.org/cosmicvoids/vide_public.git
synced 2025-07-04 15:21:11 +00:00
Re-Implemented jozov2 using priority_queues (far less complicated and maybe better performing)
This commit is contained in:
Guilhem Lavaux
parent
9a63f289d7
commit
91a222a3c3
2 changed files with 124 additions and 124 deletions
|
|
@ -139,10 +139,10 @@ int main(int argc,char **argv)
|
|||
if (z[i].np == 1)
|
||||
continue;
|
||||
|
||||
txt << format("%d %d %d %e %e %d %d %e %d %f %6.2e")
|
||||
txt << format("%d %d %d %e %e %d %d %e %d %f %6.2e %f")
|
||||
% (h+1) % i % z[i].core % p[z[i].core].dens % z[i].vol
|
||||
% z[i].np % z[i].nhl % z[i].voljoin % z[i].npjoin
|
||||
% z[i].denscontrast % prob << endl;
|
||||
% z[i].denscontrast % prob % z[i].leak << endl;
|
||||
|
||||
} /* h+1 to start from 1, not zero */
|
||||
txt.close();
|
||||
|
|
|
|||
|
|
@ -2,6 +2,7 @@
|
|||
#include <omp.h>
|
||||
#endif
|
||||
|
||||
#include <queue>
|
||||
#include <set>
|
||||
#include <iostream>
|
||||
#include <fstream>
|
||||
|
|
@ -13,6 +14,40 @@
|
|||
using namespace std;
|
||||
using boost::format;
|
||||
|
||||
struct ZoneDensityPair
|
||||
{
|
||||
int h;
|
||||
double density;
|
||||
|
||||
bool operator<(const ZoneDensityPair& p2) const
|
||||
{
|
||||
return density > p2.density;
|
||||
}
|
||||
};
|
||||
|
||||
typedef priority_queue<ZoneDensityPair> ZoneQueue;
|
||||
|
||||
static void build_process_queue(ZoneQueue& q, ZONE *z, char *inyet, int h)
|
||||
{
|
||||
ZoneDensityPair zdp;
|
||||
ZONE& z_h = z[h];
|
||||
bool interior = true;
|
||||
|
||||
assert(inyet[h] == 1);
|
||||
for (int za = 0; za < z_h.nadj; za++)
|
||||
{
|
||||
zdp.h = z_h.adj[za];
|
||||
zdp.density = z_h.slv[za];
|
||||
if (inyet[zdp.h] == 0)
|
||||
{
|
||||
q.push(zdp);
|
||||
interior = false;
|
||||
}
|
||||
}
|
||||
if (interior)
|
||||
inyet[h] = 2;
|
||||
}
|
||||
|
||||
void doWatershed(PARTICLE *p, pid_t np, ZONE *z, int numZones, float maxvol, float voltol)
|
||||
{
|
||||
int *iord;
|
||||
|
|
@ -28,36 +63,33 @@ void doWatershed(PARTICLE *p, pid_t np, ZONE *z, int numZones, float maxvol, flo
|
|||
|
||||
iord = new int[numZones];
|
||||
|
||||
findrtop(sorter, numZones, iord, numZones);
|
||||
// findrtop(sorter, numZones, iord, numZones);
|
||||
delete[] sorter;
|
||||
|
||||
#pragma omp parallel
|
||||
//#pragma omp parallel
|
||||
{
|
||||
char *inyet, *inyet2;
|
||||
int *zonelist, *zonelist2;
|
||||
int nhl;
|
||||
int *links = new int[NLINKS];
|
||||
bool *done_zones;
|
||||
int prev_ii = -1;
|
||||
|
||||
inyet = new char[numZones];
|
||||
inyet2 = new char[numZones];
|
||||
zonelist = new int[numZones];
|
||||
zonelist2 = new int[numZones];
|
||||
done_zones = new bool[numZones];
|
||||
|
||||
fill(inyet, inyet + numZones, 0);
|
||||
fill(inyet2, inyet2 + numZones, 0);
|
||||
fill(done_zones, done_zones + numZones, false);
|
||||
|
||||
nhl = 0;
|
||||
#pragma omp for schedule(dynamic,1)
|
||||
//#pragma omp for schedule(dynamic,1)
|
||||
for (int h = 0; h < numZones; h++)
|
||||
{
|
||||
int nhlcount = 0;
|
||||
float lowvol;
|
||||
float lowvol, z_cur_core_dens;
|
||||
bool beaten;
|
||||
set<int> to_process;
|
||||
priority_queue<ZoneDensityPair> to_process;
|
||||
int link0;
|
||||
|
||||
for (int hl = 0; hl < nhl; hl++)
|
||||
inyet[zonelist[hl]] = 0;
|
||||
|
|
@ -65,144 +97,113 @@ void doWatershed(PARTICLE *p, pid_t np, ZONE *z, int numZones, float maxvol, flo
|
|||
zonelist[0] = h;
|
||||
inyet[h] = 1;
|
||||
nhl = 1;
|
||||
to_process.insert(h);
|
||||
z[h].npjoin = z[h].np;
|
||||
z_cur_core_dens = p[z[h].core].dens;
|
||||
|
||||
build_process_queue(to_process, z, inyet, h);
|
||||
|
||||
do {
|
||||
/* Find the lowest-volume (highest-density) adjacency */
|
||||
int nl = 0;
|
||||
|
||||
beaten = false;
|
||||
lowvol = BIGFLT;
|
||||
|
||||
set<int>::iterator iter = to_process.begin();
|
||||
while (iter != to_process.end())
|
||||
{
|
||||
int h2 = *iter;
|
||||
ZONE& z_cur = z[h2];
|
||||
bool interior = true, touched = false;
|
||||
assert(inyet[h2] == 1);
|
||||
|
||||
for (int za = 0; za < z_cur.nadj; za++)
|
||||
{
|
||||
if (inyet[z_cur.adj[za]] == 0)
|
||||
{
|
||||
interior = false;
|
||||
if (z_cur.slv[za] == lowvol)
|
||||
{
|
||||
links[nl] = z_cur.adj[za];
|
||||
touched = true;
|
||||
nl ++;
|
||||
if (nl == NLINKS)
|
||||
{
|
||||
printf("Too many links with the same rho_sl! Increase NLINKS from %d\n",nl);
|
||||
exit(0);
|
||||
}
|
||||
}
|
||||
else
|
||||
if (z_cur.slv[za] < lowvol)
|
||||
{
|
||||
lowvol = z_cur.slv[za];
|
||||
links[0] = z_cur.adj[za];
|
||||
nl = 1;
|
||||
touched = true;
|
||||
}
|
||||
}
|
||||
}
|
||||
++iter;
|
||||
if (interior || !touched)
|
||||
to_process.erase(h2);
|
||||
if (interior)
|
||||
inyet[h2] = 2; /* No bordering exter. zones */
|
||||
}
|
||||
|
||||
if (nl == 0)
|
||||
if (to_process.empty())
|
||||
{
|
||||
beaten = true;
|
||||
z[h].leak = maxvol;
|
||||
continue;
|
||||
}
|
||||
|
||||
do
|
||||
{
|
||||
lowvol = to_process.top().density;
|
||||
link0 = to_process.top().h;
|
||||
to_process.pop();
|
||||
}
|
||||
while ((inyet[link0] != 0) && (!to_process.empty()));
|
||||
|
||||
if (to_process.empty())
|
||||
{
|
||||
beaten = true;
|
||||
z[h].leak = maxvol;
|
||||
continue;
|
||||
}
|
||||
|
||||
nl++;
|
||||
|
||||
if (lowvol > voltol)
|
||||
if (lowvol > voltol)
|
||||
{
|
||||
beaten = true;
|
||||
z[h].leak = lowvol;
|
||||
continue;
|
||||
}
|
||||
|
||||
for (int l = 0; l < nl; l++)
|
||||
{
|
||||
//if (!done_zones[links[l]]) /* Equivalent to p[z[links[l]].core].dens < p[z[h].core].dens) as zones are sorted. */
|
||||
if (p[z[links[l]].core].dens < p[z[h].core].dens)
|
||||
beaten = true;
|
||||
}
|
||||
|
||||
if (beaten)
|
||||
if (p[z[link0].core].dens < z_cur_core_dens)
|
||||
{
|
||||
beaten = true;
|
||||
z[h].leak = lowvol;
|
||||
continue;
|
||||
}
|
||||
|
||||
/* Add everything linked to the link(s) */
|
||||
int nhl2 = 0;
|
||||
for (int l = 0; l < nl; l++)
|
||||
{
|
||||
if (inyet2[links[l]] == 0)
|
||||
{
|
||||
zonelist2[nhl2] = links[l];
|
||||
inyet2[links[l]] = 1;
|
||||
nhl2 ++;
|
||||
bool added = true;
|
||||
while (added && !beaten)
|
||||
{
|
||||
added = false;
|
||||
for (int hl = 0; (hl < nhl2) && (!beaten); hl++)
|
||||
{
|
||||
int h2 = zonelist2[hl];
|
||||
|
||||
if (inyet2[h2] == 1) {
|
||||
bool interior = true; /* Guilty until proven innocent */
|
||||
for (int za = 0; za < z[h2].nadj; za ++) {
|
||||
int link2 = z[h2].adj[za];
|
||||
if ((inyet[link2]+inyet2[link2]) == 0) {
|
||||
interior = false;
|
||||
if (z[h2].slv[za] <= lowvol) {
|
||||
//if (!done_zones[link2]) { // Equivalent to p[z[link2].core].dens < p[z[h].core].dens)
|
||||
if (p[z[link2].core].dens < p[z[h].core].dens) {
|
||||
beaten = true;
|
||||
break;
|
||||
}
|
||||
zonelist2[nhl2] = link2;
|
||||
inyet2[link2] = 1;
|
||||
nhl2++;
|
||||
added = true;
|
||||
}
|
||||
}
|
||||
}
|
||||
if (interior)
|
||||
inyet2[h2] = 2;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
/* See if there's a beater */
|
||||
if (beaten) {
|
||||
z[h].leak = lowvol;
|
||||
} else {
|
||||
for (int h2 = 0; h2 < nhl2; h2++) {
|
||||
int new_h = zonelist2[h2];
|
||||
|
||||
zonelist2[0] = link0;
|
||||
inyet2[link0] = 1;
|
||||
nhl2=1;
|
||||
|
||||
bool added = true;
|
||||
while (added && !beaten)
|
||||
{
|
||||
added = false;
|
||||
for (int hl = 0; (hl < nhl2) && (!beaten); hl++)
|
||||
{
|
||||
int h2 = zonelist2[hl];
|
||||
|
||||
if (inyet2[h2] == 1) {
|
||||
bool interior = true; /* Guilty until proven innocent */
|
||||
|
||||
for (int za = 0; za < z[h2].nadj; za ++) {
|
||||
int link2 = z[h2].adj[za];
|
||||
|
||||
if ((inyet[link2]+inyet2[link2]) == 0) {
|
||||
interior = false;
|
||||
if (z[h2].slv[za] <= lowvol) {
|
||||
if (p[z[link2].core].dens < z_cur_core_dens) {
|
||||
beaten = true;
|
||||
break;
|
||||
}
|
||||
zonelist2[nhl2] = link2;
|
||||
inyet2[link2] = 1;
|
||||
nhl2++;
|
||||
added = true;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
if (interior)
|
||||
inyet2[h2] = 2;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
/* See if there's a beater */
|
||||
if (beaten) {
|
||||
z[h].leak = lowvol;
|
||||
} else {
|
||||
for (int h2 = 0; h2 < nhl2; h2++) {
|
||||
int new_h = zonelist2[h2];
|
||||
|
||||
zonelist[nhl] = new_h;
|
||||
assert(inyet[new_h] == 0);
|
||||
inyet[new_h] = 1;
|
||||
if (inyet2[new_h] != 2)
|
||||
to_process.insert(new_h);
|
||||
build_process_queue(to_process, z, inyet, new_h);
|
||||
nhl++;
|
||||
z[h].npjoin += z[new_h].np;
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
for (int hl = 0; hl < nhl2; hl++)
|
||||
inyet2[zonelist2[hl]] = 0;
|
||||
if (nhl/10000 > nhlcount) {
|
||||
|
|
@ -212,10 +213,12 @@ void doWatershed(PARTICLE *p, pid_t np, ZONE *z, int numZones, float maxvol, flo
|
|||
(cout << format(" %d [%d]") % nhl % to_process.size()).flush();
|
||||
nhlcount = nhl/10000;
|
||||
}
|
||||
} while((lowvol < BIGFLT) && (!beaten));
|
||||
|
||||
z[h].denscontrast = z[h].leak/p[z[h].core].dens;
|
||||
if (z[h].denscontrast < 1.) z[h].denscontrast = 1.;
|
||||
}
|
||||
while((lowvol < BIGFLT) && (!beaten));
|
||||
|
||||
z[h].denscontrast = z[h].leak/p[z[h].core].dens;
|
||||
if (z[h].denscontrast < 1.)
|
||||
z[h].denscontrast = 1.;
|
||||
|
||||
|
||||
/* Don't sort; want the core zone to be first */
|
||||
|
|
@ -233,19 +236,16 @@ void doWatershed(PARTICLE *p, pid_t np, ZONE *z, int numZones, float maxvol, flo
|
|||
z[h].zonelist[q] = zonelist[q];
|
||||
}
|
||||
|
||||
done_zones[h] = true;
|
||||
z[h].nhl = nhl;
|
||||
}
|
||||
delete[] zonelist;
|
||||
delete[] zonelist2;
|
||||
delete[] links;
|
||||
delete[] inyet;
|
||||
delete[] inyet2;
|
||||
delete[] done_zones;
|
||||
|
||||
|
||||
}
|
||||
delete[] iord;
|
||||
|
||||
|
||||
double maxdenscontrast = 0;
|
||||
#pragma omp parallel shared(maxdenscontrast)
|
||||
{
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue