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

Contents of /eyes/analyze_green_area.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.1 - (show annotations)
Thu Mar 8 15:49:21 2012 UTC (6 years, 7 months ago) by hib
Branch: MAIN
CVS Tags: HEAD
File MIME type: text/plain
more cleanup for girl3 (Jenni) and cyclone is still in process
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==0)&&(input_pixel.green>=65536-256)&&(input_pixel.blue==0)) {
944 rbuf[y*width+x] = 32768;
945 points_to_look_at[number_of_points].x=x;
946 points_to_look_at[number_of_points].y=y;
947 points_to_look_at[number_of_points].pos_tally=0.;
948 points_to_look_at[number_of_points].pos_count=1.;
949 number_of_points++;
950 }
951 else {
952 rbuf[y*width+x] = 0;
953 }
954 }
955 }
956 if (y < (ssize_t) MagickGetImageHeight(input_wand))
957 ThrowWandException(input_wand);
958 input_iterator=DestroyPixelIterator(input_iterator);
959 input_wand=DestroyMagickWand(input_wand);
960
961
962 *prbuf = rbuf;
963 *pwidth = width;
964 *pheight=height;
965 *pnumber_of_points = number_of_points;
966 *ppoints_to_look_at = points_to_look_at;
967 return (0);
968 }
969
970
971
972 int sweep_inside(struct point *points_to_look_at,
973 struct point *second_points_to_look_at,
974 unsigned short *rbuf,
975 unsigned short *rbackbuf,
976 int j,
977 int width,
978 int height,
979 int scale)
980 {
981 struct point *ht;
982 int number_of_points;
983 int i;
984 ht = NULL;
985 while (j) {
986 HASH_CLEAR(hh,ht);
987 /* first swap second points with points*/
988 {
989 struct point *temp;
990 temp = second_points_to_look_at;
991 second_points_to_look_at = points_to_look_at;
992 points_to_look_at = temp;
993 }
994
995 number_of_points = j;
996 j=0;
997
998
999 /* right - now lets figure out the value for these records */
1000 for (i=0;i<number_of_points;i++) {
1001 int cx,cy;
1002 struct point *p;
1003 p=points_to_look_at+i;
1004 cx = p->x;
1005 cy = p->y;
1006 double last_point = p->pos_tally/p->pos_count;
1007 double point = rbuf[cy*width+cx];
1008 /* The value is going from 32769 to 65535, starting higher and sweeping down
1009 The point will be 65535 down also, whcih is slightly different from the outside
1010 sweep
1011 */
1012 unsigned short actual = rbackbuf[cy*width+cx];
1013 double number_of_points_to_get_here = (point-32768)+1.; // know it is not 0
1014 if (number_of_points_to_get_here<=0.0) {
1015 fprintf(stderr,"utoh point = scale\n");
1016 number_of_points_to_get_here=1.;
1017 }
1018 double next_point = last_point - (last_point - 32769.) / number_of_points_to_get_here;
1019
1020 int flag= 0;
1021 if (rbackbuf[cy*width+cx]==0) {
1022 rbackbuf[cy*width+cx] = next_point;
1023 flag=1;
1024 }
1025 else {
1026 if ((point > 32768)&&(actual < (int)(next_point))) {
1027 rbackbuf[cy*width+cx] = next_point;
1028 flag=1;
1029 }
1030 }
1031 if (flag) {
1032 // fprintf(stderr,"%d,%d last %lf actual %d me %lf npoints %lf next %lf\n",cx,cy,last_point,
1033 // actual,point,number_of_points_to_get_here,next_point);
1034 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);
1035 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);
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,cy+1,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+1,cy,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,cy-1,next_point,second_points_to_look_at,&j,width,height,rbuf,rbackbuf,cx,cy,flag,&ht);
1042 }
1043 } /* for all points */
1044
1045
1046 } /* while we have new points to look at */
1047
1048 }
1049
1050
1051
1052
1053
1054
1055
1056
1057
1058 int sweep_outside(struct point *points_to_look_at,
1059 struct point *second_points_to_look_at,
1060 unsigned short *rbuf,
1061 unsigned short *rbackbuf,
1062 int j,
1063 int width,
1064 int height)
1065 {
1066 int number_of_points;
1067 int i;
1068 struct point *ht;
1069 ht = NULL;
1070
1071 int cnt;
1072 cnt=0;
1073 while (j) {
1074 HASH_CLEAR(hh,ht);
1075
1076 /* first swap second points with points*/
1077 {
1078 struct point *temp;
1079 temp = second_points_to_look_at;
1080 second_points_to_look_at = points_to_look_at;
1081 points_to_look_at = temp;
1082 }
1083
1084 number_of_points = j;
1085 j=0;
1086
1087
1088 /* right - now lets figure out the value for these records */
1089 for (i=0;i<number_of_points;i++) {
1090 int cx,cy;
1091 struct point *p;
1092 p=points_to_look_at+i;
1093 cx = p->x;
1094 cy = p->y;
1095 double last_point = p->pos_tally/p->pos_count;
1096 double point = rbuf[cy*width+cx];
1097 /* The value is going from 0 to 32768, not including 32768
1098 The point will be 32768 - number of points to get here
1099 So we */
1100 unsigned short actual = rbackbuf[cy*width+cx];
1101 double number_of_points_to_get_here = 32768 - point; // know it is not 0
1102 double next_point = last_point + (32768.-last_point) / number_of_points_to_get_here;
1103
1104 int flag= 0;
1105 if (rbackbuf[cy*width+cx]==0) {
1106 rbackbuf[cy*width+cx] = next_point;
1107 flag=1;
1108 }
1109 else {
1110 if ((point < 32768)&&(actual > (int)(next_point))) {
1111 rbackbuf[cy*width+cx] = next_point;
1112 flag=1;
1113 }
1114 }
1115 if (flag) {
1116 // fprintf(stderr,"%d,%d last %lf actual %d me %lf npoints %lf next %lf\n",cx,cy,last_point,
1117 // actual,point,number_of_points_to_get_here,next_point);
1118 // 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);
1119 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);
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,cy+1,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+1,cy,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,cy-1,next_point,second_points_to_look_at,&j,width,height,rbuf,rbackbuf,cx,cy,flag,&ht);
1126 }
1127 } /* for all points */
1128 cnt++;
1129 // if (cnt>600) break; - this was used for debugging
1130
1131 } /* while we have new points to look at */
1132
1133 }
1134
1135
1136
1137
1138
1139
1140
1141
1142 int main(int argc,char **argv)
1143 {
1144 /* set the stack limit to 500 megabytes */
1145 const rlim_t stacksize = 2000*1024*1024; // 2 gigabytes
1146 struct rlimit rl;
1147 int result;
1148 result = getrlimit(RLIMIT_STACK, &rl);
1149 if (result == 0)
1150 {
1151 if (rl.rlim_cur < stacksize)
1152 {
1153 rl.rlim_cur = stacksize;
1154 result = setrlimit(RLIMIT_STACK, &rl);
1155 if (result != 0)
1156 {
1157 fprintf(stderr, "setrlimit returned result = %d\n", result);
1158 exit(-1);
1159 }
1160 }
1161 }
1162 else {
1163 fprintf(stderr, "getrlimit returned result = %d\n", result);
1164 exit(-1);
1165 }
1166
1167 /* seed1 ip_Address fname */
1168 char *iname=argv[1];
1169 char *oname=argv[2];
1170
1171 unsigned char *otherrand;
1172 int other;
1173 unsigned short *rbuf;
1174 unsigned short *rbackbuf;
1175 int number_of_points;
1176 int number_second_points;
1177 struct point *points_to_look_at;
1178 struct point *second_points_to_look_at;
1179 struct point *end_points_to_look_at;
1180 struct point *temp_points;
1181 int number_end_points;
1182 MagickWand *m_wand = NULL;
1183 PixelWand *c_wand = NULL;
1184 PixelWand *c2_wand = NULL;
1185 MagickPixelPacket pixel;
1186
1187 MagickWandGenesis();
1188 ssize_t width,height;
1189
1190
1191
1192 read_image(iname,&rbuf,&width,&height,&points_to_look_at,&number_of_points);
1193 rbackbuf = malloc(width*height*2*sizeof(unsigned short));
1194 if (rbackbuf == NULL) {fprintf(stderr,"memory\n");exit(-1);}
1195 second_points_to_look_at = malloc(width*height*sizeof(struct point));
1196 if (second_points_to_look_at == NULL) {fprintf(stderr,"memory\n");exit(-1);}
1197 end_points_to_look_at = malloc(width*height*sizeof(struct point));
1198 if (end_points_to_look_at == NULL) {fprintf(stderr,"memory\n");exit(-1);}
1199 temp_points = malloc(width*height*sizeof(struct point));
1200 if (temp_points == NULL) {fprintf(stderr,"memory\n");exit(-1);}
1201
1202
1203
1204 /*
1205 So we just find the innies and outies. The innies will have 32769 and the outies will have 32767
1206 One layer thick around all of everything.
1207 */
1208 /* OK we find an inside point. Once found, we set everytihng to 1 */
1209 {
1210 int cnt;
1211 cnt = 0;
1212 int x,y;
1213 while (cnt<100000) {
1214 x = rand()%width;
1215 y = rand()%height;
1216 if (rbuf[y*width+x] != 32768) {
1217 if (inside(rbuf,width,height,x,y)) {
1218 break;
1219 }
1220 }
1221 cnt++;
1222 }
1223 if (cnt>=100000) {
1224 fprintf(stderr,"error %s missing green\n",argv[1]);
1225 exit(-1);
1226 }
1227 fill_inside(rbuf,x,y,width,height);
1228 }
1229
1230 /* check to see that we have some outside */
1231 {int cx,cy;
1232 for (cx=0;cx<width;cx++) {
1233 for (cy=0;cy<height;cy++) {
1234 if (rbuf[cy*width+cx]==0) break;
1235 }
1236 if (cy<height) break;
1237 }
1238 if ((cx>=width)&&(cy>=height)) {
1239 fprintf(stderr,"error %s does not close up\n",argv[1]);
1240 exit(-1);
1241 }
1242 }
1243
1244
1245 #define P(x,y) (((y)<height)&&((x)<width)&&((y)>=0)&&((x)>=0)&&(rbuf[(y)*width+(x)]==32768))
1246 /* prune inside to within one of the green line and also add the outside border */
1247 int j;
1248 struct point *ht;
1249 ht = NULL;
1250 struct point *endpoint_ht;
1251 endpoint_ht = NULL;
1252 j=0;
1253 {int cx,cy;
1254 for (cx=0;cx<width;cx++) {
1255 for (cy=0;cy<height;cy++) {
1256 if (rbuf[cy*width+cx]==32769) {
1257 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)) {
1258 add_if_not_in_list(cx,cy,second_points_to_look_at,&j,&ht);
1259 }
1260 else {
1261 rbuf[cy*width+cx]=0;
1262 }
1263 }
1264 else if (rbuf[cy*width+cx]==0) { /* outside */
1265 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)) {
1266 rbuf[cy*width+cx]=32767; /* fix it */
1267 }
1268 }
1269 }
1270 }
1271 }
1272
1273
1274
1275 /* handle inside */
1276 number_end_points = 0;
1277 while (j) {
1278 /* first swap second points with points*/
1279 {
1280 struct point *temp;
1281 temp = second_points_to_look_at;
1282 second_points_to_look_at = points_to_look_at;
1283 points_to_look_at = temp;
1284 }
1285 number_of_points = j;
1286
1287
1288 HASH_CLEAR(hh,ht);
1289 j=0;
1290 while (number_of_points) {
1291 int i=(rand()>>1) % number_of_points;
1292 int cx=points_to_look_at[i].x;
1293 int cy=points_to_look_at[i].y;
1294 number_of_points--;
1295 points_to_look_at[i] = points_to_look_at[number_of_points];
1296 int x,y;
1297 int oldj=j; // used to check if we are an end point
1298 x=cx-1;
1299 y=cy;
1300 if (x>=0) {
1301 if (rbuf[y*width+x]==0) { // do we not have a value
1302 rbuf[y*width+x] = rbuf[cy*width+cx]+ 1; /* go up!*/
1303 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1304 } /* if we are ripe for marking as inside or outside */
1305 } /* if we should check this pixel */
1306 x=cx-1;
1307 y=cy-1;
1308 if ((x>=0)&&(y>=0)) {
1309 if (rbuf[y*width+x]==0) { // do we not have a value
1310 rbuf[y*width+x] = rbuf[cy*width+cx] + 1; /* go up! But diagonally */
1311 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1312 } /* if we are ripe for marking as inside or outside */
1313 } /* if we should check this pixel */
1314 x=cx;
1315 y=cy-1;
1316 if (y>=0) {
1317 if (rbuf[y*width+x]==0) { // do we not have a value
1318 rbuf[y*width+x] = rbuf[cy*width+cx]+ 1; /* go up!*/
1319 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1320 } /* if we are ripe for marking as inside or outside */
1321 } /* if we should check this pixel */
1322 x=cx+1;
1323 y=cy-1;
1324 if ((x<width)&&(y>=0)) {
1325 if (rbuf[y*width+x]==0) { // do we not have a value
1326 rbuf[y*width+x] = rbuf[cy*width+cx] + 1; /* go up! But diagonally */
1327 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1328 } /* if we are ripe for marking as inside or outside */
1329 } /* if we should check this pixel */
1330 x=cx+1;
1331 y=cy;
1332 if (x<width) {
1333 if (rbuf[y*width+x]==0) { // do we not have a value
1334 rbuf[y*width+x] = rbuf[cy*width+cx]+ 1; /* go up!*/
1335 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1336 } /* if we are ripe for marking as inside or outside */
1337 } /* if we should check this pixel */
1338 x=cx+1;
1339 y=cy+1;
1340 if ((x<width)&&(y<height)) {
1341 if (rbuf[y*width+x]==0) { // do we not have a value
1342 rbuf[y*width+x] = rbuf[cy*width+cx] + 1; /* go up! But diagonally */
1343 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1344 } /* if we are ripe for marking as inside or outside */
1345 } /* if we should check this pixel */
1346 x=cx;
1347 y=cy+1;
1348 if (y<height) {
1349 if (rbuf[y*width+x]==0) { // do we not have a value
1350 rbuf[y*width+x] = rbuf[cy*width+cx]+ 1; /* go up!*/
1351 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1352 } /* if we are ripe for marking as inside or outside */
1353 } /* if we should check this pixel */
1354 x=cx-1;
1355 y=cy+1;
1356 if ((x>=0)&&(y<height)) {
1357 if (rbuf[y*width+x]==0) { // do we not have a value
1358 rbuf[y*width+x] = rbuf[cy*width+cx] + 1; /* go up! But diagonally */
1359 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1360 } /* if we are ripe for marking as inside or outside */
1361 } /* if we should check this pixel */
1362 if (j==oldj) { // we are an end point
1363 /* well, we might be, or we might be snuffed from another
1364 point, in which case, some of our neighbors will be in the list */
1365 if (!in_list(cx+1,cy ,second_points_to_look_at,j,ht) &&
1366 !in_list(cx+1,cy+1,second_points_to_look_at,j,ht) &&
1367 !in_list(cx ,cy+1,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-1,cy ,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 ,cy-1,second_points_to_look_at,j,ht) &&
1372 !in_list(cx+1,cy-1,second_points_to_look_at,j,ht)) {
1373 end_points_to_look_at[number_end_points].x=cx;
1374 end_points_to_look_at[number_end_points].y=cy;
1375 end_points_to_look_at[number_end_points].pos_count=1; /* we will do the tally later */
1376 {
1377 struct point *c;
1378 c=end_points_to_look_at+number_end_points;
1379 c->key = c->y*30000+c->x;
1380 HASH_ADD_INT(endpoint_ht,key,c);
1381 number_end_points++;
1382 }
1383 }
1384 }
1385 } /* while we have points from this round to do */
1386 } /* while we have points to do */
1387
1388 HASH_CLEAR(hh,ht);
1389
1390
1391
1392 /* OK - well, I don't trust the end points yet - its still a bit stuffy, so I will
1393 see if there are any neightbors that are higher */
1394 int i;
1395 i=0;
1396 while (i<number_end_points) {
1397 int cx = end_points_to_look_at[i].x;
1398 int cy = end_points_to_look_at[i].y;
1399
1400 if (is_brighter(rbuf,width,height,cx,cy,1,0)) { /* 5-8 is good. 9 is too slow */
1401 number_end_points--;
1402 end_points_to_look_at[i]= end_points_to_look_at[number_end_points];
1403 }
1404 else {
1405 i++;
1406 }
1407 }
1408
1409
1410 /* now we need to do the outside - finding the end points is a bit easier, but we don't need them anyways */
1411 /* to start, collect the outside points - which are value 32767 */
1412 j=0;
1413 {
1414 int cx,cy;
1415 for (cy=0;cy<height;cy++) {
1416 for (cx=0;cx<width;cx++) {
1417 unsigned short cur=rbuf[cy*width+cx];
1418 if (cur ==32767) {
1419 second_points_to_look_at[j].x = cx;
1420 second_points_to_look_at[j].y = cy;
1421 j++;
1422 }
1423 }
1424 }
1425 }
1426
1427
1428 int level;
1429 level=32767;
1430 while (j) {
1431 level--;
1432 /* first swap second points with points*/
1433 {
1434 struct point *temp;
1435 temp = second_points_to_look_at;
1436 second_points_to_look_at = points_to_look_at;
1437 points_to_look_at = temp;
1438 }
1439 number_of_points = j;
1440
1441 HASH_CLEAR(hh,ht);
1442 j=0;
1443 while (number_of_points) {
1444 int i=(rand()>>1) % number_of_points;
1445 int cx=points_to_look_at[i].x;
1446 int cy=points_to_look_at[i].y;
1447 number_of_points--;
1448 points_to_look_at[i] = points_to_look_at[number_of_points];
1449 int x,y;
1450 int oldj=j; // used to check if we are an end point
1451 x=cx-1;
1452 y=cy;
1453 if (x>=0) {
1454 if (rbuf[y*width+x]==0) { // do we not have a value
1455 rbuf[y*width+x] = level;
1456 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1457 } /* if we are ripe for marking as inside or outside */
1458 } /* if we should check this pixel */
1459 x=cx-1;
1460 y=cy-1;
1461 if ((x>=0)&&(y>=0)) {
1462 if (rbuf[y*width+x]==0) { // do we not have a value
1463 rbuf[y*width+x] = level;
1464 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1465 } /* if we are ripe for marking as inside or outside */
1466 } /* if we should check this pixel */
1467 x=cx;
1468 y=cy-1;
1469 if (y>=0) {
1470 if (rbuf[y*width+x]==0) { // do we not have a value
1471 rbuf[y*width+x] = level;
1472 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1473 } /* if we are ripe for marking as inside or outside */
1474 } /* if we should check this pixel */
1475 x=cx+1;
1476 y=cy-1;
1477 if ((x<width)&&(y>=0)) {
1478 if (rbuf[y*width+x]==0) { // do we not have a value
1479 rbuf[y*width+x] = level;
1480 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1481 } /* if we are ripe for marking as inside or outside */
1482 } /* if we should check this pixel */
1483 x=cx+1;
1484 y=cy;
1485 if (x<width) {
1486 if (rbuf[y*width+x]==0) { // do we not have a value
1487 rbuf[y*width+x] = level;
1488 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1489 } /* if we are ripe for marking as inside or outside */
1490 } /* if we should check this pixel */
1491 x=cx+1;
1492 y=cy+1;
1493 if ((x<width)&&(y<height)) {
1494 if (rbuf[y*width+x]==0) { // do we not have a value
1495 rbuf[y*width+x] = level;
1496 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1497 } /* if we are ripe for marking as inside or outside */
1498 } /* if we should check this pixel */
1499 x=cx;
1500 y=cy+1;
1501 if (y<height) {
1502 if (rbuf[y*width+x]==0) { // do we not have a value
1503 rbuf[y*width+x] = level;
1504 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1505 } /* if we are ripe for marking as inside or outside */
1506 } /* if we should check this pixel */
1507 x=cx-1;
1508 y=cy+1;
1509 if ((x>=0)&&(y<height)) {
1510 if (rbuf[y*width+x]==0) { // do we not have a value
1511 rbuf[y*width+x] = level;
1512 add_if_not_in_list(x,y,second_points_to_look_at,&j,&ht);
1513 } /* if we are ripe for marking as inside or outside */
1514 } /* if we should check this pixel */
1515 } /* while we have points from this round to do */
1516 } /* while we have points to do */
1517
1518
1519
1520
1521
1522
1523
1524
1525 /*OK - lets plateau the inside so that we always have one object
1526 note - we use the back buffer as a scratch space*/
1527 {
1528 unsigned short max_shade_point = get_highest_shade(rbuf,width,height);
1529 unsigned short min_shade_point = 32769;
1530 unsigned short shade_point=32769;
1531 int movement=(max_shade_point-min_shade_point-1)/20;
1532 if (movement<1) movement=1;
1533 int tally;
1534 tally=number_of_shapes_at_a_shade(shade_point,rbuf,rbackbuf,width,height,points_to_look_at,second_points_to_look_at,temp_points);
1535 if (tally == 1) {
1536
1537 while (1) { // While moving throght the sample
1538 shade_point=min_shade_point;
1539 tally=1;
1540 while (1) {
1541 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);
1542 // 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);
1543 if (tally !=1) break;
1544 shade_point += movement;
1545 if (shade_point == max_shade_point) break;
1546 if (shade_point+movement > max_shade_point) shade_point = max_shade_point-movement;
1547 }
1548 if (tally==1) break; // all good
1549 min_shade_point = shade_point;
1550 max_shade_point = shade_point+movement-1;
1551 if (min_shade_point == max_shade_point) break; // also all good
1552 movement=(max_shade_point-min_shade_point-1)/20;
1553 if (movement<1) movement=1;
1554 // fprintf(stderr,"plateau new min %d max %d movement %d\n",min_shade_point,max_shade_point,movement);
1555 }
1556 /* clip it */
1557 {
1558 int cx,cy;
1559 for (cy=0;cy<height;cy++) {
1560 for (cx=0;cx<width;cx++) {
1561 if (rbuf[cy*width+cx] > shade_point) {
1562 rbuf[cy*width+cx] = shade_point;
1563 }
1564 }
1565 }
1566 }
1567 } /* if we have one shape at the first level */
1568 else {
1569 fprintf(stderr,"for %s no clipping - lowest level does not work\n",argv[1]);
1570 }
1571 }
1572
1573
1574
1575
1576
1577 /* copy to the back buffer*/
1578 {
1579 int cx,cy;
1580 for (cy=0;cy<height;cy++) {
1581 for (cx=0;cx<width;cx++) {
1582 double cur= rbuf[cy*width+cx];
1583 if (cur==32768) {
1584 rbackbuf[cy*width+cx] = cur;
1585 }
1586 else {
1587 rbackbuf[cy*width+cx] = 0.;
1588 }
1589 }
1590 }
1591 }
1592
1593
1594
1595
1596
1597 /* OK - now we need to add the end points - which are always at the edge of the screen
1598 We set the values based on the fact that anything after should be zero, and the number
1599 we have lets us know how many gradients we are from 0
1600 First build the list
1601 */
1602 j=0;
1603 for (i=0;i<height;i++) {
1604 if (rbuf[i*width+0] < 32768) {
1605 second_points_to_look_at[j].x = 0;
1606 second_points_to_look_at[j].y = i;
1607 second_points_to_look_at[j].pos_tally = 0.;
1608 second_points_to_look_at[j].pos_count = 1.;
1609 j++;
1610 }
1611 if (rbuf[i*width+(width-1)] < 32768) {
1612 second_points_to_look_at[j].x = width-1;
1613 second_points_to_look_at[j].y = i;
1614 second_points_to_look_at[j].pos_tally = 0.;
1615 second_points_to_look_at[j].pos_count = 1.;
1616 j++;
1617 }
1618 }
1619 for (i=0;i<width;i++) {
1620 if (rbuf[0*width+i] < 32768) {
1621 second_points_to_look_at[j].x = i;
1622 second_points_to_look_at[j].y = 0;
1623 second_points_to_look_at[j].pos_tally = 0.;
1624 second_points_to_look_at[j].pos_count = 1.;
1625 j++;
1626 }
1627 if (rbuf[(height-1)*width+i] < 32768) {
1628 second_points_to_look_at[j].x = i;
1629 second_points_to_look_at[j].y = height-1;
1630 second_points_to_look_at[j].pos_tally = 0.;
1631 second_points_to_look_at[j].pos_count = 1.;
1632 j++;
1633 }
1634 }
1635
1636 sweep_outside(points_to_look_at,second_points_to_look_at,rbuf,rbackbuf,j,width,height);
1637
1638
1639
1640
1641 /* clean up the inside
1642 You see, the outside all goes to zero, at the edge, well, for the most part anyways
1643 But the inside could have these multiple points. We can't deal with that
1644 because then the picture would have multiple points showing - and we need it always to be one
1645 solid single blob
1646 */
1647 #ifdef waitfdf
1648 unsigned short max=32769.0;
1649 {
1650 int cx,cy;
1651 for (cy=0;cy<height;cy++) {
1652 for (cx=0;cx<width;cx++) {
1653 unsigned short cur=rbuf[cy*width+cx];
1654 if (cur > max) max=cur;
1655 }
1656 }
1657 }
1658 /* OK so max is now equal to the highest point. Lets see how many points we have
1659 like that
1660 */
1661 j=0;
1662 {
1663 int cx,cy;
1664 for (cy=0;cy<height;cy++) {
1665 for (cx=0;cx<width;cx++) {
1666 unsigned short cur=rbuf[cy*width+cx];
1667 if (cur > max) max=cur;
1668 }
1669 }
1670 }
1671 #endif
1672
1673
1674
1675
1676 /* lets do the inside now as well */
1677
1678
1679 unsigned short size = get_highest_shade(rbuf,width,height);
1680 double scale = (size-32768);
1681 double full_over_scale = 32767.0 / scale;
1682
1683
1684
1685
1686
1687
1688 /* Cool, now we have to "kern" the inside values of rbuf from the start point. */
1689
1690 /* so this part was done above - all we have to do is apply scale and full_over_scale
1691 The code is still here but commented so you can see it logically */
1692 #ifdef done_Earlier
1693 unsigned short size=get_highest_shade(rbuf,width,height);
1694 double scale = (size-32768);
1695 double full_over_scale = 32767.0 / scale;
1696 #endif
1697 {
1698 int cx,cy;
1699 for (cy=0;cy<height;cy++) {
1700 for (cx=0;cx<width;cx++) {
1701 double cur= rbuf[cy*width+cx];
1702 if (cur>32768) {
1703 cur = ((cur -32768.) * full_over_scale) + 32768.;
1704 rbuf[cy*width+cx] = cur;
1705 rbackbuf[cy*width+cx] = cur; /* because we are not sweeping inside */
1706 }
1707 }
1708 }
1709 }
1710
1711
1712
1713 /* Cool, now we have to "kern" the outside values from the start point. */
1714 size=get_lowest_shade(rbuf,width,height);
1715 scale = 32768-(size);
1716 full_over_scale = 32767.0 / scale;
1717
1718
1719 {
1720 int cx,cy;
1721 for (cy=0;cy<height;cy++) {
1722 for (cx=0;cx<width;cx++) {
1723 double cur= rbuf[cy*width+cx];
1724 if (cur<32768) {
1725 cur = 32768 - (32768 - cur) * full_over_scale;
1726 rbuf[cy*width+cx] = cur;
1727 }
1728 }
1729 }
1730 }
1731
1732
1733
1734
1735 PixelIterator *iterator;
1736 PixelWand **pixels;
1737 m_wand = NewMagickWand();
1738 c_wand = NewPixelWand();
1739 PixelSetColor(c_wand,"black"); /* for now */
1740 MagickNewImage(m_wand,width,height,c_wand);
1741 MagickSetImageChannelDepth(m_wand, AllChannels, 16);
1742 MagickSetImageDepth(m_wand, 16);
1743
1744 iterator=NewPixelIterator(m_wand);
1745 if ((iterator == (PixelIterator *) NULL) )
1746 ThrowWandException(m_wand);
1747
1748
1749 register ssize_t x;
1750 ssize_t y;
1751 for (y=0; y < (ssize_t) height; y++)
1752 {
1753 pixels=PixelGetNextIteratorRow(iterator,&width);
1754 if ((pixels == (PixelWand **) NULL))
1755 break;
1756 for (x=0; x < (ssize_t) width; x++)
1757 {
1758 PixelGetMagickColor(pixels[x],&pixel);
1759 pixel.red=pixel.green=pixel.blue=rbackbuf[y*width+x];
1760 PixelSetMagickColor(pixels[x],&pixel);
1761 }
1762 (void) PixelSyncIterator(iterator);
1763 }
1764 if (y < (ssize_t) MagickGetImageHeight(m_wand))
1765 ThrowWandException(m_wand);
1766 iterator=DestroyPixelIterator(iterator);
1767
1768
1769 MagickGaussianBlurImage(m_wand,20.,2.);
1770
1771 /*
1772 Write the image then destroy it.
1773 */
1774 MagickBooleanType status;
1775 status=MagickWriteImages(m_wand,oname,MagickTrue);
1776 if (status == MagickFalse)
1777 ThrowWandException(m_wand);
1778 m_wand=DestroyMagickWand(m_wand);
1779 MagickWandTerminus();
1780 return(0);
1781 }
1782

  ViewVC Help
Powered by ViewVC 1.1.5