1 /*
2 # This file is part of the Astrometry.net suite.
3 # Licensed under a 3-clause BSD style license - see LICENSE
4 */
5
6 // This file gets #included in dfind.c, with
7 // DFIND2 and IMGTYPE defined appropriately.
8 // I did this so that we can handle int* and
9 // unsigned char* images using the same code.
10
11 #include "os-features.h"
12 #include "errors.h"
13 #include "log.h"
14
DFIND2(const IMGTYPE * image,int nx,int ny,int * object,int * pnobjects)15 int DFIND2(const IMGTYPE* image,
16 int nx,
17 int ny,
18 int* object,
19 int* pnobjects) {
20 int ix, iy, i;
21 int maxgroups = initial_max_groups;
22 dimage_label_t *equivs = malloc(sizeof(dimage_label_t) * maxgroups);
23 int maxlabel = 0;
24
25 /* Keep track of 'on' pixels to avoid later rescanning */
26 il* on_pixels = il_new(256);
27
28 /* Find blobs and track equivalences */
29 for (iy = 0; iy < ny; iy++) {
30 for (ix = 0; ix < nx; ix++) {
31 int thislabel, thislabelmin;
32
33 object[nx*iy+ix] = -1;
34 if (!image[nx*iy+ix])
35 continue;
36
37 /* Store location of each 'on' pixel. */
38 il_append(on_pixels, nx * iy + ix);
39
40 /* If the pixel to the left [exists and] is on, this pixel
41 joins its group. */
42 if (ix && image[nx*iy+ix-1]) {
43 /* Old group */
44 object[nx*iy+ix] = object[nx*iy+ix-1];
45
46 } else {
47 /* New blob */
48 // FIXME this part should become uf_new_group()
49 if (maxlabel >= maxgroups) {
50 maxgroups *= 2;
51 equivs = realloc(equivs, sizeof(dimage_label_t) * maxgroups);
52 assert(equivs);
53 }
54 object[nx*iy+ix] = maxlabel;
55 equivs[maxlabel] = maxlabel;
56 maxlabel++;
57
58 if (maxlabel == LABEL_MAX) {
59 logverb("Ran out of labels. Relabelling...\n");
60 maxlabel = relabel_image(on_pixels, maxlabel, equivs, object);
61 logverb("After relabelling, we need %i labels\n", maxlabel);
62 if (maxlabel == LABEL_MAX) {
63 ERROR("Ran out of labels.");
64 exit(-1);
65 }
66 }
67 }
68
69 thislabel = object[nx*iy + ix];
70
71 /* Compute minimum equivalence label for this pixel */
72 thislabelmin = collapsing_find_minlabel(thislabel, equivs);
73
74 if (iy == 0)
75 continue;
76
77 /* Check three pixels above this one which are 'neighbours' */
78 for (i = MAX(0, ix - 1); i <= MIN(ix + 1, nx - 1); i++) {
79 if (image[nx*(iy-1)+i]) {
80 int otherlabel = object[nx*(iy-1) + i];
81
82 /* Find min of the other */
83 int otherlabelmin = collapsing_find_minlabel(otherlabel, equivs);
84
85 /* Merge groups if necessary */
86 if (thislabelmin != otherlabelmin) {
87 int oldlabelmin = MAX(thislabelmin, otherlabelmin);
88 int newlabelmin = MIN(thislabelmin, otherlabelmin);
89 thislabelmin = newlabelmin;
90 equivs[oldlabelmin] = newlabelmin;
91 equivs[thislabel] = newlabelmin;
92 /* Update other pixel too */
93 object[nx*(iy-1) + i] = newlabelmin;
94 }
95 }
96 }
97 object[nx*iy + ix] = thislabelmin;
98 }
99 }
100
101 /* Re-label the groups before returning */
102 maxlabel = relabel_image(on_pixels, maxlabel, equivs, object);
103 //logverb("After final relabelling, %i labels were used.\n", maxlabel);
104 if (pnobjects)
105 *pnobjects = maxlabel;
106
107 free(equivs);
108 il_free(on_pixels);
109 return 1;
110 }
111