The Perl and Raku Conference 2025: Greenville, South Carolina - June 27-29 Learn more

/* -*- Mode: C -*- */
#define PERL_NO_GET_CONTEXT 1
#include "EXTERN.h"
#include "perl.h"
#include "XSUB.h"
#include "ppport.h"
#define X0 0
#define Y0 1
#define X1 2
#define Y1 3
#define IDOFFSET (sizeof(NV) * 4)
#define BRUTEFORCECUTOFF 40
#if 1
#define DP(f)
#define DUMP(msg, av, start)
#define my_assert(a) assert(a)
#else
#define my_assert(a) if(!(a)) _failed_assertion(aTHX_ #a, __LINE__, __FILE__)
#define DP(f) f
#define DUMP(msg, av, start) _dump(aTHX_ (msg), (av), (start))
static void
_dump(pTHX_ char *msg, AV *rects, I32 start) {
I32 end = av_len(rects) + 1;
SV **svs = AvARRAY(rects);
I32 i;
fprintf (stderr, "%s = start: %d, end: %d", msg, start, end);
for (i = start; i < end; i++) {
SV *sv = svs[i];
if (SvOK(sv)) {
STRLEN len, j;
NV* nv = (NV*)SvPV(svs[i], len);
IV* iv = (IV*)(SvPV_nolen(svs[i]) + IDOFFSET);
len = (len - IDOFFSET) / sizeof(IV);
fprintf (stderr, " [%.0f %.0f %.0f %.0f |", nv[0], nv[1], nv[2], nv[3]);
for (j = 0; j < len; j++)
fprintf(stderr, " %d", iv[j]);
fprintf(stderr, "]");
}
else
fprintf(stderr, " undef");
}
fprintf(stderr, "\n");
fflush(stderr);
}
void _failed_assertion(pTHX_ char *str, int line, char *file) {
fprintf(stderr, "assertion %s failed at %s line %d\n", str, line, file);
fflush(stderr);
exit(1);
}
#endif
static SV *
av_sure_fetch(pTHX_ AV *av, I32 i) {
SV **sv = av_fetch(av, i, 0);
my_assert(sv);
return *sv;
}
static NV
sqr (NV a) {
return a * a;
}
int
double_cmp(pTHX_ double *a, double *b) {
double fa = *a;
double fb = *b;
DP(printf("cmp(%f, %f) => %d\n", fa, fb, (fa < fb) ? -1 : ((fa > fb) ? 1 : 0)));
return (fa < fb) ? -1 : ((fa > fb) ? 1 : 0);
}
void
sort_inplace(pTHX_ double **v, int size) {
sortsv((SV**)v, size, (SVCOMPARE_t)&double_cmp);
}
static NV
find_best_cut(pTHX_ AV *rects, I32 start, I32 end, int dir, NV *bestv) {
NV **v0, **v1, **vc0, **vc1;
NV v, med, best;
int op, cl;
int i;
SV **svs;
I32 size = end - start;
my_assert(bestv);
DUMP("fbc in", rects, start);
DP(fprintf(stderr, "end: %d\n", end));
Newx(v0, size + 1, NV *);
Newx(v1, size + 1, NV *);
v0[size] = v1[size] = NULL;
vc0 = v0; vc1 = v1;
svs = AvARRAY(rects) + start;
size = end - start;
for (i = 0; i < size; i++) {
NV *nv = (NV*)SvPV_nolen(svs[i]);
if (dir == 'x') {
v0[i] = nv + X0;
v1[i] = nv + X1;
}
else {
v0[i] = nv + Y0;
v1[i] = nv + Y1;
}
}
v0[size] = v1[size] = NULL;
sort_inplace(aTHX_ v0, size);
sort_inplace(aTHX_ v1, size);
op = cl = 0;
med = 0.5 * size;
best = .24 * sqr(size);
my_assert(best >= 0);
while (*v0 && *v1) {
NV v, good;
NV l, r;
v = (**v0 <= **v1) ? **v0 : **v1;
while (*v0 && v == **v0) {
op++;
v0++;
}
while (*v1 && v == **v1) {
cl++;
v1++;
}
my_assert(op > 0 && op <= size);
my_assert(cl >= 0 && cl <= size);
l = op - med;
r = size - cl - med;
good = sqr(l) + sqr(r);
my_assert(good >= 0);
if (good < best) {
DP(fprintf(stderr, "find_best_cut l: %.2f, r: %.2f, good: %.2f\n", l, r, good));
best = good;
*bestv = v;
}
}
Safefree(vc0);
Safefree(vc1);
return best;
}
static void
_break(pTHX_ AV *rects, I32 start, AV *parts);
static void
_brute_force_merge(pTHX_ AV *rects, I32 start, AV *parts) {
I32 i, len;
for (i = start; i <= av_len(rects); i++) {
SV *svr = (AvARRAY(rects))[i];
I32 j;
for (j=0; j <= av_len(parts);) {
NV *r = (NV*) SvPV_nolen(svr);
SV *svp = (AvARRAY(parts))[j];
NV *p = (NV*) SvPV_nolen(svp);
if ((r[X1] > p[X0]) &&
(r[X0] < p[X1]) &&
(r[Y1] > p[Y0]) &&
(r[Y0] < p[Y1])) {
if ((r[X0] == p[X0]) &&
(r[Y0] == p[Y0]) &&
(r[X1] == p[X1]) &&
(r[Y1] == p[Y1])) {
SV *last;
STRLEN len;
char *pv = SvPV(svp, len);
sv_catpvn(svr, pv + IDOFFSET, len - IDOFFSET);
last = av_pop(parts);
if (last == svp)
SvREFCNT_dec(last);
else
av_store(parts, j, last);
}
else {
NV x[4], y[4];
/* sort xs and ys */
if (r[X0] < p[X0]) {
x[0] = r[X0]; x[1] = p[X0];
}
else {
x[0] = p[X0]; x[1] = r[X0];
}
if (r[X1] < p[X1]) {
x[2] = r[X1]; x[3] = p[X1];
}
else {
x[2] = p[X1]; x[3] = r[X1];
}
if (r[Y0] < p[Y0]) {
y[0] = r[Y0]; y[1] = p[Y0];
}
else {
y[0] = p[Y0]; y[1] = r[Y0];
}
if (r[Y1] < p[Y1]) {
y[2] = r[Y1]; y[3] = p[Y1];
}
else {
y[2] = p[Y1]; y[3] = r[Y1];
}
if ( ( ((x[3] - x[0]) > (y[3] - y[0])) &&
((x[0] != x[1]) || (x[2] != x[3])) ) ||
((y[0] == y[1]) && (y[2] == y[3])) ) {
NV b = ((sqr(x[1] - x[0]) + sqr(x[1] - x[3]))
< (sqr(x[2] - x[0]) + sqr(x[2] - x[3])) ? x[1] : x[2]);
if ((r[X0] < b) && (r[X1] > b)) {
SV *svcp = newSVsv(svr);
NV *cp = (NV*)SvPV_nolen(svcp);
cp[X0] = b;
av_push(rects, svcp);
/* fprintf(stderr, "[%f %f %f %f] -> [%f %f %f %f], [%f %f %f %f]\n",
r[0], r[1], r[2], r[3],
cp[0], cp[1], cp[2], cp[3],
b, r[1], r[2], r[3]); fflush(stderr); */
r[X1] = b;
}
if ((p[X0] < b) && (p[X1] > b)) {
SV *svcp = newSVsv(svp);
NV *cp = (NV*)SvPV_nolen(svcp);
cp[X0] = b;
av_push(parts, svcp);
/* fprintf(stderr, "[%f %f %f %f] -> [%f %f %f %f], [%f %f %f %f]\n",
p[0], p[1], p[2], p[3],
cp[0], cp[1], cp[2], cp[3],
p[0], p[1], b, p[3]); fflush(stderr); */
p[X1] = b;
}
}
else {
NV b = ((sqr(y[1] - y[0]) + sqr(y[1] - y[3]))
< (sqr(y[2] - y[0]) + sqr(y[2] - y[3])) ? y[1] : y[2]);
if ((r[Y0] < b) && (r[Y1] > b)) {
SV *svcp = newSVsv(svr);
NV *cp = (NV*)SvPV_nolen(svcp);
cp[Y0] = b;
av_push(rects, svcp);
/* fprintf(stderr, "[%f %f %f %f] -> [%f %f %f %f], [%f %f %f %f]\n",
r[0], r[1], r[2], r[3],
cp[0], cp[1], cp[2], cp[3],
r[0], b, r[2], r[3]); fflush(stderr); */
r[Y1] = b;
}
if ((p[Y0] < b) && (p[Y1] > b)) {
SV *svcp = newSVsv(svp);
NV *cp = (NV*)SvPV_nolen(svcp);
cp[Y0] = b;
av_push(parts, svcp);
/* fprintf(stderr, "[%f %f %f %f] -> [%f %f %f %f], [%f %f %f %f]\n",
p[0], p[1], p[2], p[3],
cp[0], cp[1], cp[2], cp[3],
p[0], p[1], p[2], b); fflush(stderr); */
p[Y1] = b;
}
}
continue;
}
}
j++;
}
}
/*
for (i = av_len(parts); i >= 0; i--)
av_push(rects, av_pop(parts));
*/
if ((len = av_len(parts)) >= 0) {
SV **svrs, **svps;
I32 start = av_len(rects) + 1;
av_extend(rects, start + len);
svrs = AvARRAY(rects) + start;
svps = AvARRAY(parts);
AvFILLp(parts) = -1;
AvFILLp(rects) = start + len;
do {
svrs[len] = svps[len];
svps[len] = &PL_sv_undef;
} while (--len >= 0);
}
}
static void
_brute_force_break(pTHX_ AV *rects, I32 start, AV *parts) {
I32 i, j, end1;
DUMP("bfb in", rects, start);
end1 = av_len(rects);
if (end1 < start + 1)
return;
if (end1 - start > 2 * BRUTEFORCECUTOFF ) {
SV **svs;
I32 end;
I32 middle = (start + end1 + 1) / 2;
_break(aTHX_ rects, middle, parts);
svs = AvARRAY(rects);
end = av_len(rects) + 1;
i = start;
j = end;
while (i < middle && j > middle) {
j--;
SV *tmp = svs[i];
svs[i] = svs[j];
svs[j] = tmp;
i++;
}
middle = start + end - middle;
_break(aTHX_ rects, middle, parts);
end = av_len(rects);
if ((end - start) > (end1 - start) * 1.2)
return _break(aTHX_ rects, start, parts);
while (end-- >= middle)
av_push(parts, av_pop(rects));
return _brute_force_merge(aTHX_ rects, start, parts);
}
while (--end1 >= start) {
SV *last, *next;
SV **svs;
DP(fprintf(stderr, "bfb: start: %d, end1: %d, end: %d\n", start, end1, av_len(rects) + 1));
svs = AvARRAY(rects);
last = svs[AvFILLp(rects)];
svs[AvFILLp(rects)--] = &PL_sv_undef;
next = svs[end1];
svs[end1] = last;
av_push(parts, next);
_brute_force_merge(aTHX_ rects, end1, parts);
}
return;
}
static void
_break(pTHX_ AV *rects, I32 start, AV *parts) {
NV bestx, bestxx, besty, bestyy, div;
int off;
I32 i, j, middle, end;
SV **svs;
DUMP("break", rects, start);
while (1) {
end = av_len(rects) + 1;
if ((end - start) <= BRUTEFORCECUTOFF)
return _brute_force_break(aTHX_ rects, start, parts);
bestx = find_best_cut(aTHX_ rects, start, end, 'x', &bestxx);
besty = ((bestx == 0) ? 1 : find_best_cut(aTHX_ rects, start, end, 'y', &bestyy));
if (bestx < besty) {
off = X0;
div = bestxx;
DP(fprintf(stderr, "cutting at x=%.0f, best=%.2f\n", bestxx, bestx));
}
else {
off = Y0;
div = bestyy;
DP(fprintf(stderr, "cutting at y=%.0f, best=%.2f\n", bestyy, besty));
}
svs = AvARRAY(rects);
i = start;
middle = end;
while (i < middle) {
SV *sv = svs[i];
NV n0 = ((NV*)SvPV_nolen(sv))[off];
if (n0 < div) {
middle--;
svs[i] = svs[middle];
svs[middle] = sv;
}
else
i++;
}
DUMP("b0", rects, start);
if (middle == start || middle == end)
return _brute_force_break(aTHX_ rects, start, parts);
_break(aTHX_ rects, middle, parts);
DUMP("b1", rects, start);
svs = AvARRAY(rects);
end = av_len(rects) + 1;
i = start;
j = end;
while (i < middle && j > middle) {
j--;
SV *tmp = svs[i];
svs[i] = svs[j];
svs[j] = tmp;
i++;
}
DUMP("b2", rects, start);
end += start - middle;
off += 2;
i = start;
middle = end;
DP(fprintf(stderr, "i: %d, middle: %d\n", i, middle));
while (i < middle) {
SV *sv = svs[i];
NV n0 = ((NV*)SvPV_nolen(sv))[off];
if (n0 > div) {
middle--;
svs[i] = svs[middle];
svs[middle] = sv;
}
else
i++;
}
DUMP("b3", rects, start);
if (middle == start)
return _brute_force_break(aTHX_ rects, start, parts);
start = middle;
}
/* _break(aTHX_ rects, middle); */
}
MODULE = Algorithm::BreakOverlappingRectangles PACKAGE = Algorithm::BreakOverlappingRectangles
PROTOTYPES: DISABLE
void
_break_rectangles(rects)
AV *rects;
CODE:
if (SvMAGICAL((SV*)rects))
Perl_croak(aTHX_ "internal error: unacceptable magic AV found");
_break(aTHX_ rects, 0, (AV*)sv_2mortal((SV*)newAV()));