#include "multigrid.h"
/* if defined, residual is calculated and loop is aborted, if no improvement
- * is reached. this does not seem to work, so it should be undefined */
-//#define RESIDUAL
+ * is reached or the target resudual change is reached. */
+#define RESIDUAL
/* just process, but do not generate neighbor list and do not solve multigrid */
//#define TEST_NO_RENDER
#define NEIGHBORS_2D 9
#define NEIGHBORS_3D 27
-struct neighbors_2d {
- int num; /* # of neighbors */
- int index[NEIGHBORS_2D]; /* index of neighbors in an array */
- double weight[NEIGHBORS_2D]; /* calculated weight of neighbors */
-};
-struct neighbors_3d {
- int num; /* # of neighbors */
- int index[NEIGHBORS_3D]; /* index of neighbors in an array */
- double weight[NEIGHBORS_3D]; /* calculated weight of neighbors */
-};
-
static struct neigh_off {
int x, y, z;
} neigh_off[NEIGHBORS_3D] = {
}
/* create an array of neighbors and use flow vectors
+ *
+ * the elements of the array consist of:
+ * int num
+ * int index[neighbors]
+ * double weight[neighbors]
+ * num specifies the number of neighbors used.
+ * neighbors is fixed for all entries in the array and depends on 2d or 3d
+ * image. since the number of neighbors depends on the dimensions of an image,
+ * they are not fixed. a structure cannot be used in this case.
+ *
* flow_x and flow_x are the vectors that point to neighbors of the next z plane
* flow_x and flow_x are the vectors that point to neighbors of the previous z plane
*/
int index;
double *flow_x_ptr, *flow_y_ptr, *flow_ix_ptr, *flow_iy_ptr;
unsigned char *neighbors, *ret_neighbors;
+ int nb_num, index_offset;
unsigned long nb_size;
- int nb_num;
int *i_index;
ASSERT_FLOW(w, h, k, flow, flow_i);
flow_iy_ptr = flow_ix_ptr + w*h*(k-1);
nb_num = (k==1) ? NEIGHBORS_2D : NEIGHBORS_3D;
- nb_size = (k==1) ? sizeof(struct neighbors_2d) : sizeof(struct neighbors_3d);
+ nb_size = sizeof(int) + (sizeof(int) + sizeof(double)) * nb_num;
+ index_offset = sizeof(int);
ret_neighbors = neighbors = calloc(w*h*k, nb_size);
// printf("alloc neighbor grid=%lu Mbytes\n", nb_size*w*h*k/1024/1024);
if (!ret_neighbors) {
- printf("failed to allocate neighbor grid\n");
+ printf("failed to allocate neighbor grid (%lu megabytes)\n", nb_size * w*h*k / 1024 / 1024 );
return NULL;
}
for (z = 0; z < k; z++) {
for (y = 0; y < h; y++) {
for (x = 0; x < w; x++) {
- if (k==1)
- i_index = ((struct neighbors_2d *)neighbors)->index;
- else
- i_index = ((struct neighbors_3d *)neighbors)->index;
+ i_index = (int *)(neighbors+index_offset);
for (i = 0, n = 0; i < nb_num; i++) {
/* xx, yy, zz is the neighbor of color source
* The first neighbor is always the point itself!
continue;
i_index[n++] = index;
}
- if (k == 1)
- ((struct neighbors_2d *)neighbors)->num = n;
- else
- ((struct neighbors_3d *)neighbors)->num = n;
+ *((int *)neighbors) = n;
neighbors += nb_size;
}
}
double sqdiff[NEIGHBORS_3D];
double log0_01 = log(0.01);
double *lum_ptr;
+ int nb_num, index_offset, weight_offset;
unsigned long nb_size;
int *i_index; double *w_index;
k = darrayGetDimensions(luminance)[2];
lum_ptr = darrayGetPr(luminance);
- nb_size = (k==1) ? sizeof(struct neighbors_2d) : sizeof(struct neighbors_3d);
+ nb_num = (k==1) ? NEIGHBORS_2D : NEIGHBORS_3D;
+ nb_size = sizeof(int) + (sizeof(int) + sizeof(double)) * nb_num;
+ index_offset = sizeof(int);
+ weight_offset = sizeof(int) + sizeof(int) * nb_num;
for (z = 0; z < k; z++) {
for (y = 0; y < h; y++) {
for (x = 0; x < w; x++, neighbors += nb_size) {
- if (k==1) {
- num = ((struct neighbors_2d *)neighbors)->num;
- i_index = ((struct neighbors_2d *)neighbors)->index;
- w_index = ((struct neighbors_2d *)neighbors)->weight;
- } else {
- num = ((struct neighbors_3d *)neighbors)->num;
- i_index = ((struct neighbors_3d *)neighbors)->index;
- w_index = ((struct neighbors_3d *)neighbors)->weight;
- }
+ num = *((int *)neighbors);
+ i_index = (int *)(neighbors+index_offset);
+ w_index = (double *)(neighbors+weight_offset);
/* Be sure that we have at least 3 neighbors.
* This applies for every corner in a plane.
* The first point in the neighbor list iss the point itself.
int x, y, z, w, h, k, i, num;
double *value_ptr, *value_line_ptr, *mark_line_ptr, *init_line_ptr;
double val;
+ int nb_num, index_offset, weight_offset;
unsigned long nb_size;
int *i_index; double *w_index;
ASSERT(k == darrayGetDimensions(mark)[2], "z-sizes of mark array is different");
value_ptr = darrayGetPr(value);
- nb_size = (k==1) ? sizeof(struct neighbors_2d) : sizeof(struct neighbors_3d);
+ nb_num = (k==1) ? NEIGHBORS_2D : NEIGHBORS_3D;
+ nb_size = sizeof(int) + (sizeof(int) + sizeof(double)) * nb_num;
+ index_offset = sizeof(int);
+ weight_offset = sizeof(int) + sizeof(int) * nb_num;
for (z = 0; z < k; z++) {
for (y = 0; y < h; y++) {
mark_line_ptr = darrayGetPr(mark) + y*w + z*w*h;
init_line_ptr = darrayGetPr(init) + y*w + z*w*h;
for (x = 0; x < w; x++, neighbors += nb_size) {
- if (k==1) {
- num = ((struct neighbors_2d *)neighbors)->num;
- i_index = ((struct neighbors_2d *)neighbors)->index;
- w_index = ((struct neighbors_2d *)neighbors)->weight;
- } else {
- num = ((struct neighbors_3d *)neighbors)->num;
- i_index = ((struct neighbors_3d *)neighbors)->index;
- w_index = ((struct neighbors_3d *)neighbors)->weight;
- }
+ num = *((int *)neighbors);
+ i_index = (int *)(neighbors+index_offset);
+ w_index = (double *)(neighbors+weight_offset);
/* if pixle is marked */
darrayCheckPr(mark, mark_line_ptr + x, 1);
if (mark_line_ptr[x] < 0.5) {
int x, y, z, w, h, k, i, num;
double *value_ptr, *mark_line_ptr, *init_line_ptr;
double val, residual = 0;
+ int nb_num, index_offset, weight_offset;
unsigned long nb_size;
int *i_index; double *w_index;
ASSERT(k == darrayGetDimensions(mark)[2], "z-sizes of mark array is different");
value_ptr = darrayGetPr(value);
- nb_size = (k==1) ? sizeof(struct neighbors_2d) : sizeof(struct neighbors_3d);
+ nb_num = (k==1) ? NEIGHBORS_2D : NEIGHBORS_3D;
+ nb_size = sizeof(int) + (sizeof(int) + sizeof(double)) * nb_num;
+ index_offset = sizeof(int);
+ weight_offset = sizeof(int) + sizeof(int) * nb_num;
for (z = 0; z < k; z++) {
for (y = 0; y < h; y++) {
init_line_ptr = darrayGetPr(init) + y*w + z*w*h;
mark_line_ptr = darrayGetPr(mark) + y*w + z*w*h;
for (x = 0; x < w; x++, neighbors += nb_size) {
- if (k==1) {
- num = ((struct neighbors_2d *)neighbors)->num;
- i_index = ((struct neighbors_2d *)neighbors)->index;
- w_index = ((struct neighbors_2d *)neighbors)->weight;
- } else {
- num = ((struct neighbors_3d *)neighbors)->num;
- i_index = ((struct neighbors_3d *)neighbors)->index;
- w_index = ((struct neighbors_3d *)neighbors)->weight;
- }
+ num = *((int *)neighbors);
+ i_index = (int *)(neighbors+index_offset);
+ w_index = (double *)(neighbors+weight_offset);
/* if pixle is marked, no residual */
darrayCheckPr(mark, mark_line_ptr + x, 1);
if (mark_line_ptr[x] < 0.5)
* before z)!
*
*/
-int solve_mg(int nlevel, int iters, int ncycle, darray_t **values, darray_t **marks, const darray_t *init, unsigned char **nb_list, darray_t **flows, int scalexyz)
+int solve_mg(int nlevel, int iters, int ncycle, darray_t **values, darray_t **marks, const darray_t *init, unsigned char **nb_list, darray_t **flows, int scalexyz, double target_residual_change)
{
#ifdef RESIDUAL
double curr_residual, prev_residual = 10000000000000.0;
for (i = 0; i < iters; i++)
smooth(level, values[level], marks[level], init, nb_list[level]);
#ifdef RESIDUAL
- curr_residual += residual(level, values[level], marks[level], init, nb_list[level]);
+ if (target_residual_change)
+ curr_residual += residual(level, values[level], marks[level], init, nb_list[level]);
#endif
if (level < scalexyz)
rc = restrict(values[level], values[level+1], flows[level]);
/* end if residual raises again */
#ifdef RESIDUAL
- if (prev_residual < curr_residual)
+//printf("residual %.8f (change=%%%.4f)\n", curr_residual, (prev_residual-curr_residual)/prev_residual*100);
+ if (target_residual_change)
+ printf(" %.1f%%", (prev_residual-curr_residual)/prev_residual*100); fflush(stdout);
+ if (curr_residual == 0.0) {
+ printf(" 0"); fflush(stdout);
+ break;
+ }
+ if (target_residual_change && (prev_residual-curr_residual)/prev_residual < target_residual_change)
break;
- curr_residual = prev_residual;
+ prev_residual = curr_residual;
#endif
}
projects / colorize.git / blobdiff