sandbox/Antoonvh/raycasterv.c
Adaptive Ray caster
This program sends pixel coordinates to stdout, reads the correspoding RGB codes from stdin
, and writes a ppm image to stderr
.
#include "utils.h"
#include "my_vertex.h"
#include "adapt_field.h"
vector u;
void refine_r (Point point, scalar s) {
fine (s, 0, 0, 0) = s[];
fine (s, 2, 0, 0) = s[1];
fine (s, 0, 2, 0) = s[0,1];
fine (s, 2, 2, 0) = s[1,1];
fine (s, 1, 0, 0) = nodata;
fine (s, 1, 1, 0) = nodata;
fine (s, 0, 1, 0) = nodata;
fine (s, 1, 2, 0) = nodata;
fine (s, 2, 1, 0) = nodata;
}
void refine_gb (Point point, scalar s) {
fine (s, 0, 0, 0) = s[];
fine (s, 2, 0, 0) = s[1];
fine (s, 0, 2, 0) = s[0,1];
fine (s, 2, 2, 0) = s[1,1];
}
double interpv (scalar s, double xp, double yp) {
Point point = locate (xp, yp);
if (point.level > 0) {
double w00 = (xp - (x - Delta/2.))*(yp - (y - Delta/2.));
double w10 = ((x + Delta/2.) - xp)*(yp - (y - Delta/2.));
double w01 = (xp - (x - Delta/2.))*((y + Delta/2.) - yp);
double w11 = ((x + Delta/2.) - xp)*((y + Delta/2.) - yp);
return (w00*s[1,1] + w10*s[0,1] + w01*s[1] + w11*s[])/sq(Delta);
}
else return 0;
}
void get_w (scalar * list, scalar w) {
scalar * wl = list_clone (list);
scalar s, ws;
for (s, ws in list,wl)
wavelet (s, ws);
for (int l = depth() - 2; l < depth(); l++)
foreach_coarse_level (l) {
double maxw = 0.;
foreach_child() {
double a = 0;
for (scalar ws in wl)
a += fabs(ws[]);
maxw = max(maxw, a);
}
w[] = maxw;
}
foreach()
w[] = coarse(w,0,0,0);
}
int maxlevel = 9;
double tol = 5;
int main (int argc, char ** argv) {
if (argc > 1)
maxlevel = atoi(argv[1]);
if (argc > 2)
tol = atof(argv[2]);
#if _OPENMP
int threads = npe();
#endif
int s = sizeof(double);
int Ni = 1 << maxlevel;
init_grid (1 << (maxlevel - 2));
scalar w[];
scalar r[], g[], b[], * rgb = {r, g, b};
for (scalar s in rgb) {
s.refine = refine_gb;
s.prolongation = refine_vert;
s.restriction = restriction_vert;
}
r.refine = refine_r;
foreach_vertex()
r[] = nodata;
int count = 1;
while (count > 0) {
count = 0;
foreach_vertex(reduction(+:count))
if (r[] == nodata)
count++;
fwrite (&count, 1, sizeof(int), stdout);
fflush (stdout);
if (count > 0) {
#if _OPENMP
omp_set_num_threads(1);
#endif
foreach_vertex() {
if (r[] == nodata) {
fwrite (&x, s, 1, stdout);
fwrite (&y, s, 1, stdout);
}
}
fflush (stdout);
foreach_vertex() {
if (r[] == nodata) {
for (scalar s in rgb) {
unsigned char c;
read (STDIN_FILENO, &c, 1);
s[] = c;
}
}
}
#if _OPENMP
omp_set_num_threads(threads);
#endif
if (depth() == maxlevel - 2)
multigrid_restriction(rgb);
get_w (rgb, w);
adapt_field (w, tol, 0, maxlevel, 99, list = rgb);
}
}
fprintf (stderr, "P6\n%d %d\n%d\n", Ni, Ni, 255);
for (int j = 0; j < Ni; j++) {
double yp = L0 - (j + 0.5)*L0/Ni;
for (int i = 0; i < Ni; i++) {
double xp = (i + 0.5)*L0/Ni;
unsigned char px[3];
px[0] = interpv (r, xp, yp);
px[1] = interpv (g, xp, yp);
px[2] = interpv (b, xp, yp);
fwrite (px, sizeof(unsigned char), 3, stderr);
}
}
fflush (stderr);
}