/[opencvs]/eyes/analyze_green_area_adrian.c
ViewVC logotype

Contents of /eyes/analyze_green_area_adrian.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.1 - (show annotations)
Mon Oct 6 06:46:37 2014 UTC (3 years, 6 months ago) by hib
Branch: MAIN
CVS Tags: HEAD
File MIME type: text/plain
Added more stuff - just in case, the old_eyes is in /b/archive/2014
1 #include <stdio.h>
2 #include <stdlib.h>
3 #include <math.h>
4 #include <string.h>
5 #include <sys/resource.h>
6 #include <wand/MagickWand.h>
7
8
9 #include "uthash.h"
10
11 /* This thing is real cool because it
12 Lets me circle an important area in a picture, and we will make this gradient
13 from 0 to 65536 which gives the gray scale gradients that allows us to "fit"
14 the original record based on its colors and size and stuff.
15 */
16
17
18 struct point {
19 int x;
20 int y;
21 double pos_tally;
22 double pos_count;
23 int key;
24 UT_hash_handle hh;
25 };
26
27 int the_bg_color = 0xf0f0f0;
28
29
30 int bits[] = {0x80,0x40,0x20,0x10,0x8,0x4,0x2,0x1};
31
32
33
34
35
36 /* simple linear search. Sorry - this should be a hash table*/
37 int in_list(int x,int y,struct point *points_to_look_at,int number_points,struct point *hash_table) {
38 int key;
39 key=y*30000+x;
40 struct point *s;
41 int r2;
42 int i;
43 HASH_FIND_INT(hash_table,&key,s);
44 if (s) r2=1; else r2=0;
45 return r2;
46 }
47
48
49
50
51
52
53 struct point *
54 in_list2(struct point *ppoint,struct point *points_to_look_at,int number_points,struct point *hash_table) {
55 int i;
56 int x,y;
57 int key;
58 x=ppoint->x;
59 y=ppoint->y;
60 key = y*30000+x;
61 struct point *s;
62 HASH_FIND_INT(hash_table,&key,s);
63 return s;
64 }
65
66
67
68
69
70
71 void add_if_not_in_list(int x,int y,struct point *points_to_look_at,int *pnumber_points,struct point **phash_table)
72 {
73 struct point *c;
74 struct point *hash_table;
75 struct point p;
76 hash_table = *phash_table;
77 p.x=x;
78 p.y=y;
79 p.key = y*30000+x;
80 if (!in_list2(&p,points_to_look_at,*pnumber_points,hash_table)) {
81 c=points_to_look_at + *pnumber_points;
82 c->x = x;
83 c->y = y;
84 c->key = p.key;
85 c->pos_tally = 0.;
86 c->pos_count = 1.;
87 (*pnumber_points)++;
88 HASH_ADD_INT(hash_table,key,c);
89 *phash_table = hash_table;
90 }
91 }
92
93
94
95
96
97 void add_if_not_in_list_sweep_inside(int x,int y,
98 double value,struct point *points_to_look_at,int *pnumber_points,
99 int width,int height,
100 unsigned short *rbuf,
101 unsigned short *rbackbuf,
102 int cx,
103 int cy,
104 int just_changed_flag,struct point **phash_table)
105 {
106 struct point *hash_table;
107 hash_table = *phash_table;
108 if ((x<0)||(x>=width)) return;
109 if ((y<0)||(y>=height)) return;
110
111 {
112 short rcpoint;
113 short rpoint;
114 rcpoint = rbuf[cy*width+cx];
115 rpoint = rbuf[y*width+x];
116
117 rcpoint = rbackbuf[cy*width+cx];
118 rpoint = rbackbuf[y*width+x];
119
120 if (rpoint) {
121 if (rpoint<=32768) return;
122 if (!just_changed_flag) { return; }
123 }
124 }
125
126 struct point *c;
127 struct point p;
128 p.x=x;
129 p.y=y;
130 c=in_list2(&p,points_to_look_at,*pnumber_points,hash_table);
131 if (c) {
132 if (value < c->pos_tally) {
133 c->pos_tally = value;
134 c->pos_count=1.;
135 }
136 }
137 else {
138 c=points_to_look_at+ *pnumber_points;
139 c->pos_tally = value; c->pos_count=1.;
140 c->x = x;
141 c->y = y;
142 c->key = y*30000+x;
143 (*pnumber_points)++;
144 HASH_ADD_INT(hash_table,key,c);
145 *phash_table = hash_table;
146 }
147 }
148
149
150
151
152 void put_all_of_this_shade_in_list(unsigned short shade,
153 unsigned short *rbuf,int width,int height,
154 struct point *points_to_look_at,int *pnumber_points,struct point **phash_table)
155 {
156 int x;
157 int y;
158 struct point *hash_table;
159 hash_table = *phash_table;
160 for (y=0;y<height;y++) {
161 for (x=0;x<width;x++) {
162 if (rbuf[y*width+x]>=shade) {
163 points_to_look_at[*pnumber_points].x=x;
164 points_to_look_at[*pnumber_points].y=y;
165 points_to_look_at[*pnumber_points].pos_tally=rbuf[y*width+x];
166 points_to_look_at[*pnumber_points].pos_count=1;
167 {
168 struct point *c;
169 c=points_to_look_at+ *pnumber_points;
170 c->key = y*30000+x;
171 HASH_ADD_INT(hash_table,key,c);
172 }
173 (*pnumber_points)++;
174 }
175 }
176 }
177 *phash_table = hash_table;
178 }
179
180
181 /* merge sort by xy - used to quickly identify the blobs
182 we do this blob analysis to make sure that there is one blob
183 as the main section - so we don't have two parts of the same
184 picture in the end.
185 */
186 void sort_by_xy(struct point *in_points,int number_points,
187 struct point *temp_points,
188 struct point *output_points)
189 {
190 int i;
191 for (i=0;i<number_points;i++) {
192 if ((in_points[i].x==0)&&(in_points[i].y==0)) {
193 fprintf(stderr,"bad %d\n",i);
194 }
195 }
196 struct point *source;
197 struct point *dest;
198 int count=number_points;
199 source=in_points;
200 dest=temp_points;
201 if (number_points==0) return;
202 while (1) {
203 int k=0;
204 int l1=0;
205 int l2=0;
206 int l3=0;
207 for (i=0;i<number_points;i++) {
208 if ((source[i].x==0)&&(source[i].y==0)) {
209 fprintf(stderr,"b2 %d\n",i);
210 }
211 }
212 while (l1<number_points) {
213 int i,j;
214 for (i=l1+1;i<number_points;i++) {
215 if (source[i-1].x>source[i].x)
216 break;
217 if ((source[i-1].x==source[i].x)&&(source[i-1].y>source[i].y))
218 break;
219 }
220 if ((k==0)&&(i==number_points)) { /* it is sorted */
221 if (source==output_points) { /* if we are in the dest list */
222 return;
223 }
224 /* OK, then move to the output list */
225 for (i=0;i<number_points;i++) {
226 output_points[i]=source[i];
227 }
228 return;
229 }
230 l2=i;
231 if (l2==number_points) {
232 l3=l2;
233 }
234 else {
235 for (i=l2+1;i<number_points;i++) {
236 if (source[i-1].x>source[i].x)
237 break;
238 if ((source[i-1].x==source[i].x)&&(source[i-1].y>source[i].y))
239 break;
240 }
241 l3=i;
242 }
243
244 /* now that we have two lists, merge them */
245 i=l1;
246 j=l2;
247 while ((i<l2)&&(j<l3)) {
248 int winner;
249 if (source[i].x>source[j].x) {
250 winner=j;
251 j++;
252 }
253 else if (source[i].x<source[j].x) {
254 winner=i;
255 i++;
256 }
257 else if (source[i].y>source[j].y) {
258 winner=j;
259 j++;
260 }
261 else {
262 winner=i;
263 i++;
264 }
265 if ((source[winner].x==0)&&(source[winner].y==0)) {
266 fprintf(stderr,"winner bad i %d j %d\n",i,j);
267 }
268 dest[k++]=source[winner];
269 }
270 while (i<l2) {
271 if ((source[i].x==0)&&(source[i].y==0)) {
272 fprintf(stderr,"k2 winner bad i %d\n",i);
273 }
274 dest[k++]=source[i++];
275 }
276 while (j<l3) {
277 if ((source[j].x==0)&&(source[j].y==0)) {
278 fprintf(stderr,"k3 winner bad j %d\n",j);
279 }
280 dest[k++]=source[j++];
281 }
282 l1=l3; /* continue to the next two */
283 } /* while we have lists to go through */
284 fprintf(stderr,"verify k = number, %d = %d\n",k,number_points);
285 source = dest;
286 if (dest==output_points) {
287 dest = temp_points;
288 }
289 else {
290 dest=output_points;
291 }
292 if (source==dest) {
293 fprintf(stderr,"houston - problem\n");
294 }
295 } /* while we have sort work to do -- forever (return is inside loop */
296 }
297
298
299
300 /* merge sort by pos_count - used to quickly identify the blobs
301 we do this blob analysis to make sure that there is one blob
302 as the main section - so we don't have two parts of the same
303 picture in the end.
304 */
305 void sort_by_pos_count(struct point *in_points,int number_points,
306 struct point *temp_points,
307 struct point *output_points)
308 {
309 struct point *source;
310 struct point *dest;
311 int count=number_points;
312 source=in_points;
313 dest=temp_points;
314 if (number_points==0) return;
315 while (1) {
316 int k=0;
317 int l1=0;
318 int l2=0;
319 int l3=0;
320 while (l1<number_points) {
321 int i,j;
322 for (i=l1+1;i<number_points;i++) {
323 if (source[i-1].pos_count>source[i].pos_count)
324 break;
325 }
326 if ((k==0)&&(i==number_points)) { /* it is sorted */
327 if (source==output_points) { /* if we are in the dest list */
328 return;
329 }
330 /* OK, then move to the output list */
331 for (i=0;i<number_points;i++) {
332 output_points[i]=source[i];
333 }
334 return;
335 }
336 l2=i;
337 if (l2==number_points) {
338 l3=l2;
339 }
340 else {
341 for (i=l2+1;i<number_points;i++) {
342 if (source[i-1].pos_count>source[i].pos_count)
343 break;
344 }
345 l3=i;
346 }
347
348 /* now that we have two lists, merge them */
349 i=l1;
350 j=l2;
351 while ((i<l2)&&(j<l3)) {
352 int winner;
353 if (source[i].pos_count>source[j].pos_count) {
354 winner=j;
355 j++;
356 }
357 else {
358 winner=i;
359 i++;
360 }
361 dest[k++]=source[winner];
362 }
363 while (i<l2) {
364 dest[k++]=source[i++];
365 }
366 while (j<l3) {
367 dest[k++]=source[j++];
368 }
369 l1=l3; /* continue to the next two */
370 } /* while we have lists to go through */
371 source = dest;
372 if (dest==output_points) {
373 dest = temp_points;
374 }
375 else {
376 dest=output_points;
377 }
378 } /* while we have sort work to do -- forever (return is inside loop */
379 }
380
381
382 int find_shape(unsigned short *rbuf,int x, int y, int width, int height) {
383 if ((x<=0)||(x>width)||(y<=0)||(y>height)) return(0);
384 if (rbuf[y*width+x]==1) {
385 rbuf[y*width+x]=2;
386 int i;
387 int t;
388 t=1;
389 for (i=1;i<5;i++) {
390 t +=find_shape(rbuf,x+i,y,width,height);
391 t +=find_shape(rbuf,x-i,y,width,height);
392 t +=find_shape(rbuf,x,y+i,width,height);
393 t +=find_shape(rbuf,x,y-i,width,height);
394 t +=find_shape(rbuf,x+i,y+i,width,height);
395 t +=find_shape(rbuf,x-i,y+i,width,height);
396 t +=find_shape(rbuf,x-i,y-i,width,height);
397 t +=find_shape(rbuf,x-i,y-i,width,height);
398 }
399 return t;
400 }
401 return 0;
402 }
403
404
405
406
407 int fill_inside(unsigned short *rbuf,int x, int y, int width, int height) {
408 if ((x<=0)||(x>width)||(y<=0)||(y>height)) return(0);
409 if (rbuf[y*width+x]==0) {
410 rbuf[y*width+x]=32769;
411 int i;
412 int t;
413 t=1;
414 for (i=1;i<2;i++) {
415 t +=fill_inside(rbuf,x+i,y,width,height);
416 t +=fill_inside(rbuf,x-i,y,width,height);
417 t +=fill_inside(rbuf,x,y+i,width,height);
418 t +=fill_inside(rbuf,x,y-i,width,height);
419 t +=fill_inside(rbuf,x+i,y+i,width,height);
420 t +=fill_inside(rbuf,x-i,y+i,width,height);
421 t +=fill_inside(rbuf,x-i,y-i,width,height);
422 t +=fill_inside(rbuf,x-i,y-i,width,height);
423 }
424 return t;
425 }
426 return 0;
427 }
428
429
430 unsigned short get_highest_shade(unsigned short *rbuf,int width,int height) {
431 /* For this, we have to kern the end points, but just the end points
432 you see, some of these points do not have any parent, and they
433 do not all start at 1. So we get an estimate of how high they are.
434 In actuality, the later sweeping will probably go over them, but
435 you cant be too sure*/
436 unsigned short size=32769.0;
437 {
438 int cx,cy;
439 for (cy=0;cy<height;cy++) {
440 for (cx=0;cx<width;cx++) {
441 unsigned short cur=rbuf[cy*width+cx];
442 if (cur > size) size=cur;
443 }
444 }
445 }
446 return size;
447 }
448
449
450 unsigned short get_lowest_shade(unsigned short *rbuf,int width,int height) {
451 /* For this, we have to kern the end points, but just the end points
452 you see, some of these points do not have any parent, and they
453 do not all start at 1. So we get an estimate of how high they are.
454 In actuality, the later sweeping will probably go over them, but
455 you cant be too sure*/
456 unsigned short size=32767;
457 {
458 int cx,cy;
459 for (cy=0;cy<height;cy++) {
460 for (cx=0;cx<width;cx++) {
461 unsigned short cur=rbuf[cy*width+cx];
462 if (cur < size) size=cur;
463 }
464 }
465 }
466 return size;
467 }
468
469
470 int number_of_shapes_at_a_shade(unsigned short shade,unsigned short *rbuf,unsigned short *rbuf2,int width,int height,
471 struct point *a,struct point *b,struct point *c)
472 {
473 int x;
474 int y;
475 int tally=0;
476 for (y=0;y<height;y++) {
477 for (x=0;x<width;x++) {
478 if (rbuf[y*width+x]>=shade) {
479 rbuf2[y*width+x]=1;
480 }
481 else
482 rbuf2[y*width+x]=0;
483 }
484 }
485 for (y=0;y<height;y++) {
486 for (x=0;x<width;x++) {
487 int t;
488 if (t=find_shape(rbuf2,x,y,width,height)>16) { // the shape needs to be bigger than a 4x4
489 tally++;
490 if (tally>1) return tally; /* don't need to do all of it */
491 }
492 }
493 }
494 return tally;
495 }
496
497
498
499
500 void add_if_not_in_list_sweep_outside(int x,int y,
501 double value,struct point *points_to_look_at,int *pnumber_points,
502 int width,int height,
503 unsigned short *rbuf,
504 unsigned short *rbackbuf,
505 int cx,
506 int cy,
507 int just_changed_flag,struct point **phash_table)
508 {
509 if ((x<0)||(x>=width)) return;
510 if ((y<0)||(y>=height)) return;
511
512 {
513 short rcpoint;
514 short rpoint;
515 rcpoint = rbuf[cy*width+cx];
516 rpoint = rbuf[y*width+x];
517
518 rcpoint = rbackbuf[cy*width+cx];
519 rpoint = rbackbuf[y*width+x];
520
521 if (rpoint>=32768) return;
522 if (rpoint) {
523 if (!just_changed_flag) { return; }
524 }
525 }
526
527 struct point *c;
528 struct point p;
529 struct point *hash_table;
530 p.x=x;
531 p.y=y;
532 hash_table = *phash_table;
533 c=in_list2(&p,points_to_look_at,*pnumber_points,hash_table);
534 if (c) {
535 if (value < c->pos_tally) {
536 c->pos_tally = value;
537 c->pos_count=1.;
538 }
539 }
540 else {
541 c=points_to_look_at+ *pnumber_points;
542 c->pos_tally = value; c->pos_count=1.;
543 c->x = x;
544 c->y = y;
545 c->key = c->y*30000+c->x;
546 HASH_ADD_INT(hash_table,key,c);
547 *phash_table = hash_table;
548 (*pnumber_points)++;
549 }
550 }
551
552
553
554
555
556 int inside(unsigned short *rbuf,ssize_t width,ssize_t height,int the_x,int the_y) {
557 /* simple check to see if inside - count the number of lines we cross.
558 If odd then we are inside, otherwise we are outside.
559 this aint perfect because a straggler line could goof it up,
560 but it is unlikely that a straggler will goof it up in all 4 directions.
561 Unless we drew a spiral or some shit. So I will check all 4 directions to see
562 if we are inside or outside, and then vote. If the vote is a tie, assume inside I guess.
563 */
564
565 int votes[7];
566 int number_votes;
567 int x,y;
568 int onoff;
569 int count;
570
571 number_votes = 0;
572
573 int ll;
574 for (ll=1;ll<=4;ll++) {
575 int l,m;
576 if (ll==1) { l=1;m=1;}
577 else if (ll==2) { l=2;m=1;}
578 else if (ll==3) { l=1;m=2;}
579 else /* if (ll==4) */ { l=2;m=2;}
580
581 /* first see if any adjacants are pre-set to inside or outside - so we can correct mistakes for
582 those damn curly queues */
583 /* to the left */
584 x=the_x-l;
585 y=the_y;
586 if (x>=0) { // left
587 int value=rbuf[y*width+x];
588 if ((value)&&(value>32768)) return(1); // out buddy is inside, so we are too
589 if ((value)&&(value<32768)) return(0); // out buddy is outside, so we are too
590 }
591 x=the_x-l;
592 y=the_y-m;
593 if ((x>=0)&&(y>-0)) { // left up
594 int value=rbuf[y*width+x];
595 if ((value)&&(value>32768)) return(1); // out buddy is inside, so we are too
596 if ((value)&&(value<32768)) return(0); // out buddy is outside, so we are too
597 }
598 x=the_x;
599 y=the_y-m;
600 if (y>=0) { // up
601 int value=rbuf[y*width+x];
602 if ((value)&&(value>32768)) return(1); // out buddy is inside, so we are too
603 if ((value)&&(value<32768)) return(0); // out buddy is outside, so we are too
604 }
605 x=the_x+l;
606 y=the_y-m;
607 if ((x<width)&&(y>=0)) { // right up
608 int value=rbuf[y*width+x];
609 if ((value)&&(value>32768)) return(1); // out buddy is inside, so we are too
610 if ((value)&&(value<32768)) return(0); // out buddy is outside, so we are too
611 }
612 x=the_x+l;
613 y=the_y;
614 if (x<width) { // right
615 int value=rbuf[y*width+x];
616 if ((value)&&(value>32768)) return(1); // out buddy is inside, so we are too
617 if ((value)&&(value<32768)) return(0); // out buddy is outside, so we are too
618 }
619 x=the_x+l;
620 y=the_y+m;
621 if ((x<width)&&(y<height)) { // right bottom
622 int value=rbuf[y*width+x];
623 if ((value)&&(value>32768)) return(1); // out buddy is inside, so we are too
624 if ((value)&&(value<32768)) return(0); // out buddy is outside, so we are too
625 }
626 x=the_x;
627 y=the_y+m;
628 if (y<height) { // bottom
629 int value=rbuf[y*width+x];
630 if ((value)&&(value>32768)) return(1); // out buddy is inside, so we are too
631 if ((value)&&(value<32768)) return(0); // out buddy is outside, so we are too
632 }
633 x=the_x-l;
634 y=the_y+m;
635 if ((x>=0)&&(y<height)) { // left bottom
636 int value=rbuf[y*width+x];
637 if ((value)&&(value>32768)) return(1); // out buddy is inside, so we are too
638 if ((value)&&(value<32768)) return(0); // out buddy is outside, so we are too
639 }
640 }
641 /* That was the easy way, because we got someone on the "inside" see!
642 The other way is to use the old even odd algorithm, but this could make a mistake
643 because of pigs tails. So we will try in up to 4 directions, and then vote.
644 */
645 if (the_x>=0) { // left line
646 y=the_y;
647 onoff=0;count=0;
648 for (x=0;x<the_x;x++) {
649 int flag = (rbuf[y*width+x]==32768);
650 if (onoff) {
651 if (!flag) {
652 onoff = !onoff;
653 }
654 }
655 else {
656 if (flag) {
657 onoff = !onoff;
658 count++;
659 }
660 }
661 } /* for each x position */
662 votes[number_votes++] = (count & 1); /* we think it is inside if we have an odd number of intersections */
663 }
664 if (the_x < width-1) { // right line
665 y=the_y;
666 onoff=0;count=0;
667 for (x=the_x+1;x<width;x++) {
668 int flag = (rbuf[y*width+x]==32768);
669 if (onoff) {
670 if (!flag) {
671 onoff = !onoff;
672 }
673 }
674 else {
675 if (flag) {
676 onoff = !onoff;
677 count++;
678 }
679 }
680 }
681 votes[number_votes++] = (count & 1); /* we think it is inside if we have an odd number of intersections */
682 }
683 if (the_y>=0) { // above line
684 x=the_x;
685 onoff=0;count=0;
686 for (y=0;y<the_y;y++) {
687 int flag = (rbuf[y*width+x]==32768);
688 if (onoff) {
689 if (!flag) {
690 onoff = !onoff;
691 }
692 }
693 else {
694 if (flag) {
695 onoff = !onoff;
696 count++;
697 }
698 }
699 }
700 votes[number_votes++] = (count & 1); /* we think it is inside if we have an odd number of intersections */
701 }
702 if (the_y<height-1) { // below line
703 x=the_x;
704 onoff=0;count=0;
705 for (y=the_y+1;y<height;y++) {
706 int flag = (rbuf[y*width+x]==32768);
707 if (onoff) {
708 if (!flag) {
709 onoff = !onoff;
710 }
711 }
712 else {
713 if (flag) {
714 onoff = !onoff;
715 count++;
716 }
717 }
718 }
719 votes[number_votes++] = (count & 1); /* we think it is inside if we have an odd number of intersections */
720 }
721 /* now for angles */
722 if ((the_x>=0)&&(the_y>=0)) { // up left line
723 if (the_x > the_y) {
724 y = 0;
725 the_x - the_y;
726 }
727 else {
728 x = 0;
729 y = the_y - the_x;
730 }
731 onoff=0;count=0;
732 while ((x<the_x)&&(y<the_y)) {
733 int flag = (rbuf[y*width+x]==32768);
734 if (onoff) {
735 if (!flag) {
736 onoff = !onoff;
737 }
738 }
739 else {
740 if (flag) {
741 onoff = !onoff;
742 count++;
743 }
744 }
745 x++;
746 y++;
747 } /* for each x position */
748 votes[number_votes++] = (count & 1); /* we think it is inside if we have an odd number of intersections */
749 }
750 if ((the_x<width)&&(the_y>=0)) { // up right line
751 x=the_x+1;
752 y=the_y-1;
753 onoff=0;count=0;
754 while ((x<width)&&(y>0)) {
755 int flag = (rbuf[y*width+x]==32768);
756 if (onoff) {
757 if (!flag) {
758 onoff = !onoff;
759 }
760 }
761 else {
762 if (flag) {
763 onoff = !onoff;
764 count++;
765 }
766 }
767 x++;
768 y--;
769 } /* for each x position */
770 votes[number_votes++] = (count & 1); /* we think it is inside if we have an odd number of intersections */
771 }
772 if ((the_x<width-1)&&(the_y<height-1)) { //lower right line
773 x=the_x+1;
774 y=the_y+1;
775 onoff=0;count=0;
776 while ((x>0)&&(y>0)) {
777 int flag = (rbuf[y*width+x]==32768);
778 if (onoff) {
779 if (!flag) {
780 onoff = !onoff;
781 }
782 }
783 else {
784 if (flag) {
785 onoff = !onoff;
786 count++;
787 }
788 }
789 x--;
790 y--;
791 } /* for each x position */
792 votes[number_votes++] = (count & 1); /* we think it is inside if we have an odd number of intersections */
793 }
794 if (number_votes<2) {
795 fprintf(stderr,"Error - at position %d,%d we only had %d votes for directions. Should have at least 2\n",x,y,number_votes);
796 fprintf(stderr,"Clean that area up and try again\n");
797 exit(-1);
798 }
799 int inside_tally;
800 int outside_tally;
801 inside_tally=0;
802 outside_tally=0;
803 int i;
804 for (i=0;i<number_votes;i++) {
805 if (votes[i]) inside_tally++;
806 else outside_tally++;
807 }
808 if (inside_tally==outside_tally) {
809 fprintf(stderr,"lucky you, We detected a pig tail at %d,%d\n",the_x,the_y);
810 fprintf(stderr,"This condition is when a line is not closed and happens to be the first line found\n");
811 fprintf(stderr,"We look in 4 directions to be sure, but there was a tie (%d in %d out), so we cant tell\n",inside_tally,outside_tally);
812 fprintf(stderr,"Clean that area up and try again\n");
813 exit(-1);
814 }
815 if (inside_tally > outside_tally+2) { /* make sure we are more quivocal before saying yes or no */
816 return (1);
817 }
818 else {
819 return(0); /* not sure */
820 }
821 }
822
823
824
825
826
827 int is_brighter(unsigned short *rbuf,int width,int height,int cx,int cy, int distance,unsigned short myshade)
828 {
829 if (!myshade) myshade = rbuf[cy*width+cx];
830
831 if (distance >1) {
832 return (is_brighter(rbuf,width,height,cx-1,cy ,distance-1,myshade)||
833 is_brighter(rbuf,width,height,cx-1,cy-1,distance-1,myshade)||
834 is_brighter(rbuf,width,height,cx ,cy-1,distance-1,myshade)||
835 is_brighter(rbuf,width,height,cx+1,cy-1,distance-1,myshade)||
836 is_brighter(rbuf,width,height,cx+1,cy ,distance-1,myshade)||
837 is_brighter(rbuf,width,height,cx+1,cy+1,distance-1,myshade)||
838 is_brighter(rbuf,width,height,cx ,cy+1,distance-1,myshade)||
839 is_brighter(rbuf,width,height,cx-1,cy+1,distance-1,myshade));
840 }
841 int x,y;
842 int flag;
843 flag=0;
844 if (cx>0) {
845 x = cx -1;
846 if (cy >0) {
847 y=cy-1;
848 flag |= (rbuf[y*width+x]>myshade);
849 }
850 if ((cy >=0)&&(cy < height)) {
851 y=cy;
852 flag |= (rbuf[y*width+x]>myshade);
853 }
854 if (cy < height-1) {
855 y=cy+1;
856 flag |= (rbuf[y*width+x]>myshade);
857 }
858 }
859 if ((cx >=0)&&(cx <width)) {
860 x = cx;
861 if (cy >0) {
862 y=cy-1;
863 flag |= (rbuf[y*width+x]>myshade);
864 }
865 if (cy < height-1) {
866 y=cy+1;
867 flag |= (rbuf[y*width+x]>myshade);
868 }
869 }
870 if (cx < width-1) {
871 x=cx+1;
872 if (cy >0) {
873 y=cy-1;
874 flag |= (rbuf[y*width+x]>myshade);
875 }
876 if ((cy >=0)&&(cy < height)) {
877 y=cy;
878 flag |= (rbuf[y*width+x]>myshade);
879 }
880 if (cy < height-1) {
881 y=cy+1;
882 flag |= (rbuf[y*width+x]>myshade);
883 }
884 }
885 return flag;
886 }
887
888
889 int read_image(char *iname,unsigned short **prbuf,ssize_t *pwidth,ssize_t *pheight,
890 struct point **ppoints_to_look_at,int *pnumber_of_points) {
891
892 MagickPixelPacket input_pixel;
893 MagickWand *input_wand = NULL;
894 PixelWand *c_wand = NULL;
895 PixelWand *c2_wand = NULL;
896 input_wand =NewMagickWand();
897 MagickReadImage(input_wand,iname);
898 struct point *points_to_look_at;
899 ssize_t width,height;
900 width = MagickGetImageWidth(input_wand);
901 height = MagickGetImageHeight(input_wand);
902 int number_of_points;
903 number_of_points=0;
904 unsigned short *rbuf;
905
906 rbuf= (unsigned short *)malloc(width*height*sizeof(unsigned short));
907 if (!rbuf) {fprintf(stderr,"memory\n"); exit(-1);}
908 points_to_look_at= (struct point *)malloc(width*height*sizeof(struct point));
909 if (!points_to_look_at) {fprintf(stderr,"memory\n"); exit(-1);}
910
911 #define ThrowWandException(wand) \
912 { \
913 char \
914 *description; \
915 \
916 ExceptionType \
917 severity; \
918 \
919 description=MagickGetException(wand,&severity); \
920 (void) fprintf(stderr,"%s %s %lu %s\n",GetMagickModule(),description); \
921 description=(char *) MagickRelinquishMemory(description); \
922 exit(-1); \
923 }
924
925 MagickBooleanType status;
926 PixelIterator *input_iterator;
927 PixelWand **input_pixels;
928 register ssize_t x;
929 ssize_t y;
930
931 input_iterator=NewPixelIterator(input_wand);
932 if ((input_iterator == (PixelIterator *) NULL) )
933 ThrowWandException(input_wand);
934
935 for (y=0; y < (ssize_t) height; y++)
936 {
937 input_pixels=PixelGetNextIteratorRow(input_iterator,&width);
938 if ((input_pixels == (PixelWand **) NULL))
939 break;
940 for (x=0; x < (ssize_t) width; x++)
941 {
942 PixelGetMagickColor(input_pixels[x],&input_pixel);
943 // if ((input_pixel.red== (112-1)*256 )&&(input_pixel.green== (223-1)*256)&&(input_pixel.blue== (14-1)*256)) {
944 if ((input_pixel.red== 28784 )&&(input_pixel.green== 57311)&&(input_pixel.blue== 3598))
945 {
946 rbuf[y*width+x] = 32768;
947 points_to_look_at[number_of_points].x=x;
948 points_to_look_at[number_of_points].y=y;
949 points_to_look_at[number_of_points].pos_tally=0.;
950 points_to_look_at[number_of_points].pos_count=1.;
951 number_of_points++;
952 }
953 else {
954 rbuf[y*width+x] = 0;
955 }
956 }
957 }
958 if (y < (ssize_t) MagickGetImageHeight(input_wand))
959 ThrowWandException(input_wand);
960 input_iterator=DestroyPixelIterator(input_iterator);
961 input_wand=DestroyMagickWand(input_wand);
962
963
964 *prbuf = rbuf;
965 *pwidth = width;
966 *pheight=height;
967 *pnumber_of_points = number_of_points;
968 *ppoints_to_look_at = points_to_look_at;
969 return (0);
970 }
971
972
973
974 int sweep_inside(struct point *points_to_look_at,
975 struct point *second_points_to_look_at,
976 unsigned short *rbuf,
977 unsigned short *rbackbuf,
978 int j,
979 int width,
980 int height,
981 int scale)
982 {
983 struct point *ht;
984 int number_of_points;
985 int i;
986 ht = NULL;
987 while (j) {
988 HASH_CLEAR(hh,ht);
989 /* first swap second points with points*/
990 {
991 struct point *temp;
992 temp = second_points_to_look_at;
993 second_points_to_look_at = points_to_look_at;
994 points_to_look_at = temp;
995 }
996
997 number_of_points = j;
998 j=0;
999
1000
1001 /* right - now lets figure out the value for these records */
1002 for (i=0;i<number_of_points;i++) {
1003 int cx,cy;
1004 struct point *p;
1005 p=points_to_look_at+i;
1006 cx = p->x;
1007 cy = p->y;
1008 double last_point = p->pos_tally/p->pos_count;
1009 double point = rbuf[cy*width+cx];
1010 /* The value is going from 32769 to 65535, starting higher and sweeping down
1011 The point will be 65535 down also, whcih is slightly different from the outside
1012 sweep
1013 */
1014 unsigned short actual = rbackbuf[cy*width+cx];
1015 double number_of_points_to_get_here = (point-32768)+1.; // know it is not 0
1016 if (number_of_points_to_get_here<=0.0) {
1017 fprintf(stderr,"utoh point = scale\n");
1018 number_of_points_to_get_here=1.;
1019 }
1020 double next_point = last_point - (last_point - 32769.) / number_of_points_to_get_here;
1021
1022 int flag= 0;
1023 if (rbackbuf[cy*width+cx]==0) {
1024 rbackbuf[cy*width+cx] = next_point;
1025 flag=1;
1026 }
1027 else {
1028 if ((point > 32768)&&(actual < (int)(next_point))) {
1029 rbackbuf[cy*width+cx] = next_point;
1030 flag=1;
1031 }
1032 }
1033 if (flag) {
1034 // fprintf(stderr,"%d,%d last %lf actual %d me %lf npoints %lf next %lf\n",cx,cy,last_point,
1035 // actual,point,number_of_points_to_get_here,next_point);
1036 add_if_not_in_list_sweep_inside(cx-1,cy-1,next_point,second_points_to_look_at,&j,width,height,rbuf,rbackbuf,cx,cy,flag,&ht);
1037 add_if_not_in_list_sweep_inside(cx-1,cy,next_point,second_points_to_look_at,&j,width,height,rbuf,rbackbuf,cx,cy,flag,&ht);
1038 add_if_not_in_list_sweep_inside(cx-1,cy+1,next_point,second_points_to_look_at,&j,width,height,rbuf,rbackbuf,cx,cy,flag,&ht);
1039 add_if_not_in_list_sweep_inside(cx,cy+1,next_point,second_points_to_look_at,&j,width,height,rbuf,rbackbuf,cx,cy,flag,&ht);
1040 add_if_not_in_list_sweep_inside(cx+1,cy+1,next_point,second_points_to_look_at,&j,width,height,rbuf,rbackbuf,cx,cy,flag,&ht);
1041 add_if_not_in_list_sweep_inside(cx+1,cy,next_point,second_points_to_look_at,&j,width,height,rbuf,rbackbuf,cx,cy,flag,&ht);
1042 add_if_not_in_list_sweep_inside(cx+1,cy-1,next_point,second_points_to_look_at,&j,width,height,rbuf,rbackbuf,cx,cy,flag,&ht);
1043 add_if_not_in_list_sweep_inside(cx,cy-1,next_point,second_points_to_look_at,&j,width,height,rbuf,rbackbuf,cx,cy,flag,&ht);
1044 }
1045 } /* for all points */
1046
1047
1048 } /* while we have new points to look at */
1049
1050 }
1051
1052
1053
1054
1055
1056
1057
1058
1059
1060 int sweep_outside(struct point *points_to_look_at,
1061 struct point *second_points_to_look_at,
1062 unsigned short *rbuf,
1063 unsigned short *rbackbuf,
1064 int j,
1065 int width,
1066 int height)
1067 {
1068 int number_of_points;
1069 int i;
1070 struct point *ht;
1071 ht = NULL;
1072
1073 int cnt;
1074 cnt=0;
1075 while (j) {
1076 HASH_CLEAR(hh,ht);
1077
1078 /* first swap second points with points*/
1079 {
1080 struct point *temp;
1081 temp = second_points_to_look_at;
1082 second_points_to_look_at = points_to_look_at;
1083 points_to_look_at = temp;
1084 }
1085
1086 number_of_points = j;
1087 j=0;
1088
1089
1090 /* right - now lets figure out the value for these records */
1091 for (i=0;i<number_of_points;i++) {
1092 int cx,cy;
1093 struct point *p;
1094 p=points_to_look_at+i;
1095 cx = p->x;
1096 cy = p->y;
1097 double last_point = p->pos_tally/p->pos_count;
1098 double point = rbuf[cy*width+cx];
1099 /* The value is going from 0 to 32768, not including 32768
1100 The point will be 32768 - number of points to get here
1101 So we */
1102 unsigned short actual = rbackbuf[cy*width+cx];
1103 double number_of_points_to_get_here = 32768 - point; // know it is not 0
1104 double next_point = last_point + (32768.-last_point) / number_of_points_to_get_here;
1105
1106 int flag= 0;
1107 if (rbackbuf[cy*width+cx]==0) {
1108 rbackbuf[cy*width+cx] = next_point;
1109 flag=1;
1110 }
1111 else {
1112 if ((point < 32768)&&(actual > (int)(next_point))) {
1113 rbackbuf[cy*width+cx] = next_point;
1114 flag=1;
1115 }
1116 }
1117 if (flag) {
1118 // fprintf(stderr,"%d,%d last %lf actual %d me %lf npoints %lf next %lf\n",cx,cy,last_point,
1119 // actual,point,number_of_points_to_get_here,next_point);
1120 // add_if_not_in_list_sweep_outside(cx-1,cy-1,next_point,second_points_to_look_at,&j,width,height,rbuf,rbackbuf,cx,cy,flag,&ht);
1121 add_if_not_in_list_sweep_outside(cx-1,cy,next_point,second_points_to_look_at,&j,width,height,rbuf,rbackbuf,cx,cy,flag,&ht);
1122 // add_if_not_in_list_sweep_outside(cx-1,cy+1,next_point,second_points_to_look_at,&j,width,height,rbuf,rbackbuf,cx,cy,flag,&ht);
1123 add_if_not_in_list_sweep_outside(cx,cy+1,next_point,second_points_to_look_at,&j,width,height,rbuf,rbackbuf,cx,cy,flag,&ht);
1124 // add_if_not_in_list_sweep_outside(cx+1,cy+1,next_point,second_points_to_look_at,&j,width,height,rbuf,rbackbuf,cx,cy,flag,&ht);
1125 add_if_not_in_list_sweep_outside(cx+1,cy,next_point,second_points_to_look_at,&j,width,height,rbuf,rbackbuf,cx,cy,flag,&ht);
1126 // add_if_not_in_list_sweep_outside(cx+1,cy-1,next_point,second_points_to_look_at,&j,width,height,rbuf,rbackbuf,cx,cy,flag,&ht);
1127 add_if_not_in_list_sweep_outside(cx,cy-1,next_point,second_points_to_look_at,&j,width,height,rbuf,rbackbuf,cx,cy,flag,&ht);
1128 }
1129 } /* for all points */
1130 cnt++;
1131 // if (cnt>600) break; - this was used for debugging
1132
1133 } /* while we have new points to look at */
1134
1135 }
1136
1137
1138
1139
1140
1141
1142
1143
1144 int main(int argc,char **argv)
1145 {
1146 /* set the stack limit to 500 megabytes */
1147 const rlim_t stacksize = 2000*1024*1024; // 2 gigabytes
1148 struct rlimit rl;
1149 int result;
1150 result = getrlimit(RLIMIT_STACK, &rl);
1151 if (result == 0)
1152 {
1153 if (rl.rlim_cur < stacksize)
1154 {
1155 rl.rlim_cur = stacksize;
1156 result = setrlimit(RLIMIT_STACK, &rl);
1157 if (result != 0)
1158 {
1159 fprintf(stderr, "setrlimit returned result = %d\n", result);
1160 exit(-1);
1161 }
1162 }
1163 }
1164 else {
1165 fprintf(stderr, "getrlimit returned result = %d\n", result);
1166 exit(-1);
1167 }
1168
1169 /* seed1 ip_Address fname */
1170 char *iname=argv[1];
1171 char *oname=argv[2];
1172
1173 unsigned char *otherrand;
1174 int other;
1175 unsigned short *rbuf;
1176 unsigned short *rbackbuf;
1177 int number_of_points;
1178 int number_second_points;
1179 struct point *points_to_look_at;
1180 struct point *second_points_to_look_at;
1181 struct point *end_points_to_look_at;
1182 struct point *temp_points;
1183 int number_end_points;
1184 MagickWand *m_wand = NULL;
1185 PixelWand *c_wand = NULL;
1186 PixelWand *c2_wand = NULL;
1187 MagickPixelPacket pixel;
1188
1189 MagickWandGenesis();
1190 ssize_t width,height;
1191
1192
1193
1194 read_image(iname,&rbuf,&width,&height,&points_to_look_at,&number_of_points);
1195 rbackbuf = malloc(width*height*2*sizeof(unsigned short));
1196 if (rbackbuf == NULL) {fprintf(stderr,"memory\n");exit(-1);}
1197 second_points_to_look_at = malloc(width*height*sizeof(struct point));
1198 if (second_points_to_look_at == NULL) {fprintf(stderr,"memory\n");exit(-1);}
1199 end_points_to_look_at = malloc(width*height*sizeof(struct point));
1200 if (end_points_to_look_at == NULL) {fprintf(stderr,"memory\n");exit(-1);}
1201 temp_points = malloc(width*height*sizeof(struct point));
1202 if (temp_points == NULL) {fprintf(stderr,"memory\n");exit(-1);}
1203
1204
1205
1206 /*
1207 So we just find the innies and outies. The innies will have 32769 and the outies will have 32767
1208 One layer thick around all of everything.
1209 */
1210 /* OK we find an inside point. Once found, we set everytihng to 1 */
1211 {
1212 int cnt;
1213 cnt = 0;
1214 int x,y;
1215 while (cnt<100000) {
1216 x = rand()%width;
1217 y = rand()%height;
1218 if (rbuf[y*width+x] != 32768) {
1219 if (inside(rbuf,width,height,x,y)) {
1220 break;
1221 }
1222 }
1223 cnt++;
1224 }
1225 if (cnt>=100000) {
1226 fprintf(stderr,"error %s missing green\n",argv[1]);
1227 exit(-1);
1228 }
1229 fill_inside(rbuf,x,y,width,height);
1230 }
1231
1232 /* check to see that we have some outside */
1233 {int cx,cy;
1234 for (cx=0;cx<width;cx++) {
1235 for (cy=0;cy<height;cy++) {
1236 if (rbuf[cy*width+cx]==0) break;
1237 }
1238 if (cy<height) break;
1239 }
1240 if ((cx>=width)&&(cy>=height)) {
1241 fprintf(stderr,"error %s does not close up\n",argv[1]);
1242 exit(-1);
1243 }
1244 }
1245
1246
1247 #define P(x,y) (((y)<height)&&((x)<width)&&((y)>=0)&&((x)>=0)&&(rbuf[(y)*width+(x)]==32768))
1248 /* prune inside to within one of the green line and also add the outside border */
1249 int j;
1250 struct point *ht;
1251 ht = NULL;
1252 struct point *endpoint_ht;
1253 endpoint_ht = NULL;
1254 j=0;
1255 {int cx,cy;
1256 for (cx=0;cx<width;cx++) {
1257 for (cy=0;cy<height;cy++) {
1258 if (rbuf[cy*width+cx]==32769) {
1259 if (P(cx-1,cy-1)||P(cx-1,cy)||P(cx-1,cy+1)||P(cx,cy+1)||P(cx+1,cy+1)||P(cx+1,cy)||P(cx+1,cy-1)||P(cx,cy-1)) {
1260 add_if_not_in_list(cx,cy,second_points_to_look_at,&j,&ht);
1261 }
1262 else {
1263 rbuf[cy*width+cx]=0;
1264 }
1265 }
1266 else if (rbuf[cy*width+cx]==0) { /* outside */
1267 if (P(cx-1,cy-1)||P(cx-1,cy)||P(cx-1,cy+1)||P(cx,cy+1)||P(cx+1,cy+1)||P(cx+1,cy)||P(cx+1,cy-1)||P(cx,cy-1)) {
1268 rbuf[cy*width+cx]=32767; /* fix it */
1269 }
1270 }
1271 }
1272 }
1273 }
1274
1275
1276
1277 /* handle inside */
1278 number_end_points = 0;
1279 while (j) {
1280 /* first swap second points with points*/
1281 {
1282 struct point *temp;
1283 temp = second_points_to_look_at;
1284 second_points_to_look_at = points_to_look_at;
1285 points_to_look_at = temp;
1286 }
1287 number_of_points = j;
1288
1289
1290 HASH_CLEAR(hh,ht);
1291 j=0;
1292 while (number_of_points) {
1293 int i=(rand()>>1) % number_of_points;
1294 int cx=points_to_look_at[i].x;
1295 int cy=points_to_look_at[i].y;
1296 number_of_points--;
1297 points_to_look_at[i] = points_to_look_at[number_of_points];
1298 int x,y;
1299 int oldj=j; // used to check if we are an end point
1300 x=cx-1;
1301 y=cy;
1302 if (x>=0) {
1303 if (rbuf[y*width+x]==0) { // do we not have a value
1304 rbuf[y*width+x] = rbuf[cy*width+cx]+ 1; /* go up!*/
1305 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1306 } /* if we are ripe for marking as inside or outside */
1307 } /* if we should check this pixel */
1308 x=cx-1;
1309 y=cy-1;
1310 if ((x>=0)&&(y>=0)) {
1311 if (rbuf[y*width+x]==0) { // do we not have a value
1312 rbuf[y*width+x] = rbuf[cy*width+cx] + 1; /* go up! But diagonally */
1313 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1314 } /* if we are ripe for marking as inside or outside */
1315 } /* if we should check this pixel */
1316 x=cx;
1317 y=cy-1;
1318 if (y>=0) {
1319 if (rbuf[y*width+x]==0) { // do we not have a value
1320 rbuf[y*width+x] = rbuf[cy*width+cx]+ 1; /* go up!*/
1321 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1322 } /* if we are ripe for marking as inside or outside */
1323 } /* if we should check this pixel */
1324 x=cx+1;
1325 y=cy-1;
1326 if ((x<width)&&(y>=0)) {
1327 if (rbuf[y*width+x]==0) { // do we not have a value
1328 rbuf[y*width+x] = rbuf[cy*width+cx] + 1; /* go up! But diagonally */
1329 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1330 } /* if we are ripe for marking as inside or outside */
1331 } /* if we should check this pixel */
1332 x=cx+1;
1333 y=cy;
1334 if (x<width) {
1335 if (rbuf[y*width+x]==0) { // do we not have a value
1336 rbuf[y*width+x] = rbuf[cy*width+cx]+ 1; /* go up!*/
1337 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1338 } /* if we are ripe for marking as inside or outside */
1339 } /* if we should check this pixel */
1340 x=cx+1;
1341 y=cy+1;
1342 if ((x<width)&&(y<height)) {
1343 if (rbuf[y*width+x]==0) { // do we not have a value
1344 rbuf[y*width+x] = rbuf[cy*width+cx] + 1; /* go up! But diagonally */
1345 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1346 } /* if we are ripe for marking as inside or outside */
1347 } /* if we should check this pixel */
1348 x=cx;
1349 y=cy+1;
1350 if (y<height) {
1351 if (rbuf[y*width+x]==0) { // do we not have a value
1352 rbuf[y*width+x] = rbuf[cy*width+cx]+ 1; /* go up!*/
1353 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1354 } /* if we are ripe for marking as inside or outside */
1355 } /* if we should check this pixel */
1356 x=cx-1;
1357 y=cy+1;
1358 if ((x>=0)&&(y<height)) {
1359 if (rbuf[y*width+x]==0) { // do we not have a value
1360 rbuf[y*width+x] = rbuf[cy*width+cx] + 1; /* go up! But diagonally */
1361 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1362 } /* if we are ripe for marking as inside or outside */
1363 } /* if we should check this pixel */
1364 if (j==oldj) { // we are an end point
1365 /* well, we might be, or we might be snuffed from another
1366 point, in which case, some of our neighbors will be in the list */
1367 if (!in_list(cx+1,cy ,second_points_to_look_at,j,ht) &&
1368 !in_list(cx+1,cy+1,second_points_to_look_at,j,ht) &&
1369 !in_list(cx ,cy+1,second_points_to_look_at,j,ht) &&
1370 !in_list(cx-1,cy+1,second_points_to_look_at,j,ht) &&
1371 !in_list(cx-1,cy ,second_points_to_look_at,j,ht) &&
1372 !in_list(cx-1,cy-1,second_points_to_look_at,j,ht) &&
1373 !in_list(cx ,cy-1,second_points_to_look_at,j,ht) &&
1374 !in_list(cx+1,cy-1,second_points_to_look_at,j,ht)) {
1375 end_points_to_look_at[number_end_points].x=cx;
1376 end_points_to_look_at[number_end_points].y=cy;
1377 end_points_to_look_at[number_end_points].pos_count=1; /* we will do the tally later */
1378 {
1379 struct point *c;
1380 c=end_points_to_look_at+number_end_points;
1381 c->key = c->y*30000+c->x;
1382 HASH_ADD_INT(endpoint_ht,key,c);
1383 number_end_points++;
1384 }
1385 }
1386 }
1387 } /* while we have points from this round to do */
1388 } /* while we have points to do */
1389
1390 HASH_CLEAR(hh,ht);
1391
1392
1393
1394 /* OK - well, I don't trust the end points yet - its still a bit stuffy, so I will
1395 see if there are any neightbors that are higher */
1396 int i;
1397 i=0;
1398 while (i<number_end_points) {
1399 int cx = end_points_to_look_at[i].x;
1400 int cy = end_points_to_look_at[i].y;
1401
1402 if (is_brighter(rbuf,width,height,cx,cy,1,0)) { /* 5-8 is good. 9 is too slow */
1403 number_end_points--;
1404 end_points_to_look_at[i]= end_points_to_look_at[number_end_points];
1405 }
1406 else {
1407 i++;
1408 }
1409 }
1410
1411
1412 /* now we need to do the outside - finding the end points is a bit easier, but we don't need them anyways */
1413 /* to start, collect the outside points - which are value 32767 */
1414 j=0;
1415 {
1416 int cx,cy;
1417 for (cy=0;cy<height;cy++) {
1418 for (cx=0;cx<width;cx++) {
1419 unsigned short cur=rbuf[cy*width+cx];
1420 if (cur ==32767) {
1421 second_points_to_look_at[j].x = cx;
1422 second_points_to_look_at[j].y = cy;
1423 j++;
1424 }
1425 }
1426 }
1427 }
1428
1429
1430 int level;
1431 level=32767;
1432 while (j) {
1433 level--;
1434 /* first swap second points with points*/
1435 {
1436 struct point *temp;
1437 temp = second_points_to_look_at;
1438 second_points_to_look_at = points_to_look_at;
1439 points_to_look_at = temp;
1440 }
1441 number_of_points = j;
1442
1443 HASH_CLEAR(hh,ht);
1444 j=0;
1445 while (number_of_points) {
1446 int i=(rand()>>1) % number_of_points;
1447 int cx=points_to_look_at[i].x;
1448 int cy=points_to_look_at[i].y;
1449 number_of_points--;
1450 points_to_look_at[i] = points_to_look_at[number_of_points];
1451 int x,y;
1452 int oldj=j; // used to check if we are an end point
1453 x=cx-1;
1454 y=cy;
1455 if (x>=0) {
1456 if (rbuf[y*width+x]==0) { // do we not have a value
1457 rbuf[y*width+x] = level;
1458 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1459 } /* if we are ripe for marking as inside or outside */
1460 } /* if we should check this pixel */
1461 x=cx-1;
1462 y=cy-1;
1463 if ((x>=0)&&(y>=0)) {
1464 if (rbuf[y*width+x]==0) { // do we not have a value
1465 rbuf[y*width+x] = level;
1466 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1467 } /* if we are ripe for marking as inside or outside */
1468 } /* if we should check this pixel */
1469 x=cx;
1470 y=cy-1;
1471 if (y>=0) {
1472 if (rbuf[y*width+x]==0) { // do we not have a value
1473 rbuf[y*width+x] = level;
1474 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1475 } /* if we are ripe for marking as inside or outside */
1476 } /* if we should check this pixel */
1477 x=cx+1;
1478 y=cy-1;
1479 if ((x<width)&&(y>=0)) {
1480 if (rbuf[y*width+x]==0) { // do we not have a value
1481 rbuf[y*width+x] = level;
1482 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1483 } /* if we are ripe for marking as inside or outside */
1484 } /* if we should check this pixel */
1485 x=cx+1;
1486 y=cy;
1487 if (x<width) {
1488 if (rbuf[y*width+x]==0) { // do we not have a value
1489 rbuf[y*width+x] = level;
1490 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1491 } /* if we are ripe for marking as inside or outside */
1492 } /* if we should check this pixel */
1493 x=cx+1;
1494 y=cy+1;
1495 if ((x<width)&&(y<height)) {
1496 if (rbuf[y*width+x]==0) { // do we not have a value
1497 rbuf[y*width+x] = level;
1498 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1499 } /* if we are ripe for marking as inside or outside */
1500 } /* if we should check this pixel */
1501 x=cx;
1502 y=cy+1;
1503 if (y<height) {
1504 if (rbuf[y*width+x]==0) { // do we not have a value
1505 rbuf[y*width+x] = level;
1506 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1507 } /* if we are ripe for marking as inside or outside */
1508 } /* if we should check this pixel */
1509 x=cx-1;
1510 y=cy+1;
1511 if ((x>=0)&&(y<height)) {
1512 if (rbuf[y*width+x]==0) { // do we not have a value
1513 rbuf[y*width+x] = level;
1514 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1515 } /* if we are ripe for marking as inside or outside */
1516 } /* if we should check this pixel */
1517 } /* while we have points from this round to do */
1518 } /* while we have points to do */
1519
1520
1521
1522
1523
1524
1525
1526
1527 /*OK - lets plateau the inside so that we always have one object
1528 note - we use the back buffer as a scratch space*/
1529 {
1530 unsigned short max_shade_point = get_highest_shade(rbuf,width,height);
1531 unsigned short min_shade_point = 32769;
1532 unsigned short shade_point=32769;
1533 int movement=(max_shade_point-min_shade_point-1)/20;
1534 if (movement<1) movement=1;
1535 int tally;
1536 tally=number_of_shapes_at_a_shade(shade_point,rbuf,rbackbuf,width,height,points_to_look_at,second_points_to_look_at,temp_points);
1537 if (tally == 1) {
1538
1539 while (1) { // While moving throght the sample
1540 shade_point=min_shade_point;
1541 tally=1;
1542 while (1) {
1543 tally=number_of_shapes_at_a_shade(shade_point+movement,rbuf,rbackbuf,width,height,points_to_look_at,second_points_to_look_at,temp_points);
1544 // fprintf(stderr,"plateau from %d to %d inc %d sp %d test %d result %d\n",min_shade_point,max_shade_point,movement,shade_point,shade_point+movement,tally);
1545 if (tally !=1) break;
1546 shade_point += movement;
1547 if (shade_point == max_shade_point) break;
1548 if (shade_point+movement > max_shade_point) shade_point = max_shade_point-movement;
1549 }
1550 if (tally==1) break; // all good
1551 min_shade_point = shade_point;
1552 max_shade_point = shade_point+movement-1;
1553 if (min_shade_point == max_shade_point) break; // also all good
1554 movement=(max_shade_point-min_shade_point-1)/20;
1555 if (movement<1) movement=1;
1556 // fprintf(stderr,"plateau new min %d max %d movement %d\n",min_shade_point,max_shade_point,movement);
1557 }
1558 /* clip it */
1559 {
1560 int cx,cy;
1561 for (cy=0;cy<height;cy++) {
1562 for (cx=0;cx<width;cx++) {
1563 if (rbuf[cy*width+cx] > shade_point) {
1564 rbuf[cy*width+cx] = shade_point;
1565 }
1566 }
1567 }
1568 }
1569 } /* if we have one shape at the first level */
1570 else {
1571 fprintf(stderr,"for %s no clipping - lowest level does not work\n",argv[1]);
1572 }
1573 }
1574
1575
1576
1577
1578
1579 /* copy to the back buffer*/
1580 {
1581 int cx,cy;
1582 for (cy=0;cy<height;cy++) {
1583 for (cx=0;cx<width;cx++) {
1584 double cur= rbuf[cy*width+cx];
1585 if (cur==32768) {
1586 rbackbuf[cy*width+cx] = cur;
1587 }
1588 else {
1589 rbackbuf[cy*width+cx] = 0.;
1590 }
1591 }
1592 }
1593 }
1594
1595
1596
1597
1598
1599 /* OK - now we need to add the end points - which are always at the edge of the screen
1600 We set the values based on the fact that anything after should be zero, and the number
1601 we have lets us know how many gradients we are from 0
1602 First build the list
1603 */
1604 j=0;
1605 for (i=0;i<height;i++) {
1606 if (rbuf[i*width+0] < 32768) {
1607 second_points_to_look_at[j].x = 0;
1608 second_points_to_look_at[j].y = i;
1609 second_points_to_look_at[j].pos_tally = 0.;
1610 second_points_to_look_at[j].pos_count = 1.;
1611 j++;
1612 }
1613 if (rbuf[i*width+(width-1)] < 32768) {
1614 second_points_to_look_at[j].x = width-1;
1615 second_points_to_look_at[j].y = i;
1616 second_points_to_look_at[j].pos_tally = 0.;
1617 second_points_to_look_at[j].pos_count = 1.;
1618 j++;
1619 }
1620 }
1621 for (i=0;i<width;i++) {
1622 if (rbuf[0*width+i] < 32768) {
1623 second_points_to_look_at[j].x = i;
1624 second_points_to_look_at[j].y = 0;
1625 second_points_to_look_at[j].pos_tally = 0.;
1626 second_points_to_look_at[j].pos_count = 1.;
1627 j++;
1628 }
1629 if (rbuf[(height-1)*width+i] < 32768) {
1630 second_points_to_look_at[j].x = i;
1631 second_points_to_look_at[j].y = height-1;
1632 second_points_to_look_at[j].pos_tally = 0.;
1633 second_points_to_look_at[j].pos_count = 1.;
1634 j++;
1635 }
1636 }
1637
1638 sweep_outside(points_to_look_at,second_points_to_look_at,rbuf,rbackbuf,j,width,height);
1639
1640
1641
1642
1643 /* clean up the inside
1644 You see, the outside all goes to zero, at the edge, well, for the most part anyways
1645 But the inside could have these multiple points. We can't deal with that
1646 because then the picture would have multiple points showing - and we need it always to be one
1647 solid single blob
1648 */
1649 #ifdef waitfdf
1650 unsigned short max=32769.0;
1651 {
1652 int cx,cy;
1653 for (cy=0;cy<height;cy++) {
1654 for (cx=0;cx<width;cx++) {
1655 unsigned short cur=rbuf[cy*width+cx];
1656 if (cur > max) max=cur;
1657 }
1658 }
1659 }
1660 /* OK so max is now equal to the highest point. Lets see how many points we have
1661 like that
1662 */
1663 j=0;
1664 {
1665 int cx,cy;
1666 for (cy=0;cy<height;cy++) {
1667 for (cx=0;cx<width;cx++) {
1668 unsigned short cur=rbuf[cy*width+cx];
1669 if (cur > max) max=cur;
1670 }
1671 }
1672 }
1673 #endif
1674
1675
1676
1677
1678 /* lets do the inside now as well */
1679
1680
1681 unsigned short size = get_highest_shade(rbuf,width,height);
1682 double scale = (size-32768);
1683 double full_over_scale = 32767.0 / scale;
1684
1685
1686
1687
1688
1689
1690 /* Cool, now we have to "kern" the inside values of rbuf from the start point. */
1691
1692 /* so this part was done above - all we have to do is apply scale and full_over_scale
1693 The code is still here but commented so you can see it logically */
1694 #ifdef done_Earlier
1695 unsigned short size=get_highest_shade(rbuf,width,height);
1696 double scale = (size-32768);
1697 double full_over_scale = 32767.0 / scale;
1698 #endif
1699 {
1700 int cx,cy;
1701 for (cy=0;cy<height;cy++) {
1702 for (cx=0;cx<width;cx++) {
1703 double cur= rbuf[cy*width+cx];
1704 if (cur>32768) {
1705 cur = ((cur -32768.) * full_over_scale) + 32768.;
1706 rbuf[cy*width+cx] = cur;
1707 rbackbuf[cy*width+cx] = cur; /* because we are not sweeping inside */
1708 }
1709 }
1710 }
1711 }
1712
1713
1714
1715 /* Cool, now we have to "kern" the outside values from the start point. */
1716 size=get_lowest_shade(rbuf,width,height);
1717 scale = 32768-(size);
1718 full_over_scale = 32767.0 / scale;
1719
1720
1721 {
1722 int cx,cy;
1723 for (cy=0;cy<height;cy++) {
1724 for (cx=0;cx<width;cx++) {
1725 double cur= rbuf[cy*width+cx];
1726 if (cur<32768) {
1727 cur = 32768 - (32768 - cur) * full_over_scale;
1728 rbuf[cy*width+cx] = cur;
1729 }
1730 }
1731 }
1732 }
1733
1734
1735
1736
1737 PixelIterator *iterator;
1738 PixelWand **pixels;
1739 m_wand = NewMagickWand();
1740 c_wand = NewPixelWand();
1741 PixelSetColor(c_wand,"black"); /* for now */
1742 MagickNewImage(m_wand,width,height,c_wand);
1743 MagickSetImageChannelDepth(m_wand, AllChannels, 16);
1744 MagickSetImageDepth(m_wand, 16);
1745
1746 iterator=NewPixelIterator(m_wand);
1747 if ((iterator == (PixelIterator *) NULL) )
1748 ThrowWandException(m_wand);
1749
1750
1751 register ssize_t x;
1752 ssize_t y;
1753 for (y=0; y < (ssize_t) height; y++)
1754 {
1755 pixels=PixelGetNextIteratorRow(iterator,&width);
1756 if ((pixels == (PixelWand **) NULL))
1757 break;
1758 for (x=0; x < (ssize_t) width; x++)
1759 {
1760 PixelGetMagickColor(pixels[x],&pixel);
1761 pixel.red=pixel.green=pixel.blue=rbackbuf[y*width+x];
1762 PixelSetMagickColor(pixels[x],&pixel);
1763 }
1764 (void) PixelSyncIterator(iterator);
1765 }
1766 if (y < (ssize_t) MagickGetImageHeight(m_wand))
1767 ThrowWandException(m_wand);
1768 iterator=DestroyPixelIterator(iterator);
1769
1770
1771 MagickGaussianBlurImage(m_wand,20.,2.);
1772
1773 /*
1774 Write the image then destroy it.
1775 */
1776 MagickBooleanType status;
1777 status=MagickWriteImages(m_wand,oname,MagickTrue);
1778 if (status == MagickFalse)
1779 ThrowWandException(m_wand);
1780 m_wand=DestroyMagickWand(m_wand);
1781 MagickWandTerminus();
1782 return(0);
1783 }
1784

  ViewVC Help
Powered by ViewVC 1.1.5