#include "../src/mark.h"
#include "../src/ppm.h"
#include "../src/yuv.h"
+#include "../src/dir_seperator.h"
#include "../lib/darray.h"
#include "../lib/colorize.h"
#include "palette.h"
+#include "timeline.h"
#define min(x,y) ((x < y) ? x : y)
GtkWidget *window;
GdkPixbuf *pixbuf;
int width, height;
+ char filename[256];
+ char folder[256];
};
static void colorize_destroy(GtkWidget *widget, gpointer priv)
GTK_STOCK_SAVE, GTK_RESPONSE_ACCEPT,
NULL);
+ gtk_file_chooser_set_current_folder(GTK_FILE_CHOOSER(dialog), cp->folder);
+ gtk_file_chooser_set_current_name(GTK_FILE_CHOOSER(dialog), cp->filename);
if (gtk_dialog_run(GTK_DIALOG(dialog)) == GTK_RESPONSE_ACCEPT) {
filename = strdup(gtk_file_chooser_get_filename(GTK_FILE_CHOOSER(dialog)));
}
void colorize_image(void)
{
- darray_t *gI = NULL, *cI = NULL, *markIm = NULL, *ntscIm = NULL, *resultIm = NULL;
+ darray_t *gI = NULL, *cI = NULL, *markIm = NULL, *ntscIm = NULL;
double *ptr, *ptr2, *ptr3;
int dims[4];
int width = img_width, height = img_height;
int features, change_bc;
int i, j;
int c;
+ int rc;
GdkPixbuf *pixbuf;
GtkWidget *image;
guchar *data;
struct colorize_priv *cp;
int rs;
int zoom_field;
+ char *p, *q;
if (!img_pixbuf)
return;
cI = NULL;
/* render u and v change */
- resultIm = colorize(ntscIm, markIm, NULL, NULL, NULL, NULL, 5, 1, 0);
- if (!resultIm) {
+ rc = colorize(ntscIm, markIm, NULL, NULL, 5, 1, 0);
+ if (rc < 0) {
printerror("No memory! Use smaller image or add more memory.");
goto error;
}
- /* if we have a change, we apply brightness+contrast from resultIm */
+ /* if we have a change, we apply brightness+contrast from ntscIm */
if (change_bc) {
- ptr = darrayGetPr(resultIm);
- ptr2 = darrayGetPr(resultIm);
+ ptr = darrayGetPr(ntscIm);
+ ptr2 = darrayGetPr(ntscIm);
for (i = 0; i < img_height; i++) {
for (j = 0; j < img_width; j++) {
/* apply contrast */
/* apply YUV to pixbuffer and display */
pixbuf = gdk_pixbuf_new(GDK_COLORSPACE_RGB, FALSE, 8, img_width, img_height * zoom_field);
- ptr = darrayGetPr(resultIm);
+ ptr = darrayGetPr(ntscIm);
yuv2rgb(ptr, ptr, width, height);
data = gdk_pixbuf_get_pixels(pixbuf);
rs = gdk_pixbuf_get_rowstride(pixbuf);
image = gtk_image_new_from_pixbuf(pixbuf);
gtk_widget_show(image);
- darrayDestroy(resultIm);
darrayDestroy(markIm);
+ markIm = NULL;
darrayDestroy(ntscIm);
+ ntscIm = NULL;
darrayDone();
/* private structure */
- cp = malloc(sizeof(struct colorize_priv));
+ cp = calloc(sizeof(struct colorize_priv), 1);
/* create window */
colorize_window = gtk_window_new(GTK_WINDOW_TOPLEVEL);
cp->pixbuf = pixbuf;
cp->width = img_width*zoom_field;
cp->height = img_height*zoom_field;
+ p = frame_list[timeline_selected].filename;
+ while((q = strchr(p, DIR_SEPERATOR)))
+ p = q + 1;
+ strcpy(cp->folder, frame_list[timeline_selected].filename);
+ cp->folder[p - frame_list[timeline_selected].filename] = '\0';
+ sprintf(cp->filename, "colorized_%s", p);
/* create menu*/
menu_bar = gtk_menu_bar_new();
root_menu = gtk_menu_item_new_with_mnemonic("_File");
gtk_widget_show(root_menu);
- menu_item = gtk_image_menu_item_new_from_stock(GTK_STOCK_SAVE, NULL);
+ menu_item = gtk_image_menu_item_new_from_stock(GTK_STOCK_SAVE_AS, NULL);
gtk_menu_shell_append(GTK_MENU_SHELL(menu), menu_item);
g_signal_connect_swapped(menu_item, "activate", G_CALLBACK(save_event), cp);
gtk_widget_show(menu_item);
gtk_container_add(GTK_CONTAINER(colorize_window), vbox);
gtk_widget_show(colorize_window);
- return;
-
error:
darrayDestroy(gI);
darrayDestroy(cI);
- darrayDestroy(resultIm);
darrayDestroy(markIm);
darrayDestroy(ntscIm);
#include "timeline.h"
#include "image.h"
-static GtkWidget *flow_window = NULL;
+static GtkWidget *flow_dialog = NULL;
static GtkWidget *enable_button;
-static GtkObject *steps_adj;
-static GtkObject *win_adj;
-static GtkObject *max_adj;
+static GtkObject *window_adj;
+static GtkObject *vector_adj;
+static GtkObject *uv_adj;
static void flow_destroy(GtkWidget *widget, gpointer priv)
{
- gtk_widget_destroy(flow_window);
- flow_window = NULL;
+ gtk_widget_destroy(flow_dialog);
+ flow_dialog = NULL;
}
static gboolean flow_key_press(GtkWidget *widget, GdkEventKey *event, gpointer priv)
static void close_event(gpointer *priv)
{
- gtk_widget_destroy(flow_window);
- flow_window = NULL;
+ gtk_widget_destroy(flow_dialog);
+ flow_dialog = NULL;
}
static void reset_event(gpointer *priv)
anything_modified = 1;
flow_default();
- gtk_adjustment_set_value(GTK_ADJUSTMENT(steps_adj), flow_steps);
- gtk_adjustment_set_value(GTK_ADJUSTMENT(win_adj), flow_win);
- gtk_adjustment_set_value(GTK_ADJUSTMENT(max_adj), flow_max);
+ gtk_adjustment_set_value(GTK_ADJUSTMENT(window_adj), flow_window);
+ gtk_adjustment_set_value(GTK_ADJUSTMENT(vector_adj), flow_view_vector);
+ gtk_adjustment_set_value(GTK_ADJUSTMENT(uv_adj), flow_view_uv);
if (flowview)
timeline_select_and_save(timeline_selected, timeline_selected);
// draw_image(0, 0, -1, -1);
anything_modified = 1;
flow_enable = GTK_TOGGLE_BUTTON(enable_button)->active;
- flow_steps = GTK_ADJUSTMENT(steps_adj)->value;
- flow_win = GTK_ADJUSTMENT(win_adj)->value;
- flow_max = GTK_ADJUSTMENT(max_adj)->value;
+ flow_window = GTK_ADJUSTMENT(window_adj)->value;
+ flow_view_vector = GTK_ADJUSTMENT(vector_adj)->value;
+ flow_view_uv = GTK_ADJUSTMENT(uv_adj)->value;
if (flowview)
timeline_select_and_save(timeline_selected, timeline_selected);
// draw_image(0, 0, -1, -1);
GtkWidget *menu;
GtkWidget *menu_item;
GtkWidget *enable_label;
- GtkWidget *steps_label;
- GtkWidget *steps_scroll;
- GtkWidget *win_label;
- GtkWidget *win_scroll;
- GtkWidget *max_label;
- GtkWidget *max_scroll;
+ GtkWidget *window_label;
+ GtkWidget *window_scroll;
+ GtkWidget *vector_label;
+ GtkWidget *vector_scroll;
+ GtkWidget *uv_label;
+ GtkWidget *uv_scroll;
GtkWidget *button_apply;
GtkWidget *button_reset;
GtkWidget *box;
- if (flow_window) {
- gtk_widget_destroy(flow_window);
- flow_window = NULL;
+ if (flow_dialog) {
+ gtk_widget_destroy(flow_dialog);
+ flow_dialog = NULL;
return;
}
return;
/* create window */
- flow_window = gtk_window_new(GTK_WINDOW_TOPLEVEL);
- gtk_window_set_default_size(GTK_WINDOW(flow_window), 300, 0);
- gtk_window_set_title(GTK_WINDOW(flow_window), "Optical Flow");
- g_signal_connect(flow_window, "destroy", G_CALLBACK(flow_destroy), flow_window);
- g_signal_connect(flow_window, "key_press_event", G_CALLBACK(flow_key_press), flow_window);
+ flow_dialog = gtk_window_new(GTK_WINDOW_TOPLEVEL);
+ gtk_window_set_default_size(GTK_WINDOW(flow_dialog), 300, 0);
+ gtk_window_set_title(GTK_WINDOW(flow_dialog), "Optical Flow");
+ g_signal_connect(flow_dialog, "destroy", G_CALLBACK(flow_destroy), flow_dialog);
+ g_signal_connect(flow_dialog, "key_press_event", G_CALLBACK(flow_key_press), flow_dialog);
/* create menu*/
menu_bar = gtk_menu_bar_new();
menu_item = gtk_image_menu_item_new_from_stock(GTK_STOCK_CLOSE, NULL);
gtk_menu_shell_append(GTK_MENU_SHELL(menu), menu_item);
- g_signal_connect_swapped(menu_item, "activate", G_CALLBACK(close_event), flow_window);
+ g_signal_connect_swapped(menu_item, "activate", G_CALLBACK(close_event), flow_dialog);
gtk_widget_show(menu_item);
gtk_menu_item_set_submenu(GTK_MENU_ITEM(root_menu), menu);
gtk_menu_shell_append(GTK_MENU_SHELL(menu_bar), root_menu);
// g_signal_connect(enable_button, "toggled", G_CALLBACK(enable_event), NULL);
gtk_widget_show(enable_button);
- steps_label = gtk_label_new("Steps:");
- gtk_widget_show(GTK_WIDGET(steps_label));
- steps_adj = gtk_adjustment_new(flow_steps, 1, 31, 1, 1, 1);
- steps_scroll = gtk_hscale_new(GTK_ADJUSTMENT(steps_adj));
- gtk_widget_show(steps_scroll);
+ window_label = gtk_label_new("Search Window:");
+ gtk_widget_show(GTK_WIDGET(window_label));
+ window_adj = gtk_adjustment_new(flow_window, 3, 61, 1, 1, 1);
+ window_scroll = gtk_hscale_new(GTK_ADJUSTMENT(window_adj));
+ gtk_widget_show(window_scroll);
- win_label = gtk_label_new("Search Window:");
- gtk_widget_show(GTK_WIDGET(win_label));
- win_adj = gtk_adjustment_new(flow_win, 3, 61, 1, 1, 1);
- win_scroll = gtk_hscale_new(GTK_ADJUSTMENT(win_adj));
- gtk_widget_show(win_scroll);
+ vector_label = gtk_label_new("Display vectors (1=disable):");
+ gtk_widget_show(GTK_WIDGET(vector_label));
+ vector_adj = gtk_adjustment_new(flow_view_vector, 1, 51, 1, 1, 1);
+ vector_scroll = gtk_hscale_new(GTK_ADJUSTMENT(vector_adj));
+ gtk_widget_show(vector_scroll);
- max_label = gtk_label_new("Maximum flow vector length:");
- gtk_widget_show(GTK_WIDGET(max_label));
- max_adj = gtk_adjustment_new(flow_max, 0, 61, 1, 1, 1);
- max_scroll = gtk_hscale_new(GTK_ADJUSTMENT(max_adj));
- gtk_widget_show(max_scroll);
+ uv_label = gtk_label_new("Display next as UV (0=disable):");
+ gtk_widget_show(GTK_WIDGET(uv_label));
+ uv_adj = gtk_adjustment_new(flow_view_uv, 0, 101, 1, 1, 1);
+ uv_scroll = gtk_hscale_new(GTK_ADJUSTMENT(uv_adj));
+ gtk_widget_show(uv_scroll);
/* buttons */
button_reset = gtk_button_new();
gtk_box_pack_start(GTK_BOX(vbox), menu_bar, FALSE, FALSE, 2);
gtk_box_pack_start(GTK_BOX(vbox), enable_label, FALSE, FALSE, 2);
gtk_box_pack_start(GTK_BOX(vbox), enable_button, FALSE, FALSE, 2);
- gtk_box_pack_start(GTK_BOX(vbox), steps_label, FALSE, FALSE, 2);
- gtk_box_pack_start(GTK_BOX(vbox), steps_scroll, FALSE, FALSE, 2);
- gtk_box_pack_start(GTK_BOX(vbox), win_label, FALSE, FALSE, 2);
- gtk_box_pack_start(GTK_BOX(vbox), win_scroll, FALSE, FALSE, 2);
- gtk_box_pack_start(GTK_BOX(vbox), max_label, FALSE, FALSE, 2);
- gtk_box_pack_start(GTK_BOX(vbox), max_scroll, FALSE, FALSE, 2);
+ gtk_box_pack_start(GTK_BOX(vbox), window_label, FALSE, FALSE, 2);
+ gtk_box_pack_start(GTK_BOX(vbox), window_scroll, FALSE, FALSE, 2);
+ gtk_box_pack_start(GTK_BOX(vbox), vector_label, FALSE, FALSE, 2);
+ gtk_box_pack_start(GTK_BOX(vbox), vector_scroll, FALSE, FALSE, 2);
+ gtk_box_pack_start(GTK_BOX(vbox), uv_label, FALSE, FALSE, 2);
+ gtk_box_pack_start(GTK_BOX(vbox), uv_scroll, FALSE, FALSE, 2);
gtk_box_pack_start(GTK_BOX(vbox), hbox, FALSE, FALSE, 20);
gtk_widget_show(vbox);
- gtk_container_add(GTK_CONTAINER(flow_window), vbox);
- gtk_widget_show(flow_window);
+ gtk_container_add(GTK_CONTAINER(flow_dialog), vbox);
+ gtk_widget_show(flow_dialog);
}
no_mark:
/* create pixbuf */
- create_or_reset_pixbuf();
+ create_or_reset_pixbuf(1);
}
-void create_or_reset_pixbuf(void)
+void create_or_reset_pixbuf(int resize)
{
+ int w, h, dw = 0, dh = 0;
+ GdkScreen *screen = gtk_window_get_screen(GTK_WINDOW(main_window));
+ GdkRectangle max;
+
+ /* calculate the maximum windows size, so that the border will not exceed the screen size */
+ gdk_window_get_frame_extents(gtk_widget_get_window(main_window), &max);
+ max.width -= main_window->allocation.width;
+ max.height -= main_window->allocation.height;
+ if (max.width < 0)
+ max.width = 0;
+ if (max.height < 0)
+ max.height = 0;
+ max.width = gdk_screen_get_width(screen)-max.width;
+ max.height = gdk_screen_get_height(screen)-max.height;
+
+ if (resize && screen) {
+try_smaller:
+ w = img_scroll->allocation.width-20;
+ h = img_scroll->allocation.height-20;
+ if (w < img_width*img_scale_x/16) {
+ w = main_window->allocation.width - w + img_width*img_scale_x/16;
+ if (w > max.width) {
+ if (img_scale_x > 4) {
+ int aspect = img_scale_y / img_scale_x;
+ img_scale_x -= 4;
+ img_scale_y = img_scale_x * aspect;
+ goto try_smaller;
+ }
+ w = max.width;
+ }
+ dw = w - main_window->allocation.width;
+ }
+ if (h < img_height*img_scale_y/16) {
+ h = main_window->allocation.height - h + img_height*img_scale_y/16;
+ if (h > max.height) {
+ if (img_scale_x > 4) {
+ int aspect = img_scale_y / img_scale_x;
+ img_scale_x -= 4;
+ img_scale_y = img_scale_x * aspect;
+ goto try_smaller;
+ }
+ h = max.height;
+ }
+ dh = h - main_window->allocation.height;
+ }
+// printf("%d %d\n", gdk_screen_get_width(screen), gdk_screen_get_height(screen));
+ if (dw || dh)
+ gtk_window_resize(GTK_WINDOW(main_window), main_window->allocation.width+dw, main_window->allocation.height+dh);
+ }
+
if (img_pixbuf) {
g_object_unref(img_pixbuf);
img_pixbuf = NULL;
img_pixbuf = gdk_pixbuf_new(GDK_COLORSPACE_RGB, FALSE, 8, img_width*img_scale_x/16, img_height*img_scale_y/16);
gtk_drawing_area_size(GTK_DRAWING_AREA(img_drawing_area), img_width*img_scale_x/16, img_height*img_scale_y/16);
-// gtk_widget_set_size_request(drawing_area, img_width, img_height);
/* label */
timeline_show_name(timeline_selected);
extern int anything_modified;
void create_image(const char *filename);
-void create_or_reset_pixbuf(void);
+void create_or_reset_pixbuf(int resize);
void destroy_image(void);
void draw_image(int x, int y, int w, int h);
void paint_brush(int x, int y, int size, int paint);
#define CC_APP_TITLE "Colorize GTK"
GtkWidget *main_window;
+GtkWidget *img_scroll;
GtkWidget *img_drawing_area = NULL;
GtkWidget *timeline_drawing_area;
GtkWidget *timeline_scroll;
GtkWidget *dialog;
char *filename = NULL;
+ if (anything_modified) {
+ int ret;
+ GtkWidget *dialog = gtk_message_dialog_new (GTK_WINDOW(main_window),
+ GTK_DIALOG_DESTROY_WITH_PARENT,
+ GTK_MESSAGE_ERROR,
+ GTK_BUTTONS_OK_CANCEL,
+ "Image has been changed, really open new image?");
+ ret = gtk_dialog_run(GTK_DIALOG (dialog));
+ gtk_widget_destroy(dialog);
+ if (ret != GTK_RESPONSE_OK)
+ return;
+ }
+
dialog = gtk_file_chooser_dialog_new("Select first grey image",
GTK_WINDOW(main_window),
GTK_FILE_CHOOSER_ACTION_OPEN,
GTK_STOCK_OPEN, GTK_RESPONSE_ACCEPT,
NULL);
+ if (frame_list)
+ gtk_file_chooser_set_filename(GTK_FILE_CHOOSER(dialog), frame_list[timeline_selected].filename);
if (gtk_dialog_run(GTK_DIALOG(dialog)) == GTK_RESPONSE_ACCEPT) {
filename = strdup(gtk_file_chooser_get_filename(GTK_FILE_CHOOSER(dialog)));
}
gtk_widget_destroy(dialog);
if (filename) {
- timeline_select_and_save(timeline_selected, timeline_selected);
+// do not save timeline_select_and_save(timeline_selected, timeline_selected);
create_timeline(filename);
if (frame_list)
int aspect = img_scale_y / img_scale_x;
img_scale_x += 4;
img_scale_y = img_scale_x * aspect;
- create_or_reset_pixbuf();
+ create_or_reset_pixbuf(0);
}
}
int aspect = img_scale_y / img_scale_x;
img_scale_x -= 4;
img_scale_y = img_scale_x * aspect;
- create_or_reset_pixbuf();
+ create_or_reset_pixbuf(0);
}
}
{
img_scale_x = 16;
img_scale_y = 16;
- create_or_reset_pixbuf();
+ create_or_reset_pixbuf(0);
}
static void zoom_field_event(gpointer priv)
img_scale_y = img_scale_x*2;
else
img_scale_y = img_scale_x;
- create_or_reset_pixbuf();
+ create_or_reset_pixbuf(1);
}
static void show_highlighted_event(gpointer priv)
int main(int argc, char *argv[])
{
GtkWidget *vbox, *tool_box;
- GtkWidget *scroll;
+ GtkWidget *paned;
+ GtkWidget *pal_scroll;
GtkWidget *menu_bar;
GtkWidget *separator;
GtkTreeSelection *selection;
main_window = gtk_window_new(GTK_WINDOW_TOPLEVEL);
gtk_widget_set_size_request(GTK_WIDGET(main_window), 400, 200);
- gtk_window_set_default_size(GTK_WINDOW(main_window), 800, 600);
+ gtk_window_set_default_size(GTK_WINDOW(main_window), 800, 500);
gtk_window_set_title(GTK_WINDOW(main_window), CC_APP_TITLE " Version "
#include "../version.h"
);
}
/* add paned view to vbox (middle part of vbox) */
- GtkWidget *paned = gtk_hpaned_new ();
+ paned = gtk_hpaned_new ();
gtk_widget_show(paned);
gtk_box_pack_start(GTK_BOX(vbox), paned, TRUE, TRUE, 2);
gtk_paned_set_position(GTK_PANED (paned), 250);
gtk_widget_set_events(palette_treeview, GDK_ENTER_NOTIFY_MASK | GDK_LEAVE_NOTIFY_MASK);
gtk_widget_show(palette_treeview);
- scroll = gtk_scrolled_window_new(NULL, NULL);
- gtk_scrolled_window_set_policy(GTK_SCROLLED_WINDOW(scroll), GTK_POLICY_AUTOMATIC, GTK_POLICY_ALWAYS);
- gtk_widget_show(scroll);
- gtk_container_add(GTK_CONTAINER(scroll), palette_treeview);
+ pal_scroll = gtk_scrolled_window_new(NULL, NULL);
+ gtk_scrolled_window_set_policy(GTK_SCROLLED_WINDOW(pal_scroll), GTK_POLICY_AUTOMATIC, GTK_POLICY_ALWAYS);
+ gtk_widget_show(pal_scroll);
+ gtk_container_add(GTK_CONTAINER(pal_scroll), palette_treeview);
palette_renderer = gtk_cell_renderer_text_new();
g_object_set(G_OBJECT(palette_renderer), "editable", FALSE, NULL);
selection = gtk_tree_view_get_selection(GTK_TREE_VIEW(palette_treeview));
g_signal_connect(selection, "changed", G_CALLBACK(palette_change), NULL);
- gtk_paned_pack1 (GTK_PANED (paned), scroll, TRUE, TRUE);
+ gtk_paned_pack1 (GTK_PANED (paned), pal_scroll, FALSE, TRUE);
/* add drawing area (right of hbox) */
img_drawing_area = gtk_drawing_area_new ();
| GDK_POINTER_MOTION_HINT_MASK);
gtk_widget_show(img_drawing_area);
- scroll = gtk_scrolled_window_new(NULL, NULL);
- gtk_scrolled_window_set_policy(GTK_SCROLLED_WINDOW(scroll), GTK_POLICY_AUTOMATIC, GTK_POLICY_AUTOMATIC);
- gtk_widget_show(scroll);
- gtk_scrolled_window_add_with_viewport(GTK_SCROLLED_WINDOW(scroll), img_drawing_area);
- gtk_paned_pack2 (GTK_PANED (paned), scroll, TRUE, TRUE);
+ img_scroll = gtk_scrolled_window_new(NULL, NULL);
+ gtk_scrolled_window_set_policy(GTK_SCROLLED_WINDOW(img_scroll), GTK_POLICY_AUTOMATIC, GTK_POLICY_AUTOMATIC);
+ gtk_widget_show(img_scroll);
+ gtk_scrolled_window_add_with_viewport(GTK_SCROLLED_WINDOW(img_scroll), img_drawing_area);
+ gtk_paned_pack2 (GTK_PANED (paned), img_scroll, TRUE, TRUE);
#if 0
/* add horizontal seperation to vbox */
void main_destroy(void);
extern int highlight, preview, rendered, flowview;
extern GtkWidget *main_window;
+extern GtkWidget *img_scroll;
extern GtkWidget *img_drawing_area;
extern GtkWidget *timeline_drawing_area;
extern GtkWidget *palette_treeview;
(new_frame < timeline_frames - 1) ? frame_list[new_frame + 1].filename : NULL,
img_grey_buffer,
img_width, img_height,
- flow_steps,
- flow_win,
- flow_max);
+ flow_window,
+ flow_view_vector,
+ flow_view_uv);
#endif
if (timeline_frames > 1) {
libcolorize_la_SOURCES = \
darray.c \
- multigrid.c
+ multigrid.c \
+ colorize.c
--- /dev/null
+#include <stdio.h>
+#include <string.h>
+#include <stdlib.h>
+#include <math.h>
+#include "darray.h"
+#include "multigrid.h"
+#include "colorize.h"
+
+/* colorize image array using image_mark
+ *
+ * image:
+ * first dimansion = width
+ * second dimension = height
+ * thrid dimension = components (first component is luminance, other components are features, like u, v)
+ * fourth dimension = time
+ *
+ * image_mark:
+ * first dimansion = width
+ * second dimension = height
+ * fourth dimension = time
+ */
+int colorize_prepare(struct colorize *col)
+{
+ int dimensions[4];
+ int ww, hh, kk, last_kk, d, temp, z;
+ darray_t *luminance = NULL;
+ int rc = -1;
+
+ ASSERT(darrayGetNumberOfDimensions(col->image) == 4, "image array must have 4 dimensions");
+ col->w = darrayGetDimensions(col->image)[0];
+ col->h = darrayGetDimensions(col->image)[1];
+ col->comp = darrayGetDimensions(col->image)[2];
+ col->k = darrayGetDimensions(col->image)[3];
+ ASSERT(col->comp > 0, "image array requires at least one component (Y-component)");
+
+ ASSERT(darrayGetNumberOfDimensions(col->image_mark) == 3, "image_mark array must have 3 dimensions");
+ ASSERT(col->w == darrayGetDimensions(col->image_mark)[0], "image_mark and image must both have same x-size");
+ ASSERT(col->h == darrayGetDimensions(col->image_mark)[1], "image_mark and image must both have same y-size");
+ ASSERT(col->k == darrayGetDimensions(col->image_mark)[2], "image_mark and image must both have same z-size");
+
+ /* Calculate maximum depth
+ * Scale down dimensions until at least 2x2 remains
+ */
+ temp = min(col->w, col->h);
+ for (col->max_depth = 0; temp > 2; col->max_depth++)
+ temp=(temp+1)>>1;
+ ASSERT(col->max_depth >= 2, "image too small");
+
+ col->nb_list = calloc(col->max_depth, sizeof(struct neighbors *));
+ if (!col->nb_list)
+ goto error;
+ col->marks = calloc(col->max_depth, sizeof(darray_t *));
+ if (!col->marks)
+ goto error;
+ col->values = calloc(col->max_depth, sizeof(darray_t *));
+ if (!col->values)
+ goto error;
+ col->flows = calloc(col->max_depth, sizeof(darray_t *));
+ if (!col->flows)
+ goto error;
+ col->flows_i = calloc(col->max_depth, sizeof(darray_t *));
+ if (!col->flows_i)
+ goto error;
+ if (col->scalexyz) {
+ /* the flows_xy are scaled in xy in the current level */
+ col->flows_xy = calloc(col->max_depth, sizeof(darray_t *));
+ if (!col->flows_xy)
+ goto error;
+ }
+ dimensions[0] = col->w;
+ dimensions[1] = col->h;
+ dimensions[2] = col->k;
+ luminance = darrayCreate(3, dimensions);
+ if (!luminance)
+ goto error;
+ /* generate sparse multigrid */
+ last_kk = kk = col->k;
+ printf("-> prepare sparse grid");
+ for (d = 0, ww=col->w, hh=col->h; d < col->max_depth; d++, ww=(ww+1)>>1, hh=(hh+1)>>1) {
+ /* create temporary arrays
+ * copy fist array to level 0 or shrink previous level
+ */
+ printf(" #%d(w=%d h=%d k=%d)", d, ww, hh, kk); fflush(stdout);
+ dimensions[0] = ww;
+ dimensions[1] = hh;
+ dimensions[2] = last_kk;
+ col->marks[d] = darrayCreate(3, dimensions);
+ if (!col->marks[d])
+ goto error;
+ dimensions[0] = ww;
+ dimensions[1] = hh;
+ dimensions[2] = last_kk-1;
+ dimensions[3] = 2;
+ col->flows[d] = darrayCreate(4, dimensions);
+ if (!col->flows[d])
+ goto error;
+ col->flows_i[d] = darrayCreate(4, dimensions);
+ if (!col->flows_i[d])
+ goto error;
+ if (d < col->scalexyz) {
+ dimensions[0] = (ww+1)>>1;
+ dimensions[1] = (hh+1)>>1;
+ dimensions[2] = kk-1;
+ dimensions[3] = 2;
+ col->flows_xy[d] = darrayCreate(4, dimensions);
+ if (!col->flows_xy[d])
+ goto error;
+ }
+ if (d == 0) {
+ for (z = 0; z < kk; z++) {
+ darrayCheckPr(luminance, darrayGetPr(luminance) + z*col->w*col->h, col->w*col->h);
+ darrayCheckPr(col->image, darrayGetPr(col->image) + z*col->w*col->h*col->comp, col->w*col->h);
+ memcpy(darrayGetPr(luminance) + z*col->w*col->h, darrayGetPr(col->image) + z*col->w*col->h*col->comp, col->w*col->h*sizeof(double));
+ }
+ if (col->flow && col->flow_i) {
+ darrayCheckPr(col->flows[d], darrayGetPr(col->flows[d]), col->w*col->h*(kk-1)*2);
+ darrayCheckPr(col->flow, darrayGetPr(col->flow), col->w*col->h*(kk-1)*2);
+ memcpy(darrayGetPr(col->flows[d]), darrayGetPr(col->flow), col->w*col->h*(kk-1)*2*sizeof(double));
+ darrayCheckPr(col->flows_i[d], darrayGetPr(col->flows_i[d]), col->w*col->h*(kk-1)*2);
+ darrayCheckPr(col->flow_i, darrayGetPr(col->flow_i), col->w*col->h*(kk-1)*2);
+ memcpy(darrayGetPr(col->flows_i[d]), darrayGetPr(col->flow_i), col->w*col->h*(kk-1)*2*sizeof(double));
+ }
+ darrayCheckPr(col->marks[d], darrayGetPr(col->marks[d]), col->w*col->h*kk);
+ darrayCheckPr(col->image_mark, darrayGetPr(col->image_mark), col->w*col->h*kk);
+ memcpy(darrayGetPr(col->marks[d]), darrayGetPr(col->image_mark), col->w*col->h*kk*sizeof(double));
+ /* invert, so 0 = marked */
+ quantize(col->marks[d], 0.5, 0, 1);
+ } else {
+ rc = shrink_xy(luminance, NULL);
+ if (rc < 0)
+ goto error;
+ rc = shrink_xy(col->marks[d-1], col->marks[d]);
+ if (rc < 0)
+ goto error;
+ rc = shrink_xy(col->flows[d-1], col->flows[d]);
+ if (rc < 0)
+ goto error;
+ rc = shrink_xy(col->flows_i[d-1], col->flows_i[d]);
+ if (rc < 0)
+ goto error;
+ /* shrink length of flow vector on smaller grid.
+ * This is not performed when comparing with
+ * the original algorithm. (it is a bug there)
+ */
+ multiply(col->flows[d], 0.5);
+ multiply(col->flows_i[d], 0.5);
+ if (d-1 < col->scalexyz) { /* process if previous d shall be shrinked in z direction */
+ /* shrink time */
+ rc = shrink_z(luminance, NULL, col->flows[d]);
+ if (rc < 0)
+ goto error;
+ rc = shrink_z(col->marks[d], NULL, col->flows[d]);
+ if (rc < 0)
+ goto error;
+ rc = shrink_flow_z(col->flows[d], NULL);
+ if (rc < 0)
+ goto error;
+ rc = shrink_flow_z(col->flows_i[d], NULL);
+ if (rc < 0)
+ goto error;
+ }
+ /* If any of the restricted pixles is 0 (marked),
+ * the result must be 0 (marked) too.
+ */
+ quantize(col->marks[d], 0.9999, 1, 0);
+ }
+ if (d < col->scalexyz) {
+ /* create x and y shrinked flow vectors for restriction and prolonging process */
+ rc = shrink_xy(col->flows[d], col->flows_xy[d]);
+ if (rc < 0)
+ goto error;
+ multiply(col->flows_xy[d], 0.5);
+ }
+
+ /* generate values array */
+ dimensions[0] = ww;
+ dimensions[1] = hh;
+ dimensions[2] = kk;
+ col->values[d] = darrayCreate(3, dimensions);
+ if (!col->values[d])
+ goto error;
+
+ /* generate array of neighbors */
+#ifndef TEST_NO_RENDER
+ col->nb_list[d] = gen_neighbor(ww, hh, kk, col->flows[d], col->flows_i[d]);
+ if (!col->nb_list[d])
+ goto error;
+ /* weighten */
+ weighten(luminance, col->nb_list[d]);
+#endif
+
+ last_kk = kk;
+ if (d < col->scalexyz)
+ kk = (kk+1)>>1;
+ }
+ darrayDestroy(luminance);
+ luminance = NULL;
+ printf("\n");
+
+ dimensions[0] = col->w;
+ dimensions[1] = col->h;
+ dimensions[2] = col->k;
+ col->init = darrayCreate(3, dimensions);
+ if (!col->init)
+ goto error;
+
+#ifdef TEST_PYRAMID
+ dimensions[0] = col->w;
+ dimensions[1] = col->h;
+ dimensions[2] = col->k;
+ luminance = darrayCreate(3, dimensions);
+ for (z = 0; z < col->k; z++) {
+ darrayCheckPr(luminance, darrayGetPr(luminance) + z*col->w*col->h, col->w*col->h);
+ darrayCheckPr(col->image, darrayGetPr(col->image) + z*col->w*col->h*comp, col->w*col->h);
+ memcpy(darrayGetPr(luminance) + z*col->w*col->h, darrayGetPr(col->image) + z*col->w*col->h*comp, col->w*col->h*sizeof(double));
+ darrayCheckPr(col->values[0], darrayGetPr(col->values[0]) + z*col->w*col->h, col->w*col->h);
+ darrayCheckPr(col->image, darrayGetPr(col->image) + z*col->w*col->h*comp, col->w*col->h);
+ memcpy(darrayGetPr(col->values[0]) + z*col->w*col->h, darrayGetPr(col->image) + z*col->w*col->h*comp, col->w*col->h*sizeof(double));
+#ifdef TEST_FLOW
+ draw_flow(col->values[0], col->flows[0], col->flows_i[0]);
+#endif
+#ifdef TEST_MARKED
+ paste(col->marks[0], col->values[0], 0, 0, 0);
+#endif
+
+ }
+ for (d = 1; d < col->max_depth; d++) {
+ rc = shrink_xy(luminance, NULL);
+ if (rc < 0)
+ goto error;
+ if (d-1 < col->scalexyz)
+ shrink_z(luminance, NULL, col->flows_xy[d-1]);
+ paste(luminance, col->values[0], 0, 0, 0);
+#ifdef TEST_FLOW
+ draw_flow(col->values[0], col->flows_x[d], col->flows_ix[d]);
+#endif
+#ifdef TEST_MARKED
+ paste(col->marks[d], col->values[0], 0, 0, 0);
+#endif
+ }
+ darrayDestroy(luminance);
+ luminance = NULL;
+ /* copy test pyramid image as Y-component */
+ for (z = 0; z < col->k; z++) {
+ darrayCheckPr(col->image, darrayGetPr(col->image) + z*col->w*col->h*col->comp, col->w*col->h);
+ darrayCheckPr(col->values[0], darrayGetPr(col->values[0]) + z*col->w*col->h, col->w*col->h);
+ memcpy(darrayGetPr(col->image) + z*col->w*col->h*col->comp, darrayGetPr(col->values[0]) + z*col->w*col->h, col->w*col->h*sizeof(double));
+ }
+#endif
+
+ rc = 0;
+
+error:
+ darrayDestroy(luminance);
+
+ return rc;
+}
+
+int colorize_solve(struct colorize *col, int quick)
+{
+ int z, d, c, ww, hh, kk;
+ int rc = -1;
+
+ if (!quick)
+ printf("-> apply component to grid and solve");
+ for (c = 1; c < col->comp; c++) {
+ int pixles, j;
+ double *value_ptr, *init_ptr, *mark_ptr;
+
+ printf(" #%d", c); fflush(stdout);
+ /* apply component
+ * - copy component into value array (level 0), otherwise 0
+ * - use maked colors for init array, otherwise use 0
+ */
+ for (z = 0; z < col->k; z++) {
+ darrayCheckPr(col->values[0], darrayGetPr(col->values[0]) + z*col->w*col->h, col->w*col->h);
+ darrayCheckPr(col->image, darrayGetPr(col->image) + c*col->w*col->h + z*col->w*col->h*col->comp, col->w*col->h);
+ memcpy(darrayGetPr(col->values[0]) + z*col->w*col->h, darrayGetPr(col->image) + c*col->w*col->h + z*col->w*col->h*col->comp, col->w*col->h*sizeof(double));
+ }
+ for (d = 1, ww=(col->w+1)>>1, hh=(col->h+1)>>1, kk=col->k; d < col->max_depth; d++, ww=(ww+1)>>1, hh=(hh+1)>>1) {
+ if (d-1 < col->scalexyz) /* scale kk, if last array was scaled */
+ kk=(kk+1)>>1;
+ value_ptr = darrayGetPr(col->values[d]);
+ pixles = ww*hh*kk;
+ for (j = 0; j < pixles; j++) {
+ darrayCheckPr(col->values[d], value_ptr, 1);
+ *value_ptr++ = 0.0;
+ }
+ }
+ value_ptr = darrayGetPr(col->values[0]);
+ init_ptr = darrayGetPr(col->init);
+ mark_ptr = darrayGetPr(col->image_mark);
+ pixles = col->w*col->h*col->k;
+ for (j = 0; j < pixles; j++) {
+ darrayCheckPr(col->image_mark, mark_ptr, 1);
+ if (*mark_ptr++ > 0.5) {
+ darrayCheckPr(col->init, init_ptr, 1);
+ darrayCheckPr(col->values[0], value_ptr, 1);
+ *init_ptr++ = *value_ptr;
+ } else {
+ darrayCheckPr(col->init, init_ptr, 1);
+ *init_ptr++ = 0.0;
+ }
+ value_ptr++;
+ }
+
+#ifndef TEST_NO_RENDER
+ /* colorize component */
+ if (quick) {
+ rc = solve_mg(col->max_depth, col->inner_iter, 1, col->values, col->marks, col->init, col->nb_list, col->flows_xy, col->scalexyz);
+ if (rc < 0)
+ goto error;
+ } else {
+ int i;
+ for (i = 0; i < col->outer_iter; i++) {
+ rc = solve_mg(col->max_depth, col->inner_iter, 2, col->values, col->marks, col->init, col->nb_list, col->flows_xy, col->scalexyz);
+ if (rc < 0)
+ goto error;
+ rc = solve_mg((col->max_depth+1) / 2, col->inner_iter, 2, col->values, col->marks, col->init, col->nb_list, col->flows_xy, col->scalexyz);
+ if (rc < 0)
+ goto error;
+ rc = solve_mg(2, col->inner_iter, 2, col->values, col->marks, col->init, col->nb_list, col->flows_xy, col->scalexyz);
+ if (rc < 0)
+ goto error;
+ rc = solve_mg(1, col->inner_iter, 4, col->values, col->marks, col->init, col->nb_list, col->flows_xy, col->scalexyz);
+ if (rc < 0)
+ goto error;
+ }
+ }
+#endif
+
+#ifdef TEST_PYRAMID
+ for (d = 1; d < col->max_depth; d++)
+ paste(col->values[d], col->values[0], 0, 0, 0);
+#endif
+
+ /* apply colorized component */
+ for (z = 0; z < col->k; z++) {
+ darrayCheckPr(col->image, darrayGetPr(col->image) + c*col->w*col->h + z*col->w*col->h*col->comp, col->w*col->h);
+ darrayCheckPr(col->values[0], darrayGetPr(col->values[0]) + z*col->w*col->h, col->w*col->h);
+ memcpy(darrayGetPr(col->image) + z*col->w*col->h*col->comp + c*col->w*col->h, darrayGetPr(col->values[0]) + z*col->w*col->h, col->w*col->h*sizeof(double));
+ }
+ }
+ if (!quick)
+ printf("\n");
+
+#ifdef TEST_SHRINK_XY
+ for (c = 0; c < col->comp; c++) {
+ darray_t *temp;
+ int dimensions[3] = {col->w, col->h, col->k};
+ for (z = 0; z < col->k; z++)
+ memcpy(darrayGetPr(col->values[0]) + z*col->w*col->h, darrayGetPr(col->image) + z*col->w*col->h*col->comp + c*col->w*col->h, col->w*col->h*sizeof(double));
+ dimensions[0] = (col->w+1)>>1;
+ dimensions[1] = (col->h+1)>>1;
+ temp = darrayCreate(3, dimensions);
+ ASSERT(temp, "no memory");
+ rc = shrink_xy(col->values[0], temp);
+ ASSERT(rc >= 0, "no memory");
+ paste(temp, col->values[0], 0, 0, 0);
+ darrayDestroy(temp);
+ for (z = 0; z < col->k; z++)
+ memcpy(darrayGetPr(col->image) + z*col->w*col->h*col->comp + c*col->w*col->h, darrayGetPr(col->values[0]) + z*col->w*col->h, col->w*col->h*sizeof(double));
+ }
+#endif
+
+#ifdef TEST_SHRINK_Z
+ for (c = 0; c < col->comp; c++) {
+ darray_t *temp;
+ int dimensions[3] = {col->w, col->h, col->k};
+ for (z = 0; z < col->k; z++)
+ memcpy(darrayGetPr(col->values[0]) + z*col->w*col->h, darrayGetPr(col->image) + z*col->w*col->h*col->comp + c*col->w*col->h, col->w*col->h*sizeof(double));
+ dimensions[2] = (col->k+1)>>1;
+ temp = darrayCreate(3, dimensions);
+ ASSERT(temp, "no memory");
+ rc = shrink_z(col->values[0], temp, col->flows[0]);
+ ASSERT(rc >= 0, "no memory");
+ paste(temp, col->values[0], 0, 0, 0);
+ darrayDestroy(temp);
+ for (z = 0; z < col->k; z++)
+ memcpy(darrayGetPr(col->image) + z*col->w*col->h*col->comp + c*col->w*col->h, darrayGetPr(col->values[0]) + z*col->w*col->h, col->w*col->h*sizeof(double));
+ }
+#endif
+
+#ifdef TEST_EXPAND_XY
+ for (c = 0; c < col->comp; c++) {
+ darray_t *temp;
+ int dimensions[3] = {col->w, col->h, col->k};
+ for (z = 0; z < col->k; z++)
+ memcpy(darrayGetPr(col->values[0]) + z*col->w*col->h, darrayGetPr(col->image) + z*col->w*col->h*col->comp + c*col->w*col->h, col->w*col->h*sizeof(double));
+ dimensions[0] = (col->w+1)>>1;
+ dimensions[1] = (col->h+1)>>1;
+ temp = darrayCreate(3, dimensions);
+ ASSERT(temp, "no memory");
+ rc = shrink_xy(col->values[0], temp);
+ ASSERT(rc >= 0, "no memory");
+ expand_xy(temp, col->values[0]);
+ darrayDestroy(temp);
+ for (z = 0; z < col->k; z++)
+ memcpy(darrayGetPr(col->image) + z*col->w*col->h*col->comp + c*col->w*col->h, darrayGetPr(col->values[0]) + z*col->w*col->h, col->w*col->h*sizeof(double));
+ }
+#endif
+
+#ifdef TEST_EXPAND_Z
+ for (c = 0; c < col->comp; c++) {
+ darray_t *temp;
+ int dimensions[3] = {col->w, col->h, col->k};
+ for (z = 0; z < col->k; z++)
+ memcpy(darrayGetPr(col->values[0]) + z*col->w*col->h, darrayGetPr(col->image) + z*col->w*col->h*col->comp + c*col->w*col->h, col->w*col->h*sizeof(double));
+ dimensions[2] = (col->k+1)>>1;
+ temp = darrayCreate(3, dimensions);
+ ASSERT(temp, "no memory");
+ rc = shrink_z(col->values[0], temp, col->flows[0]);
+ ASSERT(rc >= 0, "no memory");
+ rc = expand_z(temp, col->values[0], col->flows[0]);
+ ASSERT(rc >= 0, "no memory");
+ darrayDestroy(temp);
+ for (z = 0; z < col->k; z++)
+ memcpy(darrayGetPr(col->image) + z*col->w*col->h*col->comp + c*col->w*col->h, darrayGetPr(col->values[0]) + z*col->w*col->h, col->w*col->h*sizeof(double));
+ }
+#endif
+
+ rc = 0;
+
+error:
+ return rc;
+}
+
+void colorize_free(struct colorize *col)
+{
+ int d;
+
+ /* free memory */
+ for (d = 0; d < col->max_depth; d++) {
+ if (col->nb_list)
+ free(col->nb_list[d]);
+ if (col->marks)
+ darrayDestroy(col->marks[d]);
+ if (col->values)
+ darrayDestroy(col->values[d]);
+ if (col->flows_xy)
+ darrayDestroy(col->flows_xy[d]);
+ if (col->flows)
+ darrayDestroy(col->flows[d]);
+ if (col->flows_i)
+ darrayDestroy(col->flows_i[d]);
+ }
+ free(col->nb_list);
+ free(col->marks);
+ free(col->values);
+ free(col->flows_xy);
+ free(col->flows);
+ free(col->flows_i);
+ darrayDestroy(col->init);
+}
+
+int colorize(const darray_t *image, const darray_t *image_mark, const darray_t *flow, const darray_t *flow_i, int inner_iter, int outer_iter, int scalexyz)
+{
+ struct colorize col;
+ int rc = -1;
+
+ memset(&col, 0, sizeof(col));
+ col.image = image;
+ col.image_mark = image_mark;
+ col.flow = flow;
+ col.flow_i = flow_i;
+ col.inner_iter = inner_iter;
+ col.outer_iter = outer_iter;
+ col.scalexyz = scalexyz;
+
+ rc = colorize_prepare(&col);
+ if (rc < 0)
+ goto error;
+
+ rc = colorize_solve(&col, 0);
+ if (rc < 0)
+ goto error;
+
+ rc = 0;
+
+error:
+ colorize_free(&col);
+
+ return rc;
+}
#ifndef ___COLORIZE_H___
#define ___COLORIZE_H___
-darray_t *colorize(const darray_t *image, const darray_t *image_mark, const darray_t *flow_x, const darray_t *flow_y, const darray_t *flow_ix, const darray_t *flow_iy, int inner_iter, int outer_iter, int scalexyz);
+
+struct colorize {
+ /* input stuff */
+ const darray_t *image;
+ const darray_t *image_mark;
+ const darray_t *flow;
+ const darray_t *flow_i;
+ int inner_iter;
+ int outer_iter;
+ int scalexyz;
+
+ /* internal stuff */
+ int w, h, k, comp, max_depth;
+ unsigned char **nb_list;
+ darray_t **values;
+ darray_t **marks;
+ darray_t *init;
+ darray_t **flows;
+ darray_t **flows_i;
+ darray_t **flows_xy;
+};
+
+int colorize_prepare(struct colorize *col);
+int colorize_solve(struct colorize *col, int quick);
+void colorize_free(struct colorize *col);
+
+int colorize(const darray_t *image, const darray_t *image_mark, const darray_t *flow, const darray_t *flow_i, int inner_iter, int outer_iter, int scalexyz);
#endif
}
}
+#ifdef DARRAY_CHECK
+void darrayCheckPr(const darray_t *array, const double *ptr, int num)
+{
+ if (ptr < array->data) {
+ printf("darrayCheckPr FALED: pointer access before array data\n");
+ abort();
+ }
+ if (((unsigned char *)(array->data) + array->size) < ((unsigned char *)ptr + num*sizeof(double))) {
+ printf("darrayCheckPr FALED: pointer access past array data\n");
+ abort();
+ }
+}
+#endif
+
#ifndef _ARRAY_H
#define _ARRAY_H
+/* enable check function for access outside array data
+ * this requres darrayCheckPr to be called before each access
+ * when disabled, the darrayCheckPr call is not compiled at all
+ */
+//#define DARRAY_CHECK
+
// n-dimensional arrays of double type elements
//
void darrayDone(void);
void darrayCheckList(void);
+#ifdef DARRAY_CHECK
+void darrayCheckPr(const darray_t *array, const double *ptr, int num);
+#else
+#define darrayCheckPr(array, ptr, num)
+#endif
+
#endif /* _ARRAY_H */
#include <stdlib.h>
#include <math.h>
#include "darray.h"
-#include "colorize.h"
+#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
/* just process, but do not generate neighbor list and do not solve multigrid */
//#define TEST_NO_RENDER
/* test expanding of X, Y and Z */
//#define TEST_EXPAND_Z
-/* print error message b, if a is false */
-#define ASSERT(a, b) \
-{ if (!(a)) { fprintf(stderr, "%s %s() Error in line %d: %s\n", __FILE__, __func__, __LINE__, (b)); abort(); } }
-
-void ASSERT_FLOW(int w, int h, int k, const darray_t *flow_x, const darray_t *flow_y, const darray_t *flow_ix, const darray_t *flow_iy)
+void ASSERT_FLOW(int w, int h, int k, const darray_t *flow, const darray_t *flow_i)
{
- if (flow_x) {
- ASSERT(darrayGetNumberOfDimensions(flow_x) == 3, "flow_x array must have 3 dimensions");
- ASSERT(w == darrayGetDimensions(flow_x)[0], "x-sizes of flow_x vectors are different");
- ASSERT(h == darrayGetDimensions(flow_x)[1], "y-sizes of flow_x vectors are different");
- ASSERT(k-1 == darrayGetDimensions(flow_x)[2], "z-sizes of flow_x vectors are not k-1");
- }
- if (flow_y) {
- ASSERT(darrayGetNumberOfDimensions(flow_y) == 3, "flow_y array must have 3 dimensions");
- ASSERT(w == darrayGetDimensions(flow_y)[0], "x-sizes of flow_y vectors are different");
- ASSERT(h == darrayGetDimensions(flow_y)[1], "y-sizes of flow_y vectors are different");
- ASSERT(k-1 == darrayGetDimensions(flow_y)[2], "z-sizes of flow_y vectors are not k-1");
+ if (flow) {
+ ASSERT(darrayGetNumberOfDimensions(flow) == 4, "flow array must have 4 dimensions");
+ ASSERT(w == darrayGetDimensions(flow)[0], "x-sizes of flow vectors are different");
+ ASSERT(h == darrayGetDimensions(flow)[1], "y-sizes of flow vectors are different");
+ ASSERT(k-1 == darrayGetDimensions(flow)[2], "z-sizes of flow vectors are not k-1");
+ ASSERT(2 == darrayGetDimensions(flow)[3], "number of flow cubes are not 2");
}
- if (flow_ix) {
- ASSERT(darrayGetNumberOfDimensions(flow_ix) == 3, "flow_ix array must have 3 dimensions");
- ASSERT(w == darrayGetDimensions(flow_ix)[0], "x-sizes of flow_ix vectors are different");
- ASSERT(h == darrayGetDimensions(flow_ix)[1], "y-sizes of flow_ix vectors are different");
- ASSERT(k-1 == darrayGetDimensions(flow_ix)[2], "z-sizes of flow_ix vectors are not k-1");
- }
- if (flow_iy) {
- ASSERT(darrayGetNumberOfDimensions(flow_iy) == 3, "flow_iy array must have 3 dimensions");
- ASSERT(w == darrayGetDimensions(flow_iy)[0], "x-sizes of flow_iy vectors are different");
- ASSERT(h == darrayGetDimensions(flow_iy)[1], "y-sizes of flow_iy vectors are different");
- ASSERT(k-1 == darrayGetDimensions(flow_iy)[2], "z-sizes of flow_iy vectors are not k-1");
+ if (flow_i) {
+ ASSERT(darrayGetNumberOfDimensions(flow_i) == 4, "flow_i array must have 4 dimensions");
+ ASSERT(w == darrayGetDimensions(flow_i)[0], "x-sizes of flow_i vectors are different");
+ ASSERT(h == darrayGetDimensions(flow_i)[1], "y-sizes of flow_i vectors are different");
+ ASSERT(k-1 == darrayGetDimensions(flow_i)[2], "z-sizes of flow_i vectors are not k-1");
+ ASSERT(2 == darrayGetDimensions(flow_i)[3], "number of flow_i cubes are not 2");
}
}
-#define min(a, b) ((a < b) ? a : b)
-#define max(a, b) ((a > b) ? a : b)
-
#define NEIGHBORS_2D 9
#define NEIGHBORS_3D 27
*
* src plane 4 -------> dst plane 2
*/
-static int shrink_z(darray_t *src, darray_t *dst, const darray_t *flow_x, const darray_t *flow_y)
+int shrink_z(darray_t *src, darray_t *dst, const darray_t *flow)
{
int x, y, z, xx, yy, w, h, ks, kd;
int i, pixles;
ASSERT(w == darrayGetDimensions(dst)[0], "x-size of source and destination must be equal");
ASSERT(h == darrayGetDimensions(dst)[1], "y-size of source and destination must be equal");
ASSERT(kd == ((ks+1)>>1), "z-size of destination must be half of the source (+an extra plane if source z-size is odd)");
- ASSERT_FLOW(w, h, ks, flow_x, flow_y, NULL, NULL);
+ ASSERT_FLOW(w, h, ks, flow, NULL);
sum = darrayCreate(3, darrayGetDimensions(dst));
if (!sum)
/* clear dst array with zeros, fill sum array with ones */
pixles = w*h*kd;
dst_ptr = darrayGetPr(dst);
- for (i = 0; i < pixles; i++)
+ for (i = 0; i < pixles; i++) {
+ darrayCheckPr(dst, dst_ptr, 1);
*dst_ptr++ = 0.0;
+ }
sum_ptr = darrayGetPr(sum);
- for (i = 0; i < pixles; i++)
+ for (i = 0; i < pixles; i++) {
+ darrayCheckPr(sum, sum_ptr, 1);
*sum_ptr++ = 1.0;
+ }
/* first copy all even planes */
for (z = 0; z < ks; z++) {
if ((z&1))
continue;
+ darrayCheckPr(src, darrayGetPr(src) + z*w*h, w*h);
+ darrayCheckPr(dst, darrayGetPr(dst) + (z/2)*w*h, w*h);
memcpy(darrayGetPr(dst) + (z/2)*w*h, darrayGetPr(src) + z*w*h, w*h*sizeof(double));
}
src_ptr = darrayGetPr(src) + z*w*h;
dst_ptr = darrayGetPr(dst) + (z/2)*w*h;
sum_ptr = darrayGetPr(sum) + (z/2)*w*h;
- flow_x_ptr = darrayGetPr(flow_x) + (z-1)*w*h;
- flow_y_ptr = darrayGetPr(flow_y) + (z-1)*w*h;
+ flow_x_ptr = darrayGetPr(flow) + (z-1)*w*h;
+ flow_y_ptr = flow_x_ptr + (ks-1)*w*h;
for (y = 0; y < h; y++) {
for (x = 0; x < w; x++) {
+ darrayCheckPr(flow, flow_x_ptr, 1);
xx = x + *flow_x_ptr++;
+ darrayCheckPr(flow, flow_y_ptr, 1);
yy = y + *flow_y_ptr++;
if (xx>=0 && yy>=0 && xx<w && yy<h) {
+ darrayCheckPr(dst, dst_ptr, 1);
+ darrayCheckPr(src, src_ptr + xx + yy*w, 1);
*dst_ptr += src_ptr[xx + yy*w];
+ darrayCheckPr(sum, sum_ptr, 1);
(*sum_ptr) += 1.0;
}
dst_ptr++;
pixles = w*h*kd;
dst_ptr = darrayGetPr(dst);
sum_ptr = darrayGetPr(sum);
- for (i = 0; i < pixles; i++)
+ for (i = 0; i < pixles; i++) {
+ darrayCheckPr(dst, dst_ptr, 1);
+ darrayCheckPr(sum, sum_ptr, 1);
*dst_ptr++ /= *sum_ptr++;
+ }
darrayDestroy(sum);
if (dst_move)
* A = 0.707
* B = 0.293
*/
-static int expand_z(darray_t *src, const darray_t *dst, const darray_t *flow_x, const darray_t *flow_y)
+static int expand_z(darray_t *src, const darray_t *dst, const darray_t *flow)
{
int x, y, z, xx, yy, xxx, yyy, w, h, ks, kd, k_;
double _x, _y;
ASSERT(w == darrayGetDimensions(dst)[0], "x-size of source and destination must be equal");
ASSERT(h == darrayGetDimensions(dst)[1], "y-size of source and destination must be equal");
ASSERT(ks == ((kd+1)>>1), "z-size of source must be half of the destination (+an extra plane if destination z-size is odd)");
- ASSERT_FLOW(w, h, kd, flow_x, flow_y, NULL, NULL);
+ ASSERT_FLOW(w, h, kd, flow, NULL);
sum = darrayCreate(3, darrayGetDimensions(dst));
if (!sum)
k_ = kd&(~1); /* even number of planes */
if (k_) {
/* copy first plane with weight 1.0 */
+ darrayCheckPr(dst, darrayGetPr(dst) + z*w*h, w*h);
+ darrayCheckPr(src, darrayGetPr(src) + (z/2)*w*h, w*h);
memcpy(darrayGetPr(dst) + z*w*h, darrayGetPr(src) + (z/2)*w*h, w*h*sizeof(double));
/* fill sum array with ones */
pixles = w*h;
sum2_ptr = sum1_ptr + w*h;
src1_ptr = darrayGetPr(src) + (z/2)*w*h;
src2_ptr = src1_ptr + w*h;
- flow_x1_ptr = darrayGetPr(flow_x) + z*w*h;
- flow_y1_ptr = darrayGetPr(flow_y) + z*w*h;
+ flow_x1_ptr = darrayGetPr(flow) + z*w*h;
+ flow_y1_ptr = flow_x1_ptr + (kd-1)*w*h;
flow_x2_ptr = flow_x1_ptr + w*h;
flow_y2_ptr = flow_y1_ptr + w*h;
for (y = 0; y < h; y++) {
for (x = 0; x < w; x++) {
/* even plane (weight A) */
+ darrayCheckPr(dst, dst2_ptr + x + y*w, 1);
+ darrayCheckPr(src, src2_ptr + x + y*w, 1);
dst2_ptr[x + y*w] += src2_ptr[x + y*w] * A;
+ darrayCheckPr(sum, sum2_ptr + x + y*w, 1);
sum2_ptr[x + y*w] += A;
/* get position via firt plane */
+ darrayCheckPr(flow, flow_x1_ptr, 1);
_x = *flow_x1_ptr++;
+ darrayCheckPr(flow, flow_y1_ptr, 1);
_y = *flow_y1_ptr++;
xx = x + round(_x);
yy = y + round(_y);
if (xx<0 || yy<0 || xx>=w || yy>=h)
continue;
/* odd plane (weight A)*/
+ darrayCheckPr(dst, dst1_ptr + xx + yy*w, 1);
+ darrayCheckPr(src, src1_ptr + x + y*w, 1);
dst1_ptr[xx + yy*w] += src1_ptr[x + y*w] * A;
+ darrayCheckPr(sum, sum1_ptr + xx + yy*w, 1);
sum1_ptr[xx + yy*w] += A;
/* follow position via second plane */
+ darrayCheckPr(flow, flow_x2_ptr + xx + yy*w, 1);
_x += flow_x2_ptr[xx + yy*w];
+ darrayCheckPr(flow, flow_y2_ptr + xx + yy*w, 1);
_y += flow_y2_ptr[xx + yy*w];
xxx = x + round(_x);
yyy = y + round(_y);
if (xxx<0 || yyy<0 || xxx>=w || yyy>=h)
continue;
/* even plane (weight B) */
+ darrayCheckPr(dst, dst2_ptr + xxx + yyy*w, 1);
+ darrayCheckPr(src, src1_ptr + x + y*w, 1);
dst2_ptr[xxx + yyy*w] += src1_ptr[x + y*w] * B;
+ darrayCheckPr(sum, sum2_ptr + xxx + yyy*w, 1);
sum2_ptr[xxx + yyy*w] += B;
/* odd plane (weigth B) */
+ darrayCheckPr(dst, dst1_ptr + xx + yy*w, 1);
+ darrayCheckPr(src, src2_ptr + xxx + yyy*w, 1);
dst1_ptr[xx + yy*w] += src2_ptr[xxx + yyy*w] * B;
+ darrayCheckPr(sum, sum1_ptr + xx + yy*w, 1);
sum1_ptr[xx + yy*w] += B;
}
}
dst1_ptr = darrayGetPr(dst) + z*w*h;
sum1_ptr = darrayGetPr(sum) + z*w*h;
src1_ptr = darrayGetPr(src) + (z/2)*w*h;
- flow_x1_ptr = darrayGetPr(flow_x) + z*w*h;
- flow_y1_ptr = darrayGetPr(flow_y) + z*w*h;
+ flow_x1_ptr = darrayGetPr(flow) + z*w*h;
+ flow_y1_ptr = flow_x1_ptr + (kd-1)*w*h;
for (y = 0; y < h; y++) {
for (x = 0; x < w; x++) {
/* get position via firt plane */
+ darrayCheckPr(flow, flow_x1_ptr, 1);
_x = *flow_x1_ptr++;
+ darrayCheckPr(flow, flow_y1_ptr, 1);
_y = *flow_y1_ptr++;
xx = x + round(_x);
yy = y + round(_y);
if (xx<0 || yy<0 || xx>=w || yy>=h)
continue;
/* border (odd) plane (weight 1.0)*/
+ darrayCheckPr(dst, dst1_ptr + xx + yy*w, 1);
+ darrayCheckPr(src, src1_ptr + x + y*w, 1);
dst1_ptr[xx + yy*w] += src1_ptr[x + y*w];
+ darrayCheckPr(sum, sum1_ptr + xx + yy*w, 1);
sum1_ptr[xx + yy*w] += 1.0;
}
}
/* copy extra plane, if destiantion size is odd */
if (z < kd) {
+ darrayCheckPr(dst, darrayGetPr(dst) + z*w*h, w*h);
+ darrayCheckPr(src, darrayGetPr(src) + (z/2)*w*h, w*h);
memcpy(darrayGetPr(dst) + z*w*h, darrayGetPr(src) + (z/2)*w*h, w*h*sizeof(double));
/* fill sum array with ones */
pixles = w*h;
sum1_ptr = darrayGetPr(sum);
for (i = 0; i < pixles; i++) {
/* only if there is a sum */
+ darrayCheckPr(sum, sum1_ptr, 1);
+ darrayCheckPr(dst, dst1_ptr, 1);
if (*sum1_ptr > 0.1)
*dst1_ptr++ /= *sum1_ptr;
else
*
* If the number of flow planes is odd, the extra plane is removed.
*/
-static int shrink_flow_z(darray_t *src_x, darray_t *src_y, darray_t *dst_x, darray_t *dst_y)
+int shrink_flow_z(darray_t *src, darray_t *dst)
{
int w, h, ks, kd, x, y, z, xx, yy;
double _x, _y;
double *src_x1_ptr, *src_x2_ptr, *src_y1_ptr, *src_y2_ptr, *dst_x_ptr, *dst_y_ptr;
- int dst_x_move = 0, dst_y_move = 0;
+ int dst_move = 0;
- ASSERT(darrayGetNumberOfDimensions(src_x) == 3, "src_x array must have 3 dimensions");
- ASSERT(darrayGetNumberOfDimensions(src_y) == 3, "src_y array must have 3 dimensions");
- w = darrayGetDimensions(src_x)[0];
- h = darrayGetDimensions(src_x)[1];
- ks = darrayGetDimensions(src_x)[2];
+ ASSERT(darrayGetNumberOfDimensions(src) == 4, "src array must have 4 dimensions");
+ w = darrayGetDimensions(src)[0];
+ h = darrayGetDimensions(src)[1];
+ ks = darrayGetDimensions(src)[2];
+ ASSERT(2 == darrayGetDimensions(src)[3], "number of cubes must be 2");
kd = ks>>1;
- ASSERT(w == darrayGetDimensions(src_y)[0], "x-size must be equal");
- ASSERT(h == darrayGetDimensions(src_y)[1], "y-size must be equal");
- ASSERT(ks == darrayGetDimensions(src_y)[2], "z-size must be equal");
- if (dst_x == NULL) {
- int dimensions[3] = {w, h, kd};
+ if (dst == NULL) {
+ int dimensions[4] = {w, h, kd, 2};
- dst_x = darrayCreate(3, dimensions);
- if (!dst_x)
- goto error;
- dst_x_move = 1;
- }
- if (dst_y == NULL) {
- int dimensions[3] = {w, h, kd};
-
- dst_y = darrayCreate(3, dimensions);
- if (!dst_y)
+ dst = darrayCreate(4, dimensions);
+ if (!dst)
goto error;
- dst_y_move = 1;
+ dst_move = 1;
}
- ASSERT(darrayGetNumberOfDimensions(dst_x) == 3, "flow_x array must have 3 dimensions");
- ASSERT(darrayGetNumberOfDimensions(dst_y) == 3, "flow_y array must have 3 dimensions");
- ASSERT(w == darrayGetDimensions(dst_x)[0], "x-size must be equal");
- ASSERT(h == darrayGetDimensions(dst_x)[1], "y-size must be equal");
- ASSERT(kd == darrayGetDimensions(dst_x)[2], "z-size of dest must be half of the planes + optional odd plane");
- ASSERT(w == darrayGetDimensions(dst_y)[0], "x-size must be equal");
- ASSERT(h == darrayGetDimensions(dst_y)[1], "y-size must be equal");
- ASSERT(kd == darrayGetDimensions(dst_y)[2], "z-size of dest must be half of the planes + optional odd plane");
+ ASSERT(darrayGetNumberOfDimensions(dst) == 4, "dst array must have 4 dimensions");
+ ASSERT(w == darrayGetDimensions(dst)[0], "x-size must be equal");
+ ASSERT(h == darrayGetDimensions(dst)[1], "y-size must be equal");
+ ASSERT(kd == darrayGetDimensions(dst)[2], "z-size of dest must be half of the planes + optional odd plane");
+ ASSERT(2 == darrayGetDimensions(dst)[3], "number of cubes must be 2");
for (z = 0; z < kd; z++) {
- src_x1_ptr = darrayGetPr(src_x) + (z*2)*w*h;
- src_y1_ptr = darrayGetPr(src_y) + (z*2)*w*h;
- dst_x_ptr = darrayGetPr(dst_x) + z*w*h;
- dst_y_ptr = darrayGetPr(dst_y) + z*w*h;
- src_x2_ptr = darrayGetPr(src_x) + (z*2+1)*w*h;
- src_y2_ptr = darrayGetPr(src_y) + (z*2+1)*w*h;
+ src_x1_ptr = darrayGetPr(src) + (z*2)*w*h;
+ src_y1_ptr = src_x1_ptr + w*h*ks;
+ src_x2_ptr = src_x1_ptr + w*h;
+ src_y2_ptr = src_y1_ptr + w*h;
+ dst_x_ptr = darrayGetPr(dst) + z*w*h;
+ dst_y_ptr = dst_x_ptr + w*h*kd;
for (y = 0; y < h; y++) {
for (x = 0; x < w; x++) {
/* get position via first plane */
+ darrayCheckPr(src, src_x1_ptr, 1);
_x = *src_x1_ptr++;
+ darrayCheckPr(src, src_y1_ptr, 1);
_y = *src_y1_ptr++;
xx = x + round(_x);
yy = y + round(_y);
if (xx<0 || yy<0 || xx>=w || yy>=h) {
/* use an outside value when outside in the first plane */
+ darrayCheckPr(dst, dst_x_ptr, 1);
*dst_x_ptr++ = _x;
+ darrayCheckPr(dst, dst_y_ptr, 1);
*dst_y_ptr++ = _y;
continue;
}
/* follow position via second plane */
+ darrayCheckPr(src, src_x2_ptr + xx + yy*w, 1);
_x += src_x2_ptr[xx + yy*w];
+ darrayCheckPr(src, src_y2_ptr + xx + yy*w, 1);
_y += src_y2_ptr[xx + yy*w];
/* store displacement result */
+ darrayCheckPr(dst, dst_x_ptr, 1);
*dst_x_ptr++ = _x;
+ darrayCheckPr(dst, dst_y_ptr, 1);
*dst_y_ptr++ = _y;
}
}
}
- if (dst_x_move)
- darrayMove(dst_x, src_x);
- if (dst_y_move)
- darrayMove(dst_y, src_y);
+ if (dst_move)
+ darrayMove(dst, src);
return 0;
error:
- if (dst_x_move)
- darrayDestroy(dst_x);
- if (dst_y_move)
- darrayDestroy(dst_y);
+ if (dst_move)
+ darrayDestroy(dst);
return -1;
}
* If a dimension of a source matrix is odd, the end of the dimension is just
* copied and not scaled.
*
- * The points in a plane are shrinked from 2x2x2 to 1x1x1 by calculating the mean
+ * The points in a plane are shrinked from 2x2 to 1x1 by calculating the mean
* of all eight values.
+ *
+ * If the array has more than two dimensions, all planes in all dimensions
+ * are processed.
*/
-static int shrink_xy(darray_t *src, darray_t *dst)
+int shrink_xy(darray_t *src, darray_t *dst)
{
- int x, xx, y, z, ws, wd, w_, hs, hd, h_, k;
+ int x, xx, y, z, ws, wd, w_, hs, hd, h_, k, dims, i;
double *s1, *s2, *d;
int dst_move = 0;
-
- ASSERT(darrayGetNumberOfDimensions(src) == 3, "Source array must have 3 dimensions");
+
+ dims = darrayGetNumberOfDimensions(src);
+ ASSERT(dims >= 2, "Source array must have at least 2 dimensions");
ws = darrayGetDimensions(src)[0];
w_ = ws>>1;
hs = darrayGetDimensions(src)[1];
h_ = hs>>1;
- k = darrayGetDimensions(src)[2];
+ k = 1;
+ for (i = 2; i < dims; i++)
+ k *= darrayGetDimensions(src)[i];
/* create dst, if not given */
if (!dst) {
- int dimensions[3] = {(ws+1)>>1, (hs+1)>>1, k};
-
- dst = darrayCreate(3, dimensions);
+ int dimensions[dims];
+
+ dimensions[0] = (ws+1)>>1;
+ dimensions[1] = (hs+1)>>1;
+ for (i = 2; i < dims; i++)
+ dimensions[i] = darrayGetDimensions(src)[i];
+
+ dst = darrayCreate(dims, dimensions);
if (!dst)
goto error;
dst_move = 1;
}
- ASSERT(darrayGetNumberOfDimensions(dst) == 3, "Destination array must have 3 dimensions");
+ ASSERT(darrayGetNumberOfDimensions(dst) == dims, "Destination array must have equal dimensions than source array");
wd = darrayGetDimensions(dst)[0];
hd = darrayGetDimensions(dst)[1];
ASSERT((ws&1) || ws == wd*2, "x-size of source must be twice the size of destination (this is even source size)");
ASSERT(!(ws&1) || ws-1 == (wd-1)*2, "x-size of source-1 must be twice the size of destination-1 (this odd source size)");
ASSERT((hs&1) || hs == hd*2, "y-size of source must be twice the size of destination (this is even source size)");
ASSERT(!(hs&1) || hs-1 == (hd-1)*2, "y-size of source-1 must be twice the size of destination-1 (this odd source size)");
- ASSERT(k == darrayGetDimensions(dst)[2], "z-size must be equal");
+ for (i = 2; i < dims; i++)
+ ASSERT(darrayGetDimensions(src)[i] == darrayGetDimensions(dst)[i], "dimension 3 and greater must be equal");
for (z = 0; z < k; z++) {
for (y = 0; y < h_; y++) {
s1 = darrayGetPr(src) + y*2*ws + z*ws*hs;
s2 = s1 + ws;
d = darrayGetPr(dst) + y*wd + z*wd*hd;
- for (x = 0, xx = 0; x < w_; x++, xx+=2)
+ for (x = 0, xx = 0; x < w_; x++, xx+=2) {
+ darrayCheckPr(dst, d + x, 1);
+ darrayCheckPr(src, s1 + xx, 2);
+ darrayCheckPr(src, s2 + xx, 2);
d[x] = (s1[xx] + s1[xx+1] + s2[xx] + s2[xx+1]) * 0.25;
+ }
/* mean of 1x2x1 */
- if (x < wd)
+ if (x < wd) {
+ darrayCheckPr(dst, d + x, 1);
+ darrayCheckPr(src, s1 + xx, 1);
+ darrayCheckPr(src, s2 + xx, 1);
d[x] = (s1[xx] + s2[xx]) * 0.5;
+ }
}
if (y < hd) {
/* mean of 2x1x1 */
s1 = darrayGetPr(src) + y*2*ws + z*ws*hs;
d = darrayGetPr(dst) + y*wd + z*wd*hd;
- for (x = 0, xx = 0; x < w_; x++, xx+=2)
+ for (x = 0, xx = 0; x < w_; x++, xx+=2) {
+ darrayCheckPr(dst, d + x, 1);
+ darrayCheckPr(src, s1 + xx, 2);
d[x] = (s1[xx] + s1[xx+1]) * 0.5;
+ }
/* mean of 1x1x1 */
- if (x < wd)
+ if (x < wd) {
+ darrayCheckPr(dst, d + x, 1);
+ darrayCheckPr(src, s1 + xx, 1);
d[x] = s1[xx];
+ }
}
}
d = darrayGetPr(dst) + (y+h_)*wd + z*wd*hd;
*d++ = *s; /* x=0 plane */
for (x = 0; x < w_-1; x++) {
+ darrayCheckPr(dst, d, 1);
+ darrayCheckPr(src, s, 2);
*d++ = A*s[0] + B*s[1];
*d++ = A*s[1] + B*s[0];
s++;
}
+ darrayCheckPr(dst, d, 1);
+ darrayCheckPr(src, s, 1);
*d++ = *s++; /* x=w-1 plane */
/* extra plane (odd fine grid) */
- if ((wd&1))
+ if ((wd&1)) {
+ darrayCheckPr(dst, d, 1);
+ darrayCheckPr(src, s, 1);
*d = *s;
+ }
}
}
/* y=h-1 plane */
s = darrayGetPr(dst) + h_*wd + z*wd*hd;
d = darrayGetPr(dst) + z*wd*hd;
- for (x = 0; x < wd; x++)
+ for (x = 0; x < wd; x++) {
+ darrayCheckPr(dst, d, 1);
+ darrayCheckPr(dst, s, 1);
*d++ = *s++; /* y=0 plane */
+ }
for (y = 0; y < h_-1; y++) {
s1 = darrayGetPr(dst) + (y+h_)*wd + z*wd*hd;
s2 = s1 + wd;
d1 = darrayGetPr(dst) + (y*2+1)*wd + z*wd*hd;
d2 = d1 + wd;
for (x = 0; x < wd; x++) {
+ darrayCheckPr(dst, d1, 1);
+ darrayCheckPr(dst, d2, 1);
+ darrayCheckPr(dst, s1, 1);
+ darrayCheckPr(dst, s2, 1);
*d1++ = A * *s1 + B * *s2;
*d2++ = A * *s2 + B * *s1;
s1++;
*
* If an element is greater theshold, the element is set to gt, otherwise le.
*/
-static void quantize(const darray_t *array, double threshold, double gt, double le)
+void quantize(const darray_t *array, double threshold, double gt, double le)
{
- int n = 1, i;
+ int n = 1, dims, i;
double *ptr;
- ASSERT(darrayGetNumberOfDimensions(array) == 3, "array must have 3 dimensions");
- n *= darrayGetDimensions(array)[0];
- n *= darrayGetDimensions(array)[1];
- n *= darrayGetDimensions(array)[2];
+ dims = darrayGetNumberOfDimensions(array);
+ for (i = 0; i < dims; i++)
+ n *= darrayGetDimensions(array)[i];
ptr = darrayGetPr(array);
for (i = 0; i < n; i++, ptr++) {
+ darrayCheckPr(array, ptr, 1);
if (*ptr > threshold)
*ptr = gt;
else
/* multiply the array content with a value
*/
-static void multiply(const darray_t *array, double factor)
+void multiply(const darray_t *array, double factor)
{
- int n = 1, i;
+ int n = 1, dims, i;
double *ptr;
- ASSERT(darrayGetNumberOfDimensions(array) == 3, "array must have 3 dimensions");
- n *= darrayGetDimensions(array)[0];
- n *= darrayGetDimensions(array)[1];
- n *= darrayGetDimensions(array)[2];
+ dims = darrayGetNumberOfDimensions(array);
+ for (i = 0; i < dims; i++)
+ n *= darrayGetDimensions(array)[i];
ptr = darrayGetPr(array);
- for (i = 0; i < n; i++, ptr++)
+ for (i = 0; i < n; i++, ptr++) {
+ darrayCheckPr(array, ptr, 1);
*ptr *= factor;
+ }
}
/* create an array of neighbors and use flow vectors
* 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
*/
-static unsigned char *gen_neighbor(int w, int h, int k, const darray_t *flow_x, const darray_t *flow_y, const darray_t *flow_ix, const darray_t *flow_iy)
+unsigned char *gen_neighbor(int w, int h, int k, const darray_t *flow, const darray_t *flow_i)
{
int x, y, z, xx, yy, zz, i, j, n;
double _x, _y;
int nb_num;
int *i_index;
- ASSERT_FLOW(w, h, k, flow_x, flow_y, flow_ix, flow_iy);
- flow_x_ptr = darrayGetPr(flow_x);
- flow_y_ptr = darrayGetPr(flow_y);
- flow_ix_ptr = darrayGetPr(flow_ix);
- flow_iy_ptr = darrayGetPr(flow_iy);
+ ASSERT_FLOW(w, h, k, flow, flow_i);
+ flow_x_ptr = darrayGetPr(flow);
+ flow_y_ptr = flow_x_ptr + w*h*(k-1);
+ flow_ix_ptr = darrayGetPr(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);
if (xx >= w || yy >= h || zz >= k)
continue;
if (zz > z) {
+ darrayCheckPr(flow, flow_x_ptr + xx + yy*w + (zz-1)*w*h, 1);
_x = flow_x_ptr[xx + yy*w + (zz-1)*w*h];
+ darrayCheckPr(flow, flow_y_ptr + xx + yy*w + (zz-1)*w*h, 1);
_y = flow_y_ptr[xx + yy*w + (zz-1)*w*h];
xx = xx + round(_x);
yy = yy + round(_y);
} else if (zz < z) {
+ darrayCheckPr(flow_i, flow_ix_ptr + xx + yy*w + zz*w*h, 1);
_x = flow_ix_ptr[xx + yy*w + zz*w*h];
+ darrayCheckPr(flow_i, flow_iy_ptr + xx + yy*w + zz*w*h, 1);
_y = flow_iy_ptr[xx + yy*w + zz*w*h];
xx = xx + round(_x);
yy = yy + round(_y);
}
/* weighten neighbors */
-static void weighten(const darray_t *luminance, unsigned char *neighbors)
+void weighten(const darray_t *luminance, unsigned char *neighbors)
{
int x, y, z, i, w, h, k, num;
double val, nval, mean, variance, sigma, sum, weight, min;
variance = 0;
for (i = 0; i < num; i++) {
/* value of neighbor pixles, including current pixle */
+ darrayCheckPr(luminance, lum_ptr + i_index[i], 1);
nval = lum_ptr[i_index[i]];
/* squared diffs between current pixle and each neighbor */
sqdiff[i] = (nval-val) * (nval-val);
w_index = ((struct neighbors_3d *)neighbors)->weight;
}
/* if pixle is marked */
+ darrayCheckPr(mark, mark_line_ptr + x, 1);
if (mark_line_ptr[x] < 0.5) {
/* always use initial color on marked plane, otherwise do not smooth */
- if (level == 0)
+ if (level == 0) {
+ darrayCheckPr(value, value_line_ptr + x, 1);
+ darrayCheckPr(init, init_line_ptr + x, 1);
value_line_ptr[x] = init_line_ptr[x];
+ }
continue;
}
/* Calculate value */
val = 0;
- for (i = 1; i < num; i++) /* sum value without point itself */
+ for (i = 1; i < num; i++) { /* sum value without point itself */
+ darrayCheckPr(value, value_ptr + i_index[i], 1);
val += w_index[i] * value_ptr[i_index[i]];
- if (level == 0)
+ }
+ if (level == 0) {
+ darrayCheckPr(value, value_line_ptr + x, 1);
+ darrayCheckPr(init, init_line_ptr + x, 1);
value_line_ptr[x] = (init_line_ptr[x] - val) / w_index[0];
- else
+ } else {
+ darrayCheckPr(value, value_line_ptr + x, 1);
value_line_ptr[x] = (-val) / w_index[0];
+ }
}
}
}
}
+#ifdef RESIDUAL
/* Calculate residual
*/
static double residual(int level, darray_t *value, const darray_t *mark, const darray_t *init, unsigned char *neighbors)
w_index = ((struct neighbors_3d *)neighbors)->weight;
}
/* if pixle is marked, no residual */
+ darrayCheckPr(mark, mark_line_ptr + x, 1);
if (mark_line_ptr[x] < 0.5)
continue;
/* Calculate value */
val = 0;
- for (i = 0; i < num; i++) /* sum value including point itself */
+ for (i = 0; i < num; i++) { /* sum value including point itself */
+ darrayCheckPr(value, value_ptr + i_index[i], 1);
val += w_index[i] * value_ptr[i_index[i]];
- if (level == 0)
+ }
+ if (level == 0) {
+ darrayCheckPr(init, init_line_ptr + x, 1);
val = init_line_ptr[x] - val;
- else
+ } else
val = -val;
residual += val * val;
}
return residual;
}
+#endif
/* Restrict
*
* different depth (z dimension), a temporary array (temp) is used to
* perform the intermediate step of scaling x,y,z.
*
- * Note that the flow vectors must be scaled xy already!
+ * Note that the flow vectors must be scaled xy already (if x and y is shrinked
+ * before z)!
*/
-static int restrict(const darray_t *value, darray_t *nvalue, const darray_t *flow_x, const darray_t *flow_y)
+static int restrict(const darray_t *value, darray_t *nvalue, const darray_t *flow)
{
int ks, wd, hd, kd;
darray_t *temp = NULL;
kd = darrayGetDimensions(nvalue)[2];
ks = darrayGetDimensions(value)[2];
if (ks != kd)
- ASSERT(flow_x && flow_y, "when scaling z, flow vectors must be given");
+ ASSERT(flow, "when scaling z, flow vectors must be given");
if (ks != kd) {
int dimensions[3] = {wd, hd, ks};
rc = shrink_xy((darray_t *)value, temp);
if (rc < 0)
goto error;
- rc = shrink_z(temp, nvalue, flow_x, flow_y);
+ rc = shrink_z(temp, nvalue, flow);
if (rc < 0)
goto error;
} else {
* different depth (z dimension), a temporary array (temp3) is used to
* perform the intermediate step of scaling x,y,z.
*
- * Note that the flow vectors must be scaled xy already!
+ * Note that the flow vectors must be scaled xy already (if x and y is shrinked
+ * before z)!
*
* First the value is shrinked to the size of the value of the next level.
* The error is calculated by substracting the values of the next level
* substracted from the value. So the changes in the next level's value are
* prolonged to the value.
*/
-static int prolong(darray_t *value, const darray_t *nvalue, const darray_t *flow_x, const darray_t *flow_y)
+static int prolong(darray_t *value, const darray_t *nvalue, const darray_t *flow)
{
int ws, hs, ks, wd, hd, kd;
int i, pixles;
hd = darrayGetDimensions(value)[1];
kd = darrayGetDimensions(value)[2];
if (ks != kd)
- ASSERT(flow_x && flow_y, "when scaling z, flow vectors must be given");
+ ASSERT(flow, "when scaling z, flow vectors must be given");
/* generate temporary arrays */
temp1 = darrayCreate(3, darrayGetDimensions(value));
rc = shrink_xy(value, temp3);
if (rc < 0)
goto error;
- rc = shrink_z(temp3, temp2, flow_x, flow_y);
+ rc = shrink_z(temp3, temp2, flow);
if (rc < 0)
goto error;
} else {
temp2_ptr = darrayGetPr(temp2);
nvalue_ptr = darrayGetPr(nvalue);
pixles = ws*hs*ks;
- for (i = 0; i < pixles; i++)
+ for (i = 0; i < pixles; i++) {
+ darrayCheckPr(temp2, temp2_ptr, 1);
+ darrayCheckPr(nvalue, nvalue_ptr, 1);
*temp2_ptr++ -= *nvalue_ptr++;
+ }
/* Expand (temp2 -> temp1) */
if (ks != kd) {
- rc = expand_z(temp2, temp3, flow_x, flow_y);
+ rc = expand_z(temp2, temp3, flow);
if (rc < 0)
goto error;
expand_xy(temp3, temp1);
temp1_ptr = darrayGetPr(temp1);
value_ptr = darrayGetPr(value);
pixles = wd*hd*kd;
- for (i = 0; i < pixles; i++)
+ for (i = 0; i < pixles; i++) {
+ darrayCheckPr(temp1, temp1_ptr, 1);
+ darrayCheckPr(value, value_ptr, 1);
*value_ptr++ -= *temp1_ptr++;
+ }
rc = 0;
for (z = 0; z < ks; z++) {
for (y = 0; y < hs; y++) {
for (x = 0; x < ws; x++) {
+ darrayCheckPr(dst, dst_ptr + (x+xoff) + (y+yoff)*wd + (z+zoff)*wd*hd, 1);
+ darrayCheckPr(src, src_ptr + x + y*ws + z*ws*hs, 1);
dst_ptr[(x+xoff) + (y+yoff)*wd + (z+zoff)*wd*hd] = src_ptr[x + y*ws + z*ws*hs];
}
}
* The flow vectors are drawn to the upper left corner of the dst, which
* may be larger in size than the flow arrays.
*/
-void draw_flow(darray_t *dst, const darray_t *flow_x, const darray_t *flow_y, const darray_t *flow_ix, const darray_t *flow_iy)
+void draw_flow(darray_t *dst, const darray_t *flow, const darray_t *flow_i)
{
int wd, hd, ws, hs, ks, x, y, z;
double *flow_x_ptr, *flow_y_ptr, *flow_ix_ptr, *flow_iy_ptr;
wd = darrayGetDimensions(dst)[0];
hd = darrayGetDimensions(dst)[1];
- ws = darrayGetDimensions(flow_x)[0];
- hs = darrayGetDimensions(flow_x)[1];
- ks = darrayGetDimensions(flow_x)[2];
- flow_x_ptr = darrayGetPr(flow_x);
- flow_y_ptr = darrayGetPr(flow_y);
- flow_ix_ptr = darrayGetPr(flow_ix);
- flow_iy_ptr = darrayGetPr(flow_iy);
+ ws = darrayGetDimensions(flow)[0];
+ hs = darrayGetDimensions(flow)[1];
+ ks = darrayGetDimensions(flow)[2];
+ flow_x_ptr = darrayGetPr(flow);
+ flow_y_ptr = flow_x_ptr + ws*hs*ks;
+ flow_ix_ptr = darrayGetPr(flow_i);
+ flow_iy_ptr = flow_ix_ptr + ws*hs*ks;
for (z = 0; z <= ks; z++) {
for (y = step/2; y < hs-step/2+1; y+=step) {
/* Solve multigrid
*
- * Note that the flow vectors must be scaled xy already!
+ * Note that the flow vectors must be scaled xy already (if x and y is shrinked
+ * before z)!
+ *
*/
-static 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_xy_x, darray_t **flows_xy_y, 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)
{
+#ifdef RESIDUAL
double curr_residual, prev_residual = 10000000000000.0;
+#endif
int cycle, level, i;
int rc = -1;
for (cycle = 0; cycle < ncycle; cycle++) {
/* go down */
+#ifdef RESIDUAL
curr_residual = 0;
+#endif
for (level = 0; level < nlevel-1; level++) {
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]);
+#endif
if (level < scalexyz)
- rc = restrict(values[level], values[level+1], flows_xy_x[level], flows_xy_y[level]);
+ rc = restrict(values[level], values[level+1], flows[level]);
else
- rc = restrict(values[level], values[level+1], NULL, NULL);
+ rc = restrict(values[level], values[level+1], NULL);
if (rc < 0)
goto error;
}
/* go up */
for (level = nlevel-2; level >= 0; level--) {
if (level < scalexyz)
- rc = prolong(values[level], values[level+1], flows_xy_x[level], flows_xy_y[level]);
+ rc = prolong(values[level], values[level+1], flows[level]);
else
- rc = prolong(values[level], values[level+1], NULL, NULL);
+ rc = prolong(values[level], values[level+1], NULL);
if (rc < 0)
goto error;
for (i = 0; i < iters; i++)
}
/* end if residual raises again */
-// if (prev_residual < curr_residual)
-// break;
+#ifdef RESIDUAL
+ if (prev_residual < curr_residual)
+ break;
curr_residual = prev_residual;
+#endif
}
return 0;
return rc;
}
-/* colorize image array using image_mark
- *
- * image:
- * first dimansion = width
- * second dimension = height
- * thrid dimension = components (first component is luminance, other components are features, like u, v)
- * fourth dimension = time
- *
- * image_mark:
- * first dimansion = width
- * second dimension = height
- * fourth dimension = time
- */
-darray_t *colorize(const darray_t *image, const darray_t *image_mark, const darray_t *flow_x, const darray_t *flow_y, const darray_t *flow_ix, const darray_t *flow_iy, int inner_iter, int outer_iter, int scalexyz)
-{
- unsigned char **nb_list = NULL;
- darray_t **values = NULL, **marks = NULL, *init = NULL, **flows_x = NULL, **flows_y = NULL, **flows_ix = NULL, **flows_iy = NULL, **flows_xy_x = NULL, **flows_xy_y = NULL;
- darray_t *result = NULL, *luminance = NULL;
- int dimensions[4];
- int w, h, k, ww, hh, kk, last_kk, comp, d, max_depth, temp, z, c;
- int rc = -1;
-
- ASSERT(darrayGetNumberOfDimensions(image) == 4, "image array must have 4 dimensions");
- w = darrayGetDimensions(image)[0];
- h = darrayGetDimensions(image)[1];
- comp = darrayGetDimensions(image)[2];
- k = darrayGetDimensions(image)[3];
- ASSERT(comp > 0, "image array requires at least one component (Y-component)");
-
- ASSERT(darrayGetNumberOfDimensions(image_mark) == 3, "image_mark array must have 3 dimensions");
- ASSERT(w == darrayGetDimensions(image_mark)[0], "image_mark and image must both have same x-size");
- ASSERT(h == darrayGetDimensions(image_mark)[1], "image_mark and image must both have same y-size");
- ASSERT(k == darrayGetDimensions(image_mark)[2], "image_mark and image must both have same z-size");
-
- /* Calculate maximum depth
- * Scale down dimensions until at least 2x2 remains
- */
- temp = min(w, h);
- for (max_depth = 0; temp > 2; max_depth++)
- temp=(temp+1)>>1;
- ASSERT(max_depth >= 2, "image too small");
-
- nb_list = calloc(max_depth, sizeof(struct neighbors *));
- if (!nb_list)
- goto error;
- marks = calloc(max_depth, sizeof(darray_t *));
- if (!marks)
- goto error;
- ASSERT(marks, "no memory");
- values = calloc(max_depth, sizeof(darray_t *));
- if (!values)
- goto error;
- ASSERT(values, "no memory");
- flows_x = calloc(max_depth, sizeof(darray_t *));
- if (!flows_x)
- goto error;
- ASSERT(flows_x, "no memory");
- flows_y = calloc(max_depth, sizeof(darray_t *));
- if (!flows_y)
- goto error;
- ASSERT(flows_y, "no memory");
- flows_ix = calloc(max_depth, sizeof(darray_t *));
- if (!flows_ix)
- goto error;
- ASSERT(flows_ix, "no memory");
- flows_iy = calloc(max_depth, sizeof(darray_t *));
- if (!flows_iy)
- goto error;
- ASSERT(flows_iy, "no memory");
- if (scalexyz) {
- /* the flows_xy are scaled in xy in the current level */
- flows_xy_x = calloc(max_depth, sizeof(darray_t *));
- if (!flows_xy_x)
- goto error;
- flows_xy_y = calloc(max_depth, sizeof(darray_t *));
- if (!flows_xy_y)
- goto error;
- }
- dimensions[0] = w;
- dimensions[1] = h;
- dimensions[2] = k;
- luminance = darrayCreate(3, dimensions);
- if (!luminance)
- goto error;
- /* generate sparse multigrid */
- last_kk = kk = k;
- printf("-> prepare sparse grid");
- for (d = 0, ww=w, hh=h; d < max_depth; d++, ww=(ww+1)>>1, hh=(hh+1)>>1) {
- /* create temporary arrays
- * copy fist array to level 0 or shrink previous level
- */
- printf(" #%d(w=%d h=%d k=%d)", d, ww, hh, kk); fflush(stdout);
- dimensions[0] = ww;
- dimensions[1] = hh;
- dimensions[2] = last_kk;
- marks[d] = darrayCreate(3, dimensions);
- if (!marks[d])
- goto error;
- dimensions[0] = ww;
- dimensions[1] = hh;
- dimensions[2] = last_kk-1;
- flows_x[d] = darrayCreate(3, dimensions);
- if (!flows_x[d])
- goto error;
- flows_y[d] = darrayCreate(3, dimensions);
- if (!flows_y[d])
- goto error;
- flows_ix[d] = darrayCreate(3, dimensions);
- if (!flows_ix[d])
- goto error;
- flows_iy[d] = darrayCreate(3, dimensions);
- if (!flows_iy[d])
- goto error;
- if (d < scalexyz) {
- dimensions[0] = (ww+1)>>1;
- dimensions[1] = (hh+1)>>1;
- dimensions[2] = kk-1;
- flows_xy_x[d] = darrayCreate(3, dimensions);
- if (!flows_xy_x[d])
- goto error;
- flows_xy_y[d] = darrayCreate(3, dimensions);
- if (!flows_xy_y[d])
- goto error;
- }
- if (d == 0) {
- for (z = 0; z < kk; z++)
- memcpy(darrayGetPr(luminance) + z*w*h, darrayGetPr(image) + z*w*h*comp, w*h*sizeof(double));
- if (flow_x && flow_y && flow_ix && flow_iy) {
- memcpy(darrayGetPr(flows_x[d]), darrayGetPr(flow_x), w*h*(kk-1)*sizeof(double));
- memcpy(darrayGetPr(flows_y[d]), darrayGetPr(flow_y), w*h*(kk-1)*sizeof(double));
- memcpy(darrayGetPr(flows_ix[d]), darrayGetPr(flow_ix), w*h*(kk-1)*sizeof(double));
- memcpy(darrayGetPr(flows_iy[d]), darrayGetPr(flow_iy), w*h*(kk-1)*sizeof(double));
- }
- memcpy(darrayGetPr(marks[d]), darrayGetPr(image_mark), w*h*kk*sizeof(double));
- /* invert, so 0 = marked */
- quantize(marks[d], 0.5, 0, 1);
- } else {
- rc = shrink_xy(luminance, NULL);
- if (rc < 0)
- goto error;
- rc = shrink_xy(marks[d-1], marks[d]);
- if (rc < 0)
- goto error;
- rc = shrink_xy(flows_x[d-1], flows_x[d]);
- if (rc < 0)
- goto error;
- rc = shrink_xy(flows_y[d-1], flows_y[d]);
- if (rc < 0)
- goto error;
- rc = shrink_xy(flows_ix[d-1], flows_ix[d]);
- if (rc < 0)
- goto error;
- rc = shrink_xy(flows_iy[d-1], flows_iy[d]);
- if (rc < 0)
- goto error;
- /* shrink length of flow vector on smaller grid.
- * This is not performed when comparing with
- * the original algorithm. (it is a bug there)
- */
- multiply(flows_x[d], 0.5);
- multiply(flows_y[d], 0.5);
- multiply(flows_ix[d], 0.5);
- multiply(flows_iy[d], 0.5);
- if (d-1 < scalexyz) { /* process if previous d shall be shrinked in z direction */
- /* shrink time */
- rc = shrink_z(luminance, NULL, flows_x[d], flows_y[d]);
- if (rc < 0)
- goto error;
- rc = shrink_z(marks[d], NULL, flows_x[d], flows_y[d]);
- if (rc < 0)
- goto error;
- rc = shrink_flow_z(flows_x[d], flows_y[d], NULL, NULL);
- if (rc < 0)
- goto error;
- rc = shrink_flow_z(flows_ix[d], flows_iy[d], NULL, NULL);
- if (rc < 0)
- goto error;
- }
- /* If any of the restricted pixles is 0 (marked),
- * the result must be 0 (marked) too.
- */
- quantize(marks[d], 0.9999, 1, 0);
- }
- if (d < scalexyz) {
- /* create x and y shrinked flow vectors for restriction and prolonging process */
- rc = shrink_xy(flows_x[d], flows_xy_x[d]);
- if (rc < 0)
- goto error;
- rc = shrink_xy(flows_y[d], flows_xy_y[d]);
- if (rc < 0)
- goto error;
- multiply(flows_xy_x[d], 0.5);
- multiply(flows_xy_y[d], 0.5);
- }
-
- /* generate values array */
- dimensions[0] = ww;
- dimensions[1] = hh;
- dimensions[2] = kk;
- values[d] = darrayCreate(3, dimensions);
- if (!values[d])
- goto error;
-
- /* generate array of neighbors */
-#ifndef TEST_NO_RENDER
- nb_list[d] = gen_neighbor(ww, hh, kk, flows_x[d], flows_y[d], flows_ix[d], flows_iy[d]);
- if (!nb_list[d])
- goto error;
- /* weighten */
- weighten(luminance, nb_list[d]);
-#endif
-
- last_kk = kk;
- if (d < scalexyz)
- kk = (kk+1)>>1;
- }
- darrayDestroy(luminance);
- luminance = NULL;
- printf("\n");
-
- dimensions[0] = w;
- dimensions[1] = h;
- dimensions[2] = k;
- init = darrayCreate(3, dimensions);
- if (!init)
- goto error;
-
- /* result array for all components */
- dimensions[0] = w;
- dimensions[1] = h;
- dimensions[2] = comp;
- dimensions[3] = k;
- result = darrayCreate(4, dimensions);
- if (!result)
- goto error;
-
-#ifdef TEST_PYRAMID
- dimensions[0] = w;
- dimensions[1] = h;
- dimensions[2] = k;
- luminance = darrayCreate(3, dimensions);
- for (z = 0; z < k; z++) {
- memcpy(darrayGetPr(luminance) + z*w*h, darrayGetPr(image) + z*w*h*comp, w*h*sizeof(double));
- memcpy(darrayGetPr(values[0]) + z*w*h, darrayGetPr(image) + z*w*h*comp, w*h*sizeof(double));
-#ifdef TEST_FLOW
- draw_flow(values[0], flows_x[0], flows_y[0], flows_ix[0], flows_iy[0]);
-#endif
-#ifdef TEST_MARKED
- paste(marks[0], values[0], 0, 0, 0);
-#endif
-
- }
- for (d = 1; d < max_depth; d++) {
- rc = shrink_xy(luminance, NULL);
- if (rc < 0)
- goto error;
- if (d-1 < scalexyz)
- shrink_z(luminance, NULL, flows_xy_x[d-1], flows_xy_y[d-1]);
- paste(luminance, values[0], 0, 0, 0);
-#ifdef TEST_FLOW
- draw_flow(values[0], flows_x[d], flows_y[d], flows_ix[d], flows_iy[d]);
-#endif
-#ifdef TEST_MARKED
- paste(marks[d], values[0], 0, 0, 0);
-#endif
- }
- darrayDestroy(luminance);
- luminance = NULL;
- /* copy test pyramid image as Y-component */
- for (z = 0; z < k; z++)
- memcpy(darrayGetPr(result) + z*w*h*comp, darrayGetPr(values[0]) + z*w*h, w*h*sizeof(double));
-#else
- /* copy Y-component */
- for (z = 0; z < k; z++)
- memcpy(darrayGetPr(result) + z*w*h*comp, darrayGetPr(image) + z*w*h*comp, w*h*sizeof(double));
-#endif
-
- printf("-> apply component to grid and solve");
- for (c = 1; c < comp; c++) {
- int pixles, j;
- double *value_ptr, *init_ptr, *mark_ptr;
-
- printf(" #%d", c); fflush(stdout);
- /* apply component
- * - copy component into value array (level 0), otherwise 0
- * - use maked colors for init array, otherwise use 0
- */
- for (z = 0; z < k; z++)
- memcpy(darrayGetPr(values[0]) + z*w*h, darrayGetPr(image) + c*w*h + z*w*h*comp, w*h*sizeof(double));
- for (d = 1, ww=(w+1)>>1, hh=(h+1)>>1, kk=k; d < max_depth; d++, ww=(ww+1)>>1, hh=(hh+1)>>1) {
- if (d-1 < scalexyz) /* scale kk, if last array was scaled */
- kk=(kk+1)>>1;
- value_ptr = darrayGetPr(values[d]);
- pixles = ww*hh*kk;
- for (j = 0; j < pixles; j++)
- *value_ptr++ = 0.0;
- }
- value_ptr = darrayGetPr(values[0]);
- init_ptr = darrayGetPr(init);
- mark_ptr = darrayGetPr(image_mark);
- pixles = w*h*k;
- for (j = 0; j < pixles; j++) {
- if (*mark_ptr++ > 0.5)
- *init_ptr++ = *value_ptr;
- else
- *init_ptr++ = 0.0;
- value_ptr++;
- }
-
-#ifndef TEST_NO_RENDER
- /* colorize component */
- int i;
- for (i = 0; i < outer_iter; i++) {
- rc = solve_mg(max_depth, inner_iter, 2, values, marks, init, nb_list, flows_xy_x, flows_xy_y, scalexyz);
- if (rc < 0)
- goto error;
- rc = solve_mg((max_depth+1) / 2, inner_iter, 2, values, marks, init, nb_list, flows_xy_x, flows_xy_y, scalexyz);
- if (rc < 0)
- goto error;
- rc = solve_mg(2, inner_iter, 2, values, marks, init, nb_list, flows_xy_x, flows_xy_y, scalexyz);
- if (rc < 0)
- goto error;
- rc = solve_mg(1, inner_iter, 4, values, marks, init, nb_list, flows_xy_x, flows_xy_y, scalexyz);
- if (rc < 0)
- goto error;
- }
-#endif
-
-#ifdef TEST_PYRAMID
- for (d = 1; d < max_depth; d++)
- paste(values[d], values[0], 0, 0, 0);
-#endif
-
- /* apply colorized component */
- for (z = 0; z < k; z++)
- memcpy(darrayGetPr(result) + z*w*h*comp + c*w*h, darrayGetPr(values[0]) + z*w*h, w*h*sizeof(double));
- }
- printf("\n");
-
-#ifdef TEST_SHRINK_XY
- for (c = 0; c < comp; c++) {
- darray_t *temp;
- int dimensions[3] = {w, h, k};
- for (z = 0; z < k; z++)
- memcpy(darrayGetPr(values[0]) + z*w*h, darrayGetPr(result) + z*w*h*comp + c*w*h, w*h*sizeof(double));
- dimensions[0] = (w+1)>>1;
- dimensions[1] = (h+1)>>1;
- temp = darrayCreate(3, dimensions);
- ASSERT(temp, "no memory");
- rc = shrink_xy(values[0], temp);
- ASSERT(rc >= 0, "no memory");
- paste(temp, values[0], 0, 0, 0);
- darrayDestroy(temp);
- for (z = 0; z < k; z++)
- memcpy(darrayGetPr(result) + z*w*h*comp + c*w*h, darrayGetPr(values[0]) + z*w*h, w*h*sizeof(double));
- }
-#endif
-
-#ifdef TEST_SHRINK_Z
- for (c = 0; c < comp; c++) {
- darray_t *temp;
- int dimensions[3] = {w, h, k};
- for (z = 0; z < k; z++)
- memcpy(darrayGetPr(values[0]) + z*w*h, darrayGetPr(result) + z*w*h*comp + c*w*h, w*h*sizeof(double));
- dimensions[2] = (k+1)>>1;
- temp = darrayCreate(3, dimensions);
- ASSERT(temp, "no memory");
- rc = shrink_z(values[0], temp, flows_x[0], flows_y[0]);
- ASSERT(rc >= 0, "no memory");
- paste(temp, values[0], 0, 0, 0);
- darrayDestroy(temp);
- for (z = 0; z < k; z++)
- memcpy(darrayGetPr(result) + z*w*h*comp + c*w*h, darrayGetPr(values[0]) + z*w*h, w*h*sizeof(double));
- }
-#endif
-
-#ifdef TEST_EXPAND_XY
- for (c = 0; c < comp; c++) {
- darray_t *temp;
- int dimensions[3] = {w, h, k};
- for (z = 0; z < k; z++)
- memcpy(darrayGetPr(values[0]) + z*w*h, darrayGetPr(result) + z*w*h*comp + c*w*h, w*h*sizeof(double));
- dimensions[0] = (w+1)>>1;
- dimensions[1] = (h+1)>>1;
- temp = darrayCreate(3, dimensions);
- ASSERT(temp, "no memory");
- rc = shrink_xy(values[0], temp);
- ASSERT(rc >= 0, "no memory");
- expand_xy(temp, values[0]);
- darrayDestroy(temp);
- for (z = 0; z < k; z++)
- memcpy(darrayGetPr(result) + z*w*h*comp + c*w*h, darrayGetPr(values[0]) + z*w*h, w*h*sizeof(double));
- }
-#endif
-
-#ifdef TEST_EXPAND_Z
- for (c = 0; c < comp; c++) {
- darray_t *temp;
- int dimensions[3] = {w, h, k};
- for (z = 0; z < k; z++)
- memcpy(darrayGetPr(values[0]) + z*w*h, darrayGetPr(result) + z*w*h*comp + c*w*h, w*h*sizeof(double));
- dimensions[2] = (k+1)>>1;
- temp = darrayCreate(3, dimensions);
- ASSERT(temp, "no memory");
- rc = shrink_z(values[0], temp, flows_x[0], flows_y[0]);
- ASSERT(rc >= 0, "no memory");
- rc = expand_z(temp, values[0], flows_x[0], flows_y[0]);
- ASSERT(rc >= 0, "no memory");
- darrayDestroy(temp);
- for (z = 0; z < k; z++)
- memcpy(darrayGetPr(result) + z*w*h*comp + c*w*h, darrayGetPr(values[0]) + z*w*h, w*h*sizeof(double));
- }
-#endif
-
- rc = 0;
-
-error:
- /* free memory */
- for (d = 0; d < max_depth; d++) {
- if (nb_list)
- free(nb_list[d]);
- if (marks)
- darrayDestroy(marks[d]);
- if (values)
- darrayDestroy(values[d]);
- if (flows_xy_x)
- darrayDestroy(flows_xy_x[d]);
- if (flows_xy_y)
- darrayDestroy(flows_xy_y[d]);
- if (flows_x)
- darrayDestroy(flows_x[d]);
- if (flows_y)
- darrayDestroy(flows_y[d]);
- if (flows_ix)
- darrayDestroy(flows_ix[d]);
- if (flows_iy)
- darrayDestroy(flows_iy[d]);
- }
- free(nb_list);
- free(marks);
- free(values);
- free(flows_xy_x);
- free(flows_xy_y);
- free(flows_x);
- free(flows_y);
- free(flows_ix);
- free(flows_iy);
- darrayDestroy(init);
- darrayDestroy(luminance);
- if (rc < 0)
- darrayDestroy(result);
-
- return (rc < 0) ? NULL : result;
-}
--- /dev/null
+#ifndef ___MULTIGRID_H___
+#define ___MULTIGRID_H___
+/* print error message b, if a is false */
+#define ASSERT(a, b) \
+{ if (!(a)) { fprintf(stderr, "%s %s() Error in line %d: %s\n", __FILE__, __func__, __LINE__, (b)); abort(); } }
+
+#define min(a, b) ((a < b) ? a : b)
+#define max(a, b) ((a > b) ? a : b)
+
+int shrink_z(darray_t *src, darray_t *dst, const darray_t *flow);
+int shrink_flow_z(darray_t *src, darray_t *dst);
+int shrink_xy(darray_t *src, darray_t *dst);
+void quantize(const darray_t *array, double threshold, double gt, double le);
+void multiply(const darray_t *array, double factor);
+unsigned char *gen_neighbor(int w, int h, int k, const darray_t *flow, const darray_t *flow_i);
+void weighten(const darray_t *luminance, unsigned char *neighbors);
+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);
+#endif
+
int main(int argc, char *argv[])
{
- darray_t *gI = NULL, *cI = NULL, *markIm = NULL, *ntscIm = NULL, *resultIm;
- darray_t *dx = NULL, *dy = NULL, *idx = NULL, *idy = NULL;
-#ifdef WITH_OPENCV
- int *flow_map_x = NULL, *flow_map_y = NULL;
-#endif
+ darray_t *gI = NULL, *cI = NULL, *markIm = NULL, *ntscIm = NULL;
+ darray_t *flow = NULL, *flow_i = NULL;
double *ptr, *ptr2, *ptr3;
double diff, sum;
int dims[4];
- int width = 0, height = 0, w, h, k = 1, index = 0, count;
+ int w = 0, h = 0, load_w, load_h, k = 1, index = 0, z;
unsigned char *img_buffer, *img_buffer_all = NULL, *mark_buffer = NULL;
- int rc, i, ii, j;
+ int rc, x, y, i, ii;
time_t start, end;
const char* filename;
char first_filename[256];
time(&start);
change_bc = 0;
- for (count = 0; count < k; count++) {
+ for (z = 0; z < k; z++) {
if (sequence) {
- filename = sequence[count].filename;
+ filename = sequence[z].filename;
/* first_filename is set by load_sequence */
} else {
filename = argv[1];
}
#endif
// load original image and convert their RGB components to double RGB array
- rc = load_img(-1, &img_buffer, &w, &h, filename, index + count);
+ rc = load_img(-1, &img_buffer, &load_w, &load_h, filename, index + z);
if (rc) {
fprintf(stderr, "Failed to load grey image\n");
return 0;
}
- scale_img(img_buffer, w, h, scale);
- if (count == 0) {
- width = w / scale;
- height = h / scale;
+ scale_img(img_buffer, load_w, load_h, scale);
+ if (z == 0) {
+ w = load_w / scale;
+ h = load_h / scale;
}
- if (w/scale != width || h/scale != height) {
+ if (load_w/scale != w || load_h/scale != h) {
fprintf(stderr, "Error: All input images must have equal dimenstions.\n");
return 0;
}
// now we know the dimensions, so we can create input arrays
if (!gI) {
- dims[0] = width; dims[1] = height; dims[2] = 3; dims[3] = k;
+ dims[0] = w; dims[1] = h; dims[2] = 3; dims[3] = k;
gI = darrayCreate(4, dims);
}
if (!gI) {
exit (0);
}
if (!cI) {
- dims[0] = width; dims[1] = height; dims[2] = 3; dims[3] = k;
+ dims[0] = w; dims[1] = h; dims[2] = 3; dims[3] = k;
cI = darrayCreate(4, dims);
}
if (!cI) {
printf("failed to create marked image array\n");
exit (0);
}
- ptr = darrayGetPr(gI) + width*height*3*count;
- img2array(img_buffer, width, height, ptr, width, height);
+ ptr = darrayGetPr(gI) + w*h*3*z;
+ img2array(img_buffer, w, h, ptr, w, h);
#ifdef WITH_OPENCV
if (k > 1 && flow_enable) {
if (!img_buffer_all)
- img_buffer_all = malloc(width*height*3*k);
+ img_buffer_all = malloc(w*h*3*k);
if (!img_buffer_all) {
printf("failed to create grey image array\n");
exit (0);
}
- memcpy(img_buffer_all + width*height*3*count, img_buffer, width*height*3);
+ memcpy(img_buffer_all + w*h*3*z, img_buffer, w*h*3);
}
#endif
free(img_buffer);
char name[256];
unsigned char c;
// load marked mask and convert their RGB components to double YUV array
- memcpy(darrayGetPr(cI) + width*height*3*count, darrayGetPr(gI) + width*height*3*count, width*height*3 * sizeof(double));
+ memcpy(darrayGetPr(cI) + w*h*3*z, darrayGetPr(gI) + w*h*3*z, w*h*3 * sizeof(double));
/* add extra memory for unscaled data to prevent buffer overflow */
if (!mark_buffer)
- mark_buffer = (unsigned char *)malloc(width*height*k + w*h);
+ mark_buffer = (unsigned char *)malloc(w*h*k + load_w*load_h);
if (!mark_buffer) {
printf("no memory for mark buffer\n");
exit (0);
printf("failed to load palette for file: '%s'\n", filename);
exit (0);
}
- sprintf(name, filename, index + count);
+ sprintf(name, filename, index + z);
/* always load full unscaled image, then scale down */
- if (load_marked(mark_buffer + width*height*count, w, h, name) == 0) {
- scale_mark(mark_buffer + width*height*count, w, h, scale);
- ptr = darrayGetPr(cI) + width*height*3*count;
- for (i = 0; i < height; i++) {
- for (j = 0; j < width; j++) {
+ if (load_marked(mark_buffer + w*h*z, load_w, load_h, name) == 0) {
+ scale_mark(mark_buffer + w*h*z, load_w, load_h, scale);
+ ptr = darrayGetPr(cI) + w*h*3*z;
+ for (y = 0; y < h; y++) {
+ for (x = 0; x < w; x++) {
/* do not apply mask on index 0 */
- c = mark_buffer[i*width+j + width*height*count];
+ c = mark_buffer[y*w+x + w*h*z];
if (c == 0)
continue;
/* check for any brightness/contrast change */
/* do not apply white pixles, this meas: keep original color */
if (mark_palette[c-1].r == 255 && mark_palette[c-1].g == 255 && mark_palette[c-1].b == 255)
continue;
- ptr[i*width+j] = mark_palette[c-1].r / 255.0F;
- ptr[i*width+j + width*height] = mark_palette[c-1].g / 255.0F;
- ptr[i*width+j + width*height*2] = mark_palette[c-1].b / 255.0F;
+ ptr[y*w+x] = mark_palette[c-1].r / 255.0F;
+ ptr[y*w+x + w*h] = mark_palette[c-1].g / 255.0F;
+ ptr[y*w+x + w*h*2] = mark_palette[c-1].b / 255.0F;
}
}
} else
- memset(mark_buffer + width*height*count, 0, width*height);
+ memset(mark_buffer + w*h*z, 0, w*h);
} else {
// load marked image and convert their RGB components to double YUV array
- rc = load_img(-1, &img_buffer, &w, &h, argv[2], index + count);
+ rc = load_img(-1, &img_buffer, &load_w, &load_h, argv[2], index + z);
if (!rc) {
- if (w/scale != width || h/scale != height) {
+ if (load_w/scale != w || load_h/scale != h) {
fprintf(stderr, "Error: All input images must have equal dimenstions.\n");
return 0;
}
- ptr = darrayGetPr(cI) + width*height*3*count;
- img2array(img_buffer, width, height, ptr, width, height);
+ ptr = darrayGetPr(cI) + w*h*3*z;
+ img2array(img_buffer, w, h, ptr, w, h);
free(img_buffer);
} else {
fprintf(stderr, "Failed to load marked image, omitting...\n");
- memcpy(darrayGetPr(cI) + width*height*3*count, darrayGetPr(gI) + width*height*3*count, width*height*3 * sizeof(double));
+ memcpy(darrayGetPr(cI) + w*h*3*z, darrayGetPr(gI) + w*h*3*z, w*h*3 * sizeof(double));
}
}
}
/* create color mask and ntsc arrays for the colorization process */
- dims[0] = width; dims[1] = height; dims[2] = k;
+ dims[0] = w; dims[1] = h; dims[2] = k;
markIm = darrayCreate(3, dims);
if (!markIm) {
printf("failed to create mark array\n");
exit (0);
}
features = (change_bc) ? 4 : 2;
- dims[0] = width; dims[1] = height; dims[2] = features+1; dims[3] = k;
+ dims[0] = w; dims[1] = h; dims[2] = features+1; dims[3] = k;
ntscIm = darrayCreate(4, dims);
if (!ntscIm) {
printf("failed to create ntsc array\n");
exit (0);
}
- for (count = 0; count < k; count++) {
+ for (z = 0; z < k; z++) {
if (sequence || !strcmp(argv[2], "marked")) {
unsigned char c;
- ptr = darrayGetPr(markIm) + width*height*count;
+ ptr = darrayGetPr(markIm) + w*h*z;
// use marked mask to fill markIm
- for (i = 0; i < height; i++) {
- for (j = 0; j < width; j++) {
- if (j < width && i < height) {
+ for (y = 0; y < h; y++) {
+ for (x = 0; x < w; x++) {
+ if (x < w && y < h) {
/* do not apply mask on index 0 */
- c = mark_buffer[i*width+j + width*height*count];
+ c = mark_buffer[y*w+x + w*h*z];
} else
c = 0;
if (c)
- ptr[i*width+j] = 1.0F;
+ ptr[y*w+x] = 1.0F;
else
- ptr[i*width+j] = 0.0F;
+ ptr[y*w+x] = 0.0F;
}
}
} else {
// - sum all components to get grey image
// - apply threshold (pixle is 1F, if the absolute difference is > 0.01F)
// original code: markIm=(sum(abs(gI-cI),3)>0.01);
- ptr = darrayGetPr(gI) + width*height*3*count;
- ptr2 = darrayGetPr(cI) + width*height*3*count;
- ptr3 = darrayGetPr(markIm) + width*height*count;
- for (i = 0, ii = width * height; i < ii; i++) {
+ ptr = darrayGetPr(gI) + w*h*3*z;
+ ptr2 = darrayGetPr(cI) + w*h*3*z;
+ ptr3 = darrayGetPr(markIm) + w*h*z;
+ for (i = 0, ii = w * h; i < ii; i++) {
diff = 0;
sum = ptr[i] - ptr2[i];
if (sum < 0)
}
// convert original image into YUV
- ptr = darrayGetPr(gI) + width*height*3*count;
- rgb2yuv(ptr, ptr, width, height);
+ ptr = darrayGetPr(gI) + w*h*3*z;
+ rgb2yuv(ptr, ptr, w, h);
// convert maked image into YUV
- ptr = darrayGetPr(cI) + width*height*3*count;
- rgb2yuv(ptr, ptr, width, height);
+ ptr = darrayGetPr(cI) + w*h*3*z;
+ rgb2yuv(ptr, ptr, w, h);
if (test != BC_ONLY) {
if (sequence || !strcmp(argv[2], "marked")) {
unsigned char c;
// generate NTSC image: use luminance from original image and chrominance from original or marked image
- ptr = darrayGetPr(gI) + width*height*3*count;
- ptr2 = darrayGetPr(cI) + width*height*3*count;
- ptr3 = darrayGetPr(ntscIm) + width*height*(features+1)*count;
+ ptr = darrayGetPr(gI) + w*h*3*z;
+ ptr2 = darrayGetPr(cI) + w*h*3*z;
+ ptr3 = darrayGetPr(ntscIm) + w*h*(features+1)*z;
/* use original y component */
- memcpy(ptr3, ptr, width * height * sizeof(double));
+ memcpy(ptr3, ptr, w * h * sizeof(double));
/* apply new uv components */
- for (i = 0; i < height; i++) {
- for (j = 0; j < width; j++) {
- c = mark_buffer[i*width+j + width*height*count];
+ for (y = 0; y < h; y++) {
+ for (x = 0; x < w; x++) {
+ c = mark_buffer[y*w+x + w*h*z];
if (c == 0) {
- ptr3[width * height + width * i + j] = ptr[width * height + width * i + j];
- ptr3[width * height * 2 + width * i + j] = ptr[width * height * 2 + width * i + j];
+ ptr3[w * h + w * y + x] = ptr[w * h + w * y + x];
+ ptr3[w * h * 2 + w * y + x] = ptr[w * h * 2 + w * y + x];
} else {
- ptr3[width * height + width * i + j] = ptr2[width * height + width * i + j];
- ptr3[width * height * 2 + width * i + j] = ptr2[width * height * 2 + width * i + j];
+ ptr3[w * h + w * y + x] = ptr2[w * h + w * y + x];
+ ptr3[w * h * 2 + w * y + x] = ptr2[w * h * 2 + w * y + x];
}
}
}
} else {
// generate NTSC image: use luminance from original image and chrominance from maked image
- ptr = darrayGetPr(gI) + width*height*3*count;
- ptr2 = darrayGetPr(cI) + width*height*3*count;
- ptr3 = darrayGetPr(ntscIm) + width*height*(features+1)*count;
- memcpy(ptr3, ptr, width * height * sizeof(double));
- memcpy(ptr3 + width * height, ptr2 + width * height, width * height * sizeof(double));
- memcpy(ptr3 + width * height * 2, ptr2 + width * height * 2, width * height * sizeof(double));
+ ptr = darrayGetPr(gI) + w*h*3*z;
+ ptr2 = darrayGetPr(cI) + w*h*3*z;
+ ptr3 = darrayGetPr(ntscIm) + w*h*(features+1)*z;
+ memcpy(ptr3, ptr, w * h * sizeof(double));
+ memcpy(ptr3 + w * h, ptr2 + w * h, w * h * sizeof(double));
+ memcpy(ptr3 + w * h * 2, ptr2 + w * h * 2, w * h * sizeof(double));
}
} else {
/* use grey image as result if BC_ONLY test is selected */
- ptr = darrayGetPr(gI) + width*height*3*count;
- ptr2 = darrayGetPr(ntscIm) + width*height*(features+1)*count;
- memcpy(ptr2, ptr, width * height * 3 * sizeof(double));
+ ptr = darrayGetPr(gI) + w*h*3*z;
+ ptr2 = darrayGetPr(ntscIm) + w*h*(features+1)*z;
+ memcpy(ptr2, ptr, w * h * 3 * sizeof(double));
}
if (change_bc) {
unsigned char c;
- ptr2 = darrayGetPr(ntscIm) + width*height*(features+1)*count;
+ ptr2 = darrayGetPr(ntscIm) + w*h*(features+1)*z;
/* use original y component */
- memcpy(ptr2, ptr, width * height * sizeof(double));
+ memcpy(ptr2, ptr, w * h * sizeof(double));
/* apply brightness and contrast from makred pixles to uv components of grey image */
- for (i = 0; i < height; i++) {
- for (j = 0; j < width; j++) {
+ for (y = 0; y < h; y++) {
+ for (x = 0; x < w; x++) {
/* use unchanged brightness and contrast on index 0 */
- c = mark_buffer[i*width+j + width*height*count];
+ c = mark_buffer[y*w+x + w*h*z];
if (c == 0) {
- ptr2[i*width+j + width*height*3] = 0;
- ptr2[i*width+j + width*height*4] = 0.1;
+ ptr2[y*w+x + w*h*3] = 0;
+ ptr2[y*w+x + w*h*4] = 0.1;
} else {
- ptr2[i*width+j + width*height*3] = mark_palette[c-1].bright / 10.0;
- ptr2[i*width+j + width*height*4] = mark_palette[c-1].contrast / 10.0;
+ ptr2[y*w+x + w*h*3] = mark_palette[c-1].bright / 10.0;
+ ptr2[y*w+x + w*h*4] = mark_palette[c-1].contrast / 10.0;
}
}
}
}
}
- /* create flow vectors */
- dims[0] = width; dims[1] = height; dims[2] = k - 1;
- dx = darrayCreate(3, dims);
- if (!dx) {
- printf("failed to create array\n");
- exit (0);
- }
- dy = darrayCreate(3, dims);
- if (!dy) {
- printf("failed to create array\n");
- exit (0);
- }
- idx = darrayCreate(3, dims);
- if (!idx) {
- printf("failed to create array\n");
- exit (0);
- }
- idy = darrayCreate(3, dims);
- if (!idy) {
- printf("failed to create array\n");
- exit (0);
- }
#ifdef WITH_OPENCV
if (k > 1 && flow_enable) {
- printf("Calculating optical flow for %d frames: steps=%d win=%d max=%d\n", k, flow_steps/scale, flow_win/scale, flow_max/scale);
- flow_map_x = (int *)malloc(sizeof(int)*width*height);
- if (!flow_map_x) {
- printf("failed to alloc array\n");
+ /* create flow vectors */
+ dims[0] = w; dims[1] = h; dims[2] = k - 1; dims[3] = 2;
+ flow = darrayCreate(4, dims);
+ if (!flow) {
+ printf("failed to create array\n");
exit (0);
}
- flow_map_y = (int *)malloc(sizeof(int)*width*height);
- if (!flow_map_y) {
- printf("failed to alloc array\n");
+ flow_i = darrayCreate(4, dims);
+ if (!flow_i) {
+ printf("failed to create array\n");
exit (0);
}
+ printf("Calculating optical flow for %d frames: window=%d\n", k, flow_window/scale);
}
- for (count = 0; count < k-1; count++) {
- if (flow_map_x && flow_map_y)
- create_flow_maps(NULL, img_buffer_all + width*height*3*(count+1), img_buffer_all + width*height*3*count, width, height, flow_steps/scale, flow_win/scale, flow_max/scale, NULL, NULL, flow_map_x, flow_map_y, NULL);
- ptr = darrayGetPr(dx) + width*height*count;
- if (flow_map_x) {
- for (j = 0; j < height; j++) {
- for (i = 0; i < width; i++)
- ptr[i + width * j] = flow_map_x[width*j + i];
- }
- }
- ptr = darrayGetPr(dy) + width*height*count;
- if (flow_map_y) {
- for (j = 0; j < height; j++) {
- for (i = 0; i < width; i++)
- ptr[i + width * j] = flow_map_y[width*j + i];
- }
- }
- if (flow_map_x && flow_map_y)
- create_flow_maps(img_buffer_all + width*height*3*count, NULL, img_buffer_all + width*height*3*(count+1), width, height, flow_steps/scale, flow_win/scale, flow_max/scale, flow_map_x, flow_map_y, NULL, NULL, NULL);
- ptr = darrayGetPr(idx) + width*height*count;
- if (flow_map_x) {
- for (j = 0; j < height; j++) {
- for (i = 0; i < width; i++)
- ptr[i + width * j] = flow_map_x[width*j + i];
- }
- }
- ptr = darrayGetPr(idy) + width*height*count;
- if (flow_map_y) {
- for (j = 0; j < height; j++) {
- for (i = 0; i < width; i++)
- ptr[i + width * j] = flow_map_y[width*j + i];
- }
- }
+ for (z = 0; z < k-1; z++) {
+ if (flow)
+ create_flow_maps(NULL, img_buffer_all + w*h*3*(z+1), img_buffer_all + w*h*3*z, w, h, flow_window/scale, 0, NULL, NULL, darrayGetPr(flow) + w*h*z, darrayGetPr(flow) + w*h*z + w*h*(k-1), NULL);
+ if (flow_i)
+ create_flow_maps(img_buffer_all + w*h*3*z, NULL, img_buffer_all + w*h*3*(z+1), w, h, flow_window/scale, 0, darrayGetPr(flow_i) + w*h*z, darrayGetPr(flow_i) + w*h*z + w*h*(k-1), NULL, NULL, NULL);
}
#endif
free(img_buffer_all);
if (test != FLOW_NEXT && test != FLOW_PREV && test != MARKED && test != MASK && test != MASK_COLOR && test != BC_ONLY && test != BC_IMAGE) {
printf("Colorizing %d frames, please wait...\n", k);
- resultIm = colorize(ntscIm, markIm, dx, dy, idx, idy, in_itr_num, out_itr_num, scalexyz);
- if (!resultIm) {
+ rc = colorize(ntscIm, markIm, flow, flow_i, in_itr_num, out_itr_num, scalexyz);
+ if (rc < 0) {
if (k > 1)
printf("No memory! Use smaller frames or less frames between key frames or add more memory.");
else
printf("No memory! Use smaller image or add more memory.");
exit(-1);
}
- } else
- resultIm = darrayClone(ntscIm);
+ }
- /* if we have a change, we apply brightness+contrast from bcIm */
+ /* if we have a change, we apply brightness+contrast from ntscIm */
if (change_bc) {
- for (count = 0; count < k; count++) {
- ptr = darrayGetPr(ntscIm) + width*height*(features+1)*count;
- ptr2 = darrayGetPr(resultIm) + width*height*(features+1)*count;
- for (i = 0; i < height; i++) {
- for (j = 0; j < width; j++) {
+ for (z = 0; z < k; z++) {
+ ptr = darrayGetPr(ntscIm) + w*h*(features+1)*z;
+ ptr2 = darrayGetPr(ntscIm) + w*h*(features+1)*z;
+ for (y = 0; y < h; y++) {
+ for (x = 0; x < w; x++) {
/* apply contrast */
- ptr2[width * i + j] = (ptr2[width * i + j] - 0.5) * ptr[width * height * 4 + width * i + j] * 10.0 + 0.5;
+ ptr2[w * y + x] = (ptr2[w * y + x] - 0.5) * ptr[w * h * 4 + w * y + x] * 10.0 + 0.5;
/* apply brightness */
- ptr2[width * i + j] += ptr[width * height * 3 + width * i + j] * 10.0;
- if (ptr2[width * i + j] < 0)
- ptr2[width * i + j] = 0;
- if (ptr2[width * i + j] > 1)
- ptr2[width * i + j] = 1;
+ ptr2[w * y + x] += ptr[w * h * 3 + w * y + x] * 10.0;
+ if (ptr2[w * y + x] < 0)
+ ptr2[w * y + x] = 0;
+ if (ptr2[w * y + x] > 1)
+ ptr2[w * y + x] = 1;
#if 0
#warning TEST: show brightness and contrast change as uv vectors on a grey array */
-ptr2[width * i + j] = 0.5;
-ptr2[width * height + width * i + j] = ptr[width * height + width * i + j] * 10;
-ptr2[width * height * 2 + width * i + j] = ptr[width * height * 2 + width * i + j] * 10 - 1;
+ptr2[w * y + x] = 0.5;
+ptr2[w * h + w * y + x] = ptr[w * h + w * y + x] * 10;
+ptr2[w * h * 2 + w * y + x] = ptr[w * h * 2 + w * y + x] * 10 - 1;
#endif
}
}
}
}
- darrayDestroy(ntscIm);
- ntscIm = NULL;
-
#ifdef WITH_OPENCV
if (test == FLOW_NEXT || test == FLOW_PREV) {
/* apply flow planes to result image as u and y vector */
- for (count = 0; count < k; count++) {
- ptr = darrayGetPr(resultIm) + width*height*(features+1)*count;
+ for (z = 0; z < k; z++) {
+ ptr = darrayGetPr(ntscIm) + w*h*(features+1)*z;
if (test == FLOW_NEXT) {
- ptr2 = darrayGetPr(dx) + width*height*count;
- ptr3 = darrayGetPr(dy) + width*height*count;
+ ptr2 = darrayGetPr(flow) + w*h*z;
+ ptr3 = darrayGetPr(flow) + w*h*z*(k-1);
} else {
- ptr2 = darrayGetPr(idx) + width*height*count;
- ptr3 = darrayGetPr(idy) + width*height*count;
+ ptr2 = darrayGetPr(flow_i) + w*h*z;
+ ptr3 = darrayGetPr(flow_i) + w*h*z*(k-1);
}
- for (i = 0; i < height; i++) {
- for (j = 0; j < width; j++) {
- ptr[width * i + j] = 0.5;
- if (count < k-1) {
- ptr[width * i + j + width*height] = ptr2[width * i + j] / flow_max;
- ptr[width * i + j + width*height*2] = ptr3[width * i + j] / flow_max;
+ for (y = 0; y < h; y++) {
+ for (x = 0; x < w; x++) {
+ ptr[w * y + x] = 0.5;
+ if (z < k-1) {
+ ptr[w * y + x + w*h] = ptr2[w * y + x] / 50;
+ ptr[w * y + x + w*h*2] = ptr3[w * y + x] / 50;
} else {
- ptr[width * i + j + width*height] = 0;
- ptr[width * i + j + width*height*2] = 0;
+ ptr[w * y + x + w*h] = 0;
+ ptr[w * y + x + w*h*2] = 0;
}
}
}
if (test == MASK) {
/* apply maked mask as image */
- for (count = 0; count < k; count++) {
- ptr = darrayGetPr(resultIm) + width*height*(features+1)*count;
- ptr2 = darrayGetPr(markIm) + width*height*count;
- for (i = 0; i < height; i++) {
- for (j = 0; j < width; j++) {
- ptr[width * i + j] = ptr2[width * i + j];
- ptr[width * i + j + width*height] = ptr2[width * i + j];
- ptr[width * i + j + width*height*2] = ptr2[width * i + j];
+ for (z = 0; z < k; z++) {
+ ptr = darrayGetPr(ntscIm) + w*h*(features+1)*z;
+ ptr2 = darrayGetPr(markIm) + w*h*z;
+ for (y = 0; y < h; y++) {
+ for (x = 0; x < w; x++) {
+ ptr[w * y + x] = ptr2[w * y + x];
+ ptr[w * y + x + w*h] = ptr2[w * y + x];
+ ptr[w * y + x + w*h*2] = ptr2[w * y + x];
}
}
}
if (test == MASK_COLOR) {
/* apply maked mask on grey image */
- for (count = 0; count < k; count++) {
- ptr = darrayGetPr(resultIm) + width*height*(features+1)*count;
- ptr2 = darrayGetPr(markIm) + width*height*count;
- for (i = 0; i < height; i++) {
- for (j = 0; j < width; j++) {
+ for (z = 0; z < k; z++) {
+ ptr = darrayGetPr(ntscIm) + w*h*(features+1)*z;
+ ptr2 = darrayGetPr(markIm) + w*h*z;
+ for (y = 0; y < h; y++) {
+ for (x = 0; x < w; x++) {
/* darken unmarked areas, make maked areas with uniformed brightness */
- if (ptr2[width * i + j] < 0.5F)
- ptr[width * i + j] = ptr[width * i + j] / 4;
+ if (ptr2[w * y + x] < 0.5F)
+ ptr[w * y + x] = ptr[w * y + x] / 4;
else
- ptr[width * i + j] = 0.5F;
+ ptr[w * y + x] = 0.5F;
}
}
}
if (test == BC_IMAGE) {
/* apply bc image as result image */
- for (count = 0; count < k; count++) {
- ptr = darrayGetPr(resultIm) + width*height*(features+1)*count;
- ptr2 = darrayGetPr(resultIm) + width*height*(features+1)*count;
+ for (z = 0; z < k; z++) {
+ ptr = darrayGetPr(ntscIm) + w*h*(features+1)*z;
+ ptr2 = darrayGetPr(ntscIm) + w*h*(features+1)*z;
/* use uniformed brightness as y component */
- for (i = 0; i < height; i++) {
- for (j = 0; j < width; j++) {
- ptr[width * i + j] = 0.5F;
+ for (y = 0; y < h; y++) {
+ for (x = 0; x < w; x++) {
+ ptr[w * y + x] = 0.5F;
}
}
if (change_bc) {
- ptr2 = darrayGetPr(resultIm) + width*height*(features+1)*count;
- memcpy(ptr+width*height, ptr2+width*height*3, width * height * sizeof(double));
- memcpy(ptr+width*height*2, ptr2+width*height*4, width * height * sizeof(double));
+ ptr2 = darrayGetPr(ntscIm) + w*h*(features+1)*z;
+ memcpy(ptr+w*h, ptr2+w*h*3, w * h * sizeof(double));
+ memcpy(ptr+w*h*2, ptr2+w*h*4, w * h * sizeof(double));
} else {
- memset(ptr+width*height, 0, width * height * sizeof(double));
- memset(ptr+width*height*2, 0, width * height * sizeof(double));
+ memset(ptr+w*h, 0, w * h * sizeof(double));
+ memset(ptr+w*h*2, 0, w * h * sizeof(double));
}
}
}
// save result YUV array to image with RGB components
- img_buffer = (unsigned char *)malloc(width * height * 3);
+ img_buffer = (unsigned char *)malloc(w * h * 3);
if (!img_buffer) {
fprintf(stderr, "Failed to allocate image buffer\n");
return 0;
}
- for (count = 0; count < k; count++) {
- ptr = darrayGetPr(resultIm) + width*height*(features+1)*count;
- yuv2rgb(ptr, ptr, width, height);
- array2img(ptr, width, height, img_buffer, width, height);
+ for (z = 0; z < k; z++) {
+ ptr = darrayGetPr(ntscIm) + w*h*(features+1)*z;
+ yuv2rgb(ptr, ptr, w, h);
+ array2img(ptr, w, h, img_buffer, w, h);
if (sequence) {
static char name[256], *p, *q;
- p = sequence[count].filename;
+ p = sequence[z].filename;
while((q = strchr(p, DIR_SEPERATOR)))
p = q + 1;
- strcpy(name, sequence[count].filename);
- name[p - sequence[count].filename] = '\0';
+ strcpy(name, sequence[z].filename);
+ name[p - sequence[z].filename] = '\0';
strcat(name, "colorized_");
strcat(name, p);
filename = name;
} else
filename = argv[3];
- save_img(img_buffer, width, height, filename, index + count);
+ save_img(img_buffer, w, h, filename, index + z);
}
free(img_buffer);
printf("Elapsed time: %d minutes, %d seconds\n", (int)(end-start)/60, (int)(end-start)%60);
// destroy
- darrayDestroy(resultIm);
- resultIm = NULL;
- darrayDestroy(dx);
- dx = NULL;
- darrayDestroy(dy);
- dy = NULL;
- darrayDestroy(idx);
- idx = NULL;
- darrayDestroy(idy);
- idy = NULL;
-#ifdef WITH_OPENCV
- free(flow_map_x);
- flow_map_x = NULL;
- free(flow_map_y);
- flow_map_y = NULL;
-#endif
+ darrayDestroy(ntscIm);
+ ntscIm = NULL;
+ darrayDestroy(flow);
+ flow = NULL;
+ darrayDestroy(flow_i);
+ flow_i = NULL;
darrayDestroy(markIm);
markIm = NULL;
#include <stdio.h>
#include "ppm.h"
#include "opticalflow.h"
+#include "yuv.h"
#include <opencv2/core/core_c.h>
#include <opencv2/imgproc/imgproc_c.h>
#include <opencv2/video/tracking.hpp>
#include "dir_seperator.h"
-int flow_enable = 0, flow_steps, flow_win, flow_max;
+int flow_enable = 0, flow_window, flow_view_vector, flow_view_uv;
void flow_default(void)
{
- flow_steps = 8;
- flow_win = 30;
- flow_max = 16;
+ flow_window = 15;
+ flow_view_vector = 16;
+ flow_view_uv = 0; /* disabled */
}
/* save flow setting */
return;
}
fprintf(fp, "flow_enable %d\n", flow_enable);
- fprintf(fp, "flow_steps %d\n", flow_steps);
- fprintf(fp, "flow_win %d\n", flow_win);
- fprintf(fp, "flow_max %d\n", flow_max);
+ fprintf(fp, "flow_window %d\n", flow_window);
+ fprintf(fp, "flow_view_vector %d\n", flow_view_vector);
+ fprintf(fp, "flow_view_uv %d\n", flow_view_uv);
fclose(fp);
}
sscanf(p, "%d", &flow_enable);
continue;
}
- if (!strncmp(line, "flow_steps ", 11)) {
- p = line + 11;
- sscanf(p, "%d", &flow_steps);
+ if (!strncmp(line, "flow_window ", 12)) {
+ p = line + 12;
+ sscanf(p, "%d", &flow_window);
continue;
}
- if (!strncmp(line, "flow_win ", 9)) {
- p = line + 9;
- sscanf(p, "%d", &flow_win);
+ if (!strncmp(line, "flow_view_vector ", 17)) {
+ p = line + 17;
+ sscanf(p, "%d", &flow_view_vector);
continue;
}
- if (!strncmp(line, "flow_max ", 9)) {
- p = line + 9;
- sscanf(p, "%d", &flow_max);
+ if (!strncmp(line, "flow_view_uv ", 13)) {
+ p = line + 13;
+ sscanf(p, "%d", &flow_view_uv);
continue;
}
}
return 0;
}
-/* smooth image, it is safe if image pointer is NULL pointer */
-static IplImage *smooth_image(IplImage *image, int width, int height, int dim)
-{
- IplImage *temp;
- CvSize size;
-
- if (!image)
- return NULL;
-
- size.width = width;
- size.height = height;
-
- temp = cvCreateImage(size, IPL_DEPTH_8U, 1);
- if (!temp)
- return image;
-
- /* use odd dimenstions of convolution matrix */
- cvSmooth(image, temp, CV_GAUSSIAN, dim | 1, dim | 1, 0, 0);
- cvReleaseImage(&image);
- return temp;
-}
-
-/* remove bad features and replace them by average flow */
-static int filter_max_flow(int num, char *features_found, CvPoint2D32f *corners_from, CvPoint2D32f *corners_to, int max_flow)
+/* draw flow matrix as lines */
+static void draw_lines_from_matrix(CvMat *flow, int width, int height, int steps, int color, IplImage *image)
{
- int i;
- int n = 0;
- int average_x = 0;
- int average_y = 0;
-
- for (i = 0; i < num; i++) {
- if (features_found[i]) {
- if (max_flow
- && sqrt(pow(corners_from[i].x-corners_to[i].x, 2) +
- pow(corners_from[i].y-corners_to[i].y, 2)) > max_flow) {
- features_found[i] = 0;
- continue;
- }
- average_x += corners_to[i].x - corners_from[i].x;
- average_y += corners_to[i].y - corners_from[i].y;
- n++;
- }
- }
-
- if (n > 1) {
- average_x /= n;
- average_y /= n;
- }
- /* replace bad vectors by average */
- for (i = 0; i < num; i++) {
- if (!features_found[i]) {
- corners_to[i].x = corners_from[i].x + average_x;
- corners_to[i].y = corners_from[i].y + average_y;
- }
- }
- return n;
-}
-
-/* draw motion vectors as lines */
-static void draw_lines(int num, char *features_found, CvPoint2D32f *corners_from, CvPoint2D32f *corners_to, int color, IplImage *image)
-{
- int i;
-
- for (i = 0; i < num; i++) {
- if (features_found[i]) {
- cvLine(image, cvPoint(cvRound(corners_from[i].x), cvRound(corners_from[i].y)), cvPoint(cvRound(corners_to[i].x), cvRound(corners_to[i].y)), CV_RGB(color,color,color), 1, 8, 0);
+ int flow_width = width / steps;
+ int flow_height = height / steps;
+ int i, j, x, y;
+ const float *f;
+
+ for (i = 0, y = steps/2; i < flow_height && y < height; i++, y += steps) {
+ f = (const float *)(flow->data.ptr + flow->step * y);
+ for (j = 0, x = steps/2; j < flow_width && x < width; j++, x += steps) {
+ cvLine(image,
+ cvPoint(cvRound(x), cvRound(y)),
+ cvPoint(cvRound(f[2*x] + x), cvRound(f[2*x+1] + y)),
+ CV_RGB(color,color,color),
+ 1, 8, 0);
}
}
}
-/* take vectors to generate a map
- *
- * the vectors are linear interpolated to fill all pixles */
-static void scale_flow_map(int *step2_x, int *step2_y, CvPoint2D32f *corners_from, CvPoint2D32f *corners_to, int width, int height, int flow_width, int flow_height, int steps)
+/* take flow matrix to generate a map */
+static void copy_matrix_to_map(CvMat *flow, int width, int height, double *flow_map_x, double *flow_map_y)
{
- /* step 1: interpolate horizontally */
- float step1_x[width*flow_height];
- float step1_y[width*flow_height];
- int i, j, k, pos;
- float from_x, from_y, to_x, to_y;
-
- /* interpolate horizontally */
- for (i = 0; i < flow_height; i++) { // count vertical vectors
- /* get first vector position */
- pos = (width - (steps * (flow_width - 1))) / 2;
- from_x = corners_to[i * flow_width + 0].x - corners_from[i * flow_width + 0].x;
- from_y = corners_to[i * flow_width + 0].y - corners_from[i * flow_width + 0].y;
- /* use left most vector to fill left of first vector */
- for (k = 0; k < pos; k++) { // count pixles before first vector
- step1_x[i * width + k] = from_x;
- step1_y[i * width + k] = from_y;
- }
- /* interpolate between vectors */
- for (j = 0; j < flow_width - 1; j++) { // count horizontal vectors
- to_x = corners_to[i * flow_width + j+1].x - corners_from[i * flow_width + j+1].x;
- to_y = corners_to[i * flow_width + j+1].y - corners_from[i * flow_width + j+1].y;
- for (k = 0; k < steps; k++) { // count from vector to next vector
- step1_x[i * width + pos] = from_x * (steps-k) / steps + to_x * k / steps;
- step1_y[i * width + pos] = from_y * (steps-k) / steps + to_y * k / steps;
- pos++;
- }
- /* shift vector for next loop */
- from_x = to_x;
- from_y = to_y;
- }
- /* use right most vector to fill right of last vector */
- for (k = pos; k < width; k++) { // count pixles after last vector
- step1_x[i * width + k] = from_x;
- step1_y[i * width + k] = from_y;
- pos++;
- }
- if (pos != width) {
- printf("pos != with\n");
- abort();
- }
- }
-
- /* interpolate vertically */
- for (i = 0; i < width; i++) { // count horizontal vectors
- /* get first vector position */
- pos = (height - (steps * (flow_height - 1))) / 2;
- from_x = step1_x[0 * width + i];
- from_y = step1_y[0 * width + i];
- /* use top most vector to fill above first vector */
- for (k = 0; k < pos; k++) { // count pixles above first vector
- step2_x[k * width + i] = from_x;
- step2_y[k * width + i] = from_y;
- }
- /* interpolate between vectors */
- for (j = 0; j < flow_height - 1; j++) { // count vertical vectors
- to_x = step1_x[(j+1) * width + i];
- to_y = step1_y[(j+1) * width + i];
- for (k = 0; k < steps; k++) { // count from vector to next vector
- step2_x[pos * width + i] = from_x * (steps-k) / steps + to_x * k / steps;
- step2_y[pos * width + i] = from_y * (steps-k) / steps + to_y * k / steps;
- pos++;
- }
- /* shift vector for next loop */
- from_x = to_x;
- from_y = to_y;
- }
- /* use bottom vector to fill below of last vector */
- for (k = pos; k < height; k++) { // count pixles after last vector
- step2_x[pos * width + i] = from_x;
- step2_y[pos * width + i] = from_y;
- pos++;
- }
- if (pos != height) {
- printf("pos != height\n");
- abort();
+ int x, y;
+ const float *f;
+
+ for (y = 0; y < height; y++) {
+ f = (const float *)(flow->data.ptr + flow->step * y);
+ for (x = 0; x < width; x++) {
+ *flow_map_x++ = *f++;
+ *flow_map_y++ = *f++;
}
}
}
* flow_map_y_next = flow map for delta y of next image (result is stored here)
*
*/
-void *create_flow_maps(const unsigned char *img_prev_buffer, const unsigned char *img_next_buffer, unsigned char *img_buffer, int width, int height, int steps, int win_size, int max_flow, int *flow_map_x_prev, int *flow_map_y_prev, int *flow_map_x_next, int *flow_map_y_next, void *_image_preview)
+void *create_flow_maps(const unsigned char *img_prev_buffer, const unsigned char *img_next_buffer, unsigned char *img_buffer, int width, int height, int win_size, int steps, double *flow_map_x_prev, double *flow_map_y_prev, double *flow_map_x_next, double *flow_map_y_next, void *_image_preview)
{
IplImage *image_prev = NULL, *image = NULL, *image_next = NULL, *image_preview = _image_preview;
- IplImage *pyr_a = NULL, *pyr_b = NULL;
- CvSize size, pyr_sz;
- int flow_width = width / steps;
- int flow_height = height / steps;
- CvPoint2D32f corners_input[flow_width * flow_height];
- CvPoint2D32f corners_prev[flow_width * flow_height];
- CvPoint2D32f corners_next[flow_width * flow_height];
- char features_found_prev[flow_width * flow_height];
- char features_found_next[flow_width * flow_height];
- int input_n;
+ CvMat *flow;
+ CvSize size;
int i, j;
+ const float pyrScale = 0.5;
+ const float levels = 3;
+ const float winsize = win_size;
+ const float iterations = 8;
+ const float polyN = 5;
+ const float polySigma = 1.2;
+ const int flags = 0;
+
+ flow = cvCreateMat(height, width, CV_32FC2);
+ if (!flow) {
+ fprintf(stderr, "failed to allocate opencv matrix\n");
+ return NULL;
+ }
+
/* apply images to opencv image */
size.width = width;
size.height = height;
}
}
-#if 1
- /* smooth image, it is safe if image pointer is NULL pointer */
- image_prev = smooth_image(image_prev, width, height, steps);
- image = smooth_image(image, width, height, steps);
- image_next = smooth_image(image_next, width, height, steps);
- if (image_preview)
- image_preview = smooth_image(image_preview, width, height, steps);
-#endif
-
/* optical flow and compose previous and next image with current image */
- pyr_sz.width = image->width;// + 8;
- pyr_sz.height = image->height;// / 3;
-
- pyr_a = cvCreateImage(pyr_sz, IPL_DEPTH_8U, 1 );
- pyr_b = cvCreateImage(pyr_sz, IPL_DEPTH_8U, 1 );
-
if (image) {
- /* define corners */
- input_n = 0;
- if ((image_prev || image_next)) {
- for (i = 0; i < flow_height; i++) {
- for (j = 0; j < flow_width; j++) {
- corners_input[input_n].x = (width - (steps * (flow_width - 1))) / 2 + (j * steps);
- corners_input[input_n].y = (height - (steps * (flow_height - 1))) / 2 + (i * steps);
- input_n++;
- }
- }
- }
-
/* calc flow */
if (image_prev) {
- cvCalcOpticalFlowPyrLK(image, image_prev, pyr_a, pyr_b, corners_input, corners_prev, input_n, cvSize(win_size | 1, win_size | 1), 5, features_found_prev, NULL, cvTermCriteria( CV_TERMCRIT_ITER | CV_TERMCRIT_EPS, 20, 0.3), 0);
- filter_max_flow(input_n, features_found_prev, corners_input, corners_prev, max_flow);
- if (image_preview)
- draw_lines(input_n, features_found_prev, corners_input, corners_prev, 255, image_preview);
+ cvCalcOpticalFlowFarneback(image, image_prev, flow, pyrScale, levels, winsize, iterations, polyN, polySigma, flags);
+ if (image_preview && steps > 1)
+ draw_lines_from_matrix(flow, width, height, steps, 255, image_preview);
if (flow_map_x_prev && flow_map_y_prev)
- scale_flow_map(flow_map_x_prev, flow_map_y_prev, corners_input, corners_prev, width, height, flow_width, flow_height, steps);
+ copy_matrix_to_map(flow, width, height, flow_map_x_prev, flow_map_y_prev);
}
if (image_next) {
- cvCalcOpticalFlowPyrLK(image, image_next, pyr_a, pyr_b, corners_input, corners_next, input_n, cvSize(win_size | 1, win_size | 1), 5, features_found_next, NULL, cvTermCriteria( CV_TERMCRIT_ITER | CV_TERMCRIT_EPS, 20, 0.3), 0);
- filter_max_flow(input_n, features_found_next, corners_input, corners_next, max_flow);
- if (image_preview)
- draw_lines(input_n, features_found_next, corners_input, corners_next, 192, image_preview);
+ cvCalcOpticalFlowFarneback(image, image_next, flow, pyrScale, levels, winsize, iterations, polyN, polySigma, flags);
+ if (image_preview && steps > 1)
+ draw_lines_from_matrix(flow, width, height, steps, 192, image_preview);
if (flow_map_x_next && flow_map_y_next)
- scale_flow_map(flow_map_x_next, flow_map_y_next, corners_input, corners_next, width, height, flow_width, flow_height, steps);
+ copy_matrix_to_map(flow, width, height, flow_map_x_next, flow_map_y_next);
}
}
- if (pyr_a)
- cvReleaseImage(&pyr_a);
- if (pyr_b)
- cvReleaseImage(&pyr_b);
-
if (image)
cvReleaseImage(&image);
if (image_prev)
if (image_next)
cvReleaseImage(&image_next);
+ cvReleaseMat(&flow);
+
return image_preview;
}
* filename_next = name of next image
* img_buffer = current image buffer
* width, height = dimensions of image buffer
- * steps = pixles between each feature = size of grid to track
* win_size = flow algorithm window size
- * max_flow = maximum length of allowed flow vector (to avoid errors)
+ * steps = pixles between each feature to visualize vectors (1 to disable)
+ * uv_scale = apply vector (pointing to next plane) by applying uv-color (0 to disable)
*
*/
-void create_flow_view(const char *filename_prev, const char *filename_next, unsigned char *img_buffer, int width, int height, int steps, int win_size, int max_flow)
+void create_flow_view(const char *filename_prev, const char *filename_next, unsigned char *img_buffer, int width, int height, int win_size, int steps, int uv_scale)
{
unsigned char *img_prev_buffer = NULL, *img_next_buffer = NULL;
- int *flow_map_x_prev = NULL, *flow_map_y_prev = NULL, *flow_map_x_next = NULL, *flow_map_y_next = NULL;
+ double *flow_map_x_prev = NULL, *flow_map_y_prev = NULL, *flow_map_x_next = NULL, *flow_map_y_next = NULL;
IplImage *image_preview = NULL;
CvSize size;
- unsigned char c, r, g, b;
+ unsigned char c;
+ double r, g, b;
int w, h, rc;
int i, j;
image_preview = cvCreateImage(size, IPL_DEPTH_8U, 1);
}
- flow_map_x_prev = (int *)malloc(sizeof(int)*width*height);
- flow_map_y_prev = (int *)malloc(sizeof(int)*width*height);
- flow_map_x_next = (int *)malloc(sizeof(int)*width*height);
- flow_map_y_next = (int *)malloc(sizeof(int)*width*height);
+ flow_map_x_prev = (double *)malloc(sizeof(double)*width*height);
+ flow_map_y_prev = (double *)malloc(sizeof(double)*width*height);
+ flow_map_x_next = (double *)malloc(sizeof(double)*width*height);
+ flow_map_y_next = (double *)malloc(sizeof(double)*width*height);
if (image_preview)
- image_preview = create_flow_maps(img_prev_buffer, img_next_buffer, img_buffer, width, height, steps, win_size, max_flow, flow_map_x_prev, flow_map_y_prev, flow_map_x_next, flow_map_y_next, image_preview);
+ image_preview = create_flow_maps(img_prev_buffer, img_next_buffer, img_buffer, width, height, win_size, steps, flow_map_x_prev, flow_map_y_prev, flow_map_x_next, flow_map_y_next, image_preview);
if (img_buffer && image_preview) {
/* paint vectors from image_preview to img_buffer */
for (j = 0; j < width; j++) {
c = CV_IMAGE_ELEM(image_preview, uchar, i, j);
if (c > 224) {
- r = 255;
+ r = 1;
g = 0;
b = 0;
} else if (c > 160) {
r = 0;
- g = 255;
+ g = 1;
b = 0;
} else {
-#if 0
- if (img_next_buffer) {
- r = g = b = c / 2 + 64;
- r += flow_map_x_next[width * i + j] * 64 / max_flow;
- b += flow_map_y_next[width * i + j] * 64 / max_flow;
- } else
-#endif
- {
- r = c;
- g = c;
- b = c;
+ if (!steps && uv_scale)
+ c += 64;
+ if (img_next_buffer && uv_scale) {
+ yuv2rgb_pixle(
+ (double)c / 255.0,
+ flow_map_x_next[width * i + j] / uv_scale,
+ flow_map_y_next[width * i + j] / uv_scale,
+ &r, &g, &b);
+ } else {
+ r = (double)c / 255.0;
+ g = (double)c / 255.0;
+ b = (double)c / 255.0;
}
}
- img_buffer[(i*width+j)*3 + 0] = r;
- img_buffer[(i*width+j)*3 + 1] = g;
- img_buffer[(i*width+j)*3 + 2] = b;
+ img_buffer[(i*width+j)*3 + 0] = r * 255.0;
+ img_buffer[(i*width+j)*3 + 1] = g * 255.0;
+ img_buffer[(i*width+j)*3 + 2] = b * 255.0;
}
}
}
-extern int flow_enable, flow_steps, flow_win, flow_max;
+extern int flow_enable, flow_window, flow_levels, flow_iter, flow_view_vector, flow_view_uv;
void flow_default(void);
void save_flow(const char *filename);
int load_flow(const char *filename);
-void *create_flow_maps(const unsigned char *img_prev_buffer, const unsigned char *img_next_buffer, unsigned char *img_buffer, int width, int height, int steps, int win_size, int max_flow, int *flow_map_x_prev, int *flow_map_y_prev, int *flow_map_x_next, int *flow_map_y_next, void *_image_preview);
+void *create_flow_maps(const unsigned char *img_prev_buffer, const unsigned char *img_next_buffer, unsigned char *img_buffer, int width, int height, int win_size, int steps, double *flow_map_x_prev, double *flow_map_y_prev, double *flow_map_x_next, double *flow_map_y_next, void *_image_preview);
-void create_flow_view(const char *filename_prev, const char *filename_next, unsigned char *img_buffer, int width, int height, int steps, int win_size, int max_flow);
+void create_flow_view(const char *filename_prev, const char *filename_next, unsigned char *img_buffer, int width, int height, int win_size, int steps, int uv_scale);