1 /*
2 * li_recognizer.c
3 *
4 * Copyright 2000 Compaq Computer Corporation.
5 * Copying or modifying this code for any purpose is permitted,
6 * provided that this copyright notice is preserved in its entirety
7 * in all copies or modifications.
8 * COMPAQ COMPUTER CORPORATION MAKES NO WARRANTIES, EXPRESSED OR
9 * IMPLIED, AS TO THE USEFULNESS OR CORRECTNESS OF THIS CODE OR
10 *
11 *
12 * Adapted from cmu_recognizer.c by Jay Kistler.
13 *
14 * Where is the CMU copyright???? Gotta track it down - Jim Gettys
15 *
16 * Credit to Dean Rubine, Jim Kempf, and Ari Rapkin.
17 */
18
19 #include <u.h>
20 #include <libc.h>
21 #include <stdio.h>
22 #include <draw.h>
23 #include <scribble.h>
24 #include "scribbleimpl.h"
25
26 #include "hre_internal.h"
27 #include "li_recognizer_internal.h"
28
29 int lidebug = 0;
30
31 #define LI_MAGIC 0xACCBADDD
32
33 #define CHECK_LI_MAGIC(_a) \
34 ((_a) != nil && ((li_recognizer*)(_a))->li_magic == LI_MAGIC)
35
36
37 static void lialg_initialize(rClassifier *);
38 static int lialg_read_classifier_digest(rClassifier *);
39 static int lialg_canonicalize_examples(rClassifier *);
40 static char *lialg_recognize_stroke(rClassifier *, point_list *);
41 static void lialg_compute_lpf_parameters(void);
42
43
44 char* li_err_msg = nil;
45
46 #define bcopy(s1,s2,n) memmove(s2,s1,n)
47
48 /*Freeing classifier*/
49
50 static void
51 free_rClassifier(rClassifier* rc);
52
53 /*
54 * Point List Support
55 */
56
57 static point_list*
add_example(point_list * l,int npts,pen_point * pts)58 add_example(point_list* l,int npts,pen_point* pts)
59 {
60 pen_point* lpts = mallocz(npts*sizeof(pen_point), 1);
61 point_list *p = malloc(sizeof(point_list));
62
63 p->npts = npts;
64 p->pts = lpts;
65 p->next = l; /*Order doesn't matter, so we stick on end.*/
66
67 /*Copy points.*/
68
69 bcopy(pts, lpts, npts * sizeof(pen_point));
70
71 return(p);
72 }
73
74
75 static void
delete_examples(point_list * l)76 delete_examples(point_list* l)
77 {
78 point_list* p;
79
80 for( ; l != nil; l = p ) {
81 p = l->next;
82 free(l->pts);
83 free(l);
84 }
85 }
86
87 /*
88 * Local functions
89 */
90
91 /*
92 * recognize_internal-Form Vector, use Classifier to classify, return char.
93 */
94
95 static char*
recognize_internal(rClassifier * rec,Stroke * str,int *)96 recognize_internal(rClassifier* rec, Stroke* str, int*)
97 {
98 char *res;
99 point_list *stroke;
100
101 stroke = add_example(nil, str->npts, str->pts);
102 if (stroke == nil) return(nil);
103
104 res = lialg_recognize_stroke(rec, stroke);
105
106 delete_examples(stroke);
107 return(res);
108 }
109
110 /*
111 * file_path-Construct pathname, check for proper extension.
112 */
113
114 static int
file_path(char * dir,char * filename,char * pathname)115 file_path(char* dir,char* filename,char* pathname)
116 {
117 char* dot;
118
119 /*Check for proper extension on file name.*/
120
121 dot = strrchr(filename,'.');
122
123 if( dot == nil ) {
124 return(-1);
125 }
126
127 /*Determine whether a gesture or character classifier.*/
128
129 if( strcmp(dot,LI_CLASSIFIER_EXTENSION) != 0 ) {
130 return(-1);
131 }
132
133 /*Concatenate directory and filename into pathname.*/
134
135 strcpy(pathname,dir);
136 strcat(pathname,"/");
137 strcat(pathname,filename);
138
139 return(0);
140 }
141
142 /*read_classifier_points-Read points so classifier can be extended.*/
143
144 static int
read_classifier_points(FILE * fd,int nclss,point_list ** ex,char ** cnames)145 read_classifier_points(FILE* fd,int nclss,point_list** ex,char** cnames)
146 {
147 int i,j,k;
148 char buf[BUFSIZ];
149 int nex = 0;
150 char* names[MAXSCLASSES];
151 point_list* examples[MAXSCLASSES];
152 pen_point* pts;
153 int npts;
154
155 /*Initialize*/
156
157 for( i = 0; i < MAXSCLASSES; i++ ) {
158 names[i] = nil;
159 examples[i] = nil;
160 }
161
162 /*Go thru classes.*/
163
164 for( k = 0; k < nclss; k++ ) {
165
166 /*Read class name and number of examples.*/
167
168 if( fscanf(fd,"%d %s",&nex,buf) != 2 )
169 goto unallocate;
170
171 /*Save class name*/
172
173 names[k] = strdup(buf);
174
175 /*Read examples.*/
176
177 for( i = 0; i < nex; i++ ) {
178
179 /*Read number of points.*/
180
181 if( fscanf(fd,"%d",&npts) != 1 )
182 goto unallocate; /*Boy would I like exceptions!*/
183
184 /*Allocate array for points.*/
185
186 if( (pts = mallocz(npts*sizeof(pen_point), 1)) == nil )
187 goto unallocate;
188
189 /*Read in points.*/
190
191 for( j = 0; j < npts; j++ ) {
192 int x,y;
193 if( fscanf(fd,"%d %d",&x,&y) != 2 ) {
194 delete_pen_point_array(pts);
195 goto unallocate;
196 }
197 pts[j].Point = Pt(x, y);
198 }
199
200 /*Add example*/
201
202 if( (examples[k] = add_example(examples[k],npts,pts)) == nil ) {
203 delete_pen_point_array(pts);
204 goto unallocate;
205 }
206
207 delete_pen_point_array(pts);
208
209 }
210 }
211
212 /* ari -- end of list of classes */
213 /* fprint(2, "]\n"); */
214
215 /*Transfer to recognizer.*/
216
217 bcopy(examples,ex,sizeof(examples));
218 bcopy(names,cnames,sizeof(names));
219
220 return(0);
221
222 /*Error. Deallocate memory and return.*/
223
224 unallocate:
225 for( ; k >= 0; k-- ) {
226 delete_examples(examples[k]);
227 free(names[k]);
228 }
229
230 return(-1);
231 }
232
233 /*read_classifier-Read a classifier file.*/
234
read_classifier(FILE * fd,rClassifier * rc)235 static int read_classifier(FILE* fd,rClassifier* rc)
236 {
237
238 li_err_msg = nil;
239
240 /*Read in classifier file.*/
241
242 if(fscanf(fd, "%d", &rc->nclasses) != 1)
243 return -1;
244
245 /*Read in the example points, so classifier can be extended.*/
246
247 if( read_classifier_points(fd,rc->nclasses,rc->ex,rc->cnames) != 0 )
248 return -1;
249
250 return(0);
251 }
252
253 /*
254 * Extension Functions
255 */
256
257 /* getClasses and clearState are by Ari */
258
259 static int
recognizer_getClasses(recognizer r,char *** list,int * nc)260 recognizer_getClasses (recognizer r, char ***list, int *nc)
261 {
262 int i, nclasses;
263 li_recognizer* rec;
264 char **ret;
265
266 rec = (li_recognizer*)r->recognizer_specific;
267
268 /*Check for LI recognizer.*/
269
270 if( !CHECK_LI_MAGIC(rec) ) {
271 li_err_msg = "Not a LI recognizer";
272 return(-1);
273 }
274
275 *nc = nclasses = rec->li_rc.nclasses;
276 ret = malloc(nclasses*sizeof(char*));
277
278 for (i = 0; i < nclasses; i++)
279 ret[i] = rec->li_rc.cnames[i]; /* only the 1st char of the cname */
280 *list = ret;
281 return 0;
282 }
283
284 static int
recognizer_clearState(recognizer)285 recognizer_clearState (recognizer)
286 {
287 /*This operation isn't supported by the LI recognizer.*/
288
289 li_err_msg = "Clearing state is not supported by the LI recognizer";
290
291 return(-1);
292 }
293
isa_li(recognizer r)294 static bool isa_li(recognizer r)
295 { return(CHECK_LI_MAGIC(r)); }
296
297 static int
recognizer_train(recognizer,rc *,uint,Stroke *,rec_element *,bool)298 recognizer_train(recognizer, rc*, uint, Stroke*, rec_element*, bool)
299 {
300 /*This operation isn't supported by the LI recognizer.*/
301
302 li_err_msg = "Training is not supported by the LI recognizer";
303
304 return(-1);
305 }
306
307 int
li_recognizer_get_example(recognizer r,int class,int instance,char ** name,pen_point ** points,int * npts)308 li_recognizer_get_example (recognizer r,
309 int class,
310 int instance,
311 char **name,
312 pen_point **points,
313 int *npts)
314 {
315 li_recognizer *rec = (li_recognizer*)r->recognizer_specific;
316 int nclasses = rec->li_rc.nclasses;
317 point_list *pl;
318
319 if( !CHECK_LI_MAGIC(rec) ) {
320 li_err_msg = "Not a LI recognizer";
321 return(-1);
322 }
323 if (class > nclasses)
324 return -1;
325 pl = rec->li_rc.canonex[class];
326 while (instance && pl)
327 {
328 pl = pl->next;
329 instance--;
330 }
331 if (!pl)
332 return -1;
333 *name = rec->li_rc.cnames[class];
334 *points = pl->pts;
335 *npts = pl->npts;
336 return pl->npts; /* I hope [sjm] */
337 }
338
339 /*
340 * API Functions
341 */
342
343
344 /*li_recognizer_load-Load a classifier file.*/
345
li_recognizer_load(recognizer r,char * dir,char * filename)346 static int li_recognizer_load(recognizer r, char* dir, char* filename)
347 {
348 FILE *fd;
349 char* pathname;
350 li_recognizer* rec;
351 rClassifier* rc;
352
353 rec = (li_recognizer*)r->recognizer_specific;
354
355 /*Make sure recognizer's OK*/
356
357 if( !CHECK_LI_MAGIC(rec) ) {
358 li_err_msg = "Not a LI recognizer";
359 return(-1);
360 }
361
362 rc = &(rec->li_rc);
363
364 /*Check parameters.*/
365
366 if( filename == nil ) {
367 li_err_msg = "Invalid parameters";
368 return(-1);
369 }
370
371 /*We let the directory be null.*/
372
373 if( dir == nil || (int)strlen(dir) <= 0 ) {
374 dir = ".";
375 }
376
377 if(0)fprint(2, "dir = %s, filename = %s\n",
378 dir, filename);
379
380 /*Make full pathname and check filename*/
381
382 pathname = malloc((strlen(dir) + strlen(filename) + 2)*sizeof(char));
383 if( file_path(dir, filename, pathname) == -1 ) {
384 free(pathname);
385 li_err_msg = "Not a LI recognizer classifier file";
386 return -1;
387 }
388
389 /* Try to short-circuit the full classifier-file processing. */
390 rc->file_name = pathname;
391 if (lialg_read_classifier_digest(rc) == 0)
392 return(0);
393 rc->file_name = nil;
394
395 /*Open the file*/
396
397 if( (fd = fopen(pathname,"r")) == nil ) {
398 li_err_msg = "Can't open classifier file";
399 if(0)fprint(2, "Can't open %s.\n", pathname);
400 free(pathname);
401 return(-1);
402 }
403
404 /*If rClassifier is OK, then delete it first.*/
405
406 if( rc->file_name != nil ) {
407 free_rClassifier(rc);
408 }
409
410 /*Read classifier.*/
411
412 if( read_classifier(fd,rc) < 0 ) {
413 free(pathname);
414 return(-1);
415 }
416
417 /*Close file.*/
418
419 fclose(fd);
420
421 /*Add classifier name.*/
422
423 rc->file_name = pathname;
424
425 /* Canonicalize examples. */
426 if (lialg_canonicalize_examples(rc) != 0) {
427 free(pathname);
428 rc->file_name = nil;
429 return -1;
430 }
431
432 return(0);
433 }
434
435 /*li_recognizer_save-Save a classifier file.*/
436
li_recognizer_save(recognizer,char *,char *)437 static int li_recognizer_save(recognizer, char*, char*)
438 {
439 /*This operation isn't supported by the LI recognizer.*/
440
441 li_err_msg = "Saving is not supported by the LI recognizer";
442 return -1;
443 }
444
445 static wordset
li_recognizer_load_dictionary(recognizer,char *,char *)446 li_recognizer_load_dictionary(recognizer, char*, char*)
447 {
448 /*This operation isn't supported by the LI recognizer.*/
449
450 li_err_msg = "Dictionaries are not supported by the LI recognizer";
451 return nil;
452 }
453
454 static int
li_recognizer_save_dictionary(recognizer,char *,char *,wordset)455 li_recognizer_save_dictionary(recognizer, char*, char*, wordset)
456 {
457 /*This operation isn't supported by the LI recognizer.*/
458 li_err_msg = "Dictionaries are not supported by the LI recognizer";
459 return -1;
460 }
461
462 static int
li_recognizer_free_dictionary(recognizer,wordset)463 li_recognizer_free_dictionary(recognizer, wordset)
464 {
465 /*This operation isn't supported by the LI recognizer.*/
466
467 li_err_msg = "Dictionaries are not supported by the LI recognizer";
468
469 return -1;
470
471 }
472
473 static int
li_recognizer_add_to_dictionary(recognizer,letterset *,wordset)474 li_recognizer_add_to_dictionary(recognizer, letterset*, wordset)
475 {
476 /*This operation isn't supported by the LI recognizer.*/
477 li_err_msg = "Dictionaries are not supported by the LI recognizer";
478 return -1;
479 }
480
481 static int
li_recognizer_delete_from_dictionary(recognizer,letterset *,wordset)482 li_recognizer_delete_from_dictionary(recognizer, letterset*, wordset)
483 {
484 /*This operation isn't supported by the LI recognizer.*/
485 li_err_msg = "Dictionaries are not supported by the LI recognizer";
486 return -1;
487 }
488
489 static char*
li_recognizer_error(recognizer rec)490 li_recognizer_error(recognizer rec)
491 {
492 char* ret = li_err_msg;
493
494 /*Check for LI recognizer.*/
495
496 if( !CHECK_LI_MAGIC(rec->recognizer_specific) ) {
497 li_err_msg = "Not a LI recognizer";
498 return nil;
499 }
500 li_err_msg = nil;
501 return ret;
502 }
503
504 static int
li_recognizer_clear(recognizer r,bool)505 li_recognizer_clear(recognizer r, bool)
506 {
507 li_recognizer* rec;
508
509 rec = (li_recognizer*)r->recognizer_specific;
510 /*Check for LI recognizer.*/
511 if( !CHECK_LI_MAGIC(rec) ) {
512 li_err_msg = "Not a LI recognizer";
513 return 0;
514 }
515 return 0;
516 }
517
518 static int
li_recognizer_set_context(recognizer,rc *)519 li_recognizer_set_context(recognizer, rc*)
520 {
521 /*This operation isn't supported by the LI recognizer.*/
522 li_err_msg = "Contexts are not supported by the LI recognizer";
523 return -1;
524 }
525
526 static rc*
li_recognizer_get_context(recognizer)527 li_recognizer_get_context(recognizer)
528 {
529 /*This operation isn't supported by the LI recognizer.*/
530 li_err_msg = "Contexts are not supported by the LI recognizer";
531 return nil;
532 }
533
534 static int
li_recognizer_get_buffer(recognizer,uint *,Stroke **)535 li_recognizer_get_buffer(recognizer, uint*, Stroke**)
536 {
537 /*This operation isn't supported by the LI recognizer.*/
538 li_err_msg = "Buffer get/set are not supported by the LI recognizer";
539 return -1;
540 }
541
542 static int
li_recognizer_set_buffer(recognizer,uint,Stroke *)543 li_recognizer_set_buffer(recognizer, uint, Stroke*)
544 {
545 /*This operation isn't supported by the LI recognizer.*/
546 li_err_msg = "Buffer get/set are not supported by the LI recognizer";
547 return -1;
548 }
549
550 static int
li_recognizer_translate(recognizer r,uint ncs,Stroke * tps,bool,int * nret,rec_alternative ** ret)551 li_recognizer_translate(recognizer r, uint ncs, Stroke* tps, bool, int* nret, rec_alternative** ret)
552 {
553 char* clss;
554 li_recognizer* rec;
555 int conf;
556 rClassifier* rc;
557
558 rec = (li_recognizer*)r->recognizer_specific;
559
560 *nret = 0;
561 *ret = nil;
562
563 /*Check for LI recognizer.*/
564
565 if( !CHECK_LI_MAGIC(rec) ) {
566 li_err_msg = "Not a LI recognizer";
567 return(-1);
568 }
569
570 rc = &(rec->li_rc);
571
572 /*Check for valid parameters.*/
573 if (ncs < 1) {
574 li_err_msg = "Invalid parameters: ncs";
575 return(-1);
576 }
577 if( tps == nil) {
578 li_err_msg = "Invalid parameters: tps";
579 return(-1);
580 }
581 if( nret == nil) {
582 li_err_msg = "Invalid parameters: nret";
583 return(-1);
584 }
585 if( ret == nil) {
586 li_err_msg = "Invalid parameters: ret";
587 return(-1);
588 }
589
590 /*
591 * Go through the stroke array and recognize. Since this is a single
592 * stroke recognizer, each stroke is treated as a separate
593 * character or gesture. We allow only characters or gestures
594 * to be recognized at one time, since otherwise, handling
595 * the display of segmentation would be difficult.
596 */
597 clss = recognize_internal(rc,tps,&conf);
598 if (clss == nil) {
599 *nret = 1;
600 return(0);
601 }
602
603 /*Return values.*/
604 *nret = 1;
605 return(*clss);
606 }
607
608
609 static rec_fn*
li_recognizer_get_extension_functions(recognizer rec)610 li_recognizer_get_extension_functions(recognizer rec)
611 {
612 rec_fn* ret;
613
614 /*Check for LI recognizer.*/
615
616 if( !CHECK_LI_MAGIC(rec->recognizer_specific) ) {
617 li_err_msg = "Not a LI recognizer";
618 return(nil);
619 }
620
621 ret = make_rec_fn_array(LI_NUM_EX_FNS);
622
623 /* ari -- clearState & getClasses are mine */
624 ret[LI_GET_CLASSES] = (rec_fn)recognizer_getClasses;
625 ret[LI_CLEAR] = (rec_fn)recognizer_clearState;
626 ret[LI_ISA_LI] = (rec_fn)isa_li;
627 ret[LI_TRAIN] = (rec_fn)recognizer_train;
628
629 return(ret);
630 }
631
632 static char**
li_recognizer_get_gesture_names(recognizer)633 li_recognizer_get_gesture_names(recognizer)
634 {
635 /*This operation isn't supported by the LI recognizer.*/
636 li_err_msg = "Gestures are not supported by the LI recognizer";
637 return nil;
638 }
639
640 static xgesture
li_recognizer_set_gesture_action(recognizer,char *,xgesture,void *)641 li_recognizer_set_gesture_action(recognizer, char*, xgesture, void*)
642 {
643 /*This operation isn't supported by the LI recognizer.*/
644 li_err_msg = "Gestures are not supported by the LI recognizer";
645 return nil;
646 }
647
648 /*
649 * Exported Functions
650 */
651
652 /*RECOGNIZER_INITIALIZE-Initialize the recognizer.*/
653
654 /* note from ari: this expands via pre-processor to
655 *
656 * recognizer __recognizer_internal_initialize(rec_info* ri)
657 */
658
RECOGNIZER_INITIALIZE(ri)659 RECOGNIZER_INITIALIZE(ri)
660 {
661 recognizer r;
662 li_recognizer* rec;
663 int i;
664
665 /*Check that locale matches.*/
666
667 if( strcmp(ri->ri_locale,LI_SUPPORTED_LOCALE) != 0 ) {
668 li_err_msg = "Not a supported locale";
669 /* fprint(2, "Locale error.\n");*/
670 return(nil);
671 }
672
673 /*
674 * Check that character sets match. Note that this is only approximate,
675 * since the classifier file will have more information.
676 */
677
678 if( ri->ri_subset != nil ) {
679 for(i = 0; ri->ri_subset[i] != nil; i++ ) {
680
681 if( strcmp(ri->ri_subset[i],UPPERCASE) != 0 &&
682 strcmp(ri->ri_subset[i],LOWERCASE) != 0 &&
683 strcmp(ri->ri_subset[i],DIGITS) != 0 &&
684 strcmp(ri->ri_subset[i],GESTURE) != 0 ) {
685 li_err_msg = "Not a supported character set";
686 /* fprint(2, "charset error.\n"); */
687
688 return(nil);
689 }
690 }
691 }
692
693 /* ari */
694 r = make_recognizer(ri);
695 /* fprint(2, "past make_recognizer.\n"); */
696
697 if( r == nil ) {
698 li_err_msg = "Can't allocate storage";
699
700 return(nil);
701 }
702
703 /*Make a LI recognizer structure.*/
704
705
706 /* rec = (li_recognizer*)safe_malloc(sizeof(li_recognizer))) == nil ); */
707
708 rec = malloc(sizeof(li_recognizer));
709
710 r->recognizer_specific = rec;
711
712 rec->li_rc.file_name = nil;
713 rec->li_rc.nclasses = 0;
714
715 /*Initialize the recognizer struct.*/
716
717 r->recognizer_load_state = li_recognizer_load;
718 r->recognizer_save_state = li_recognizer_save;
719 r->recognizer_load_dictionary = li_recognizer_load_dictionary;
720 r->recognizer_save_dictionary = li_recognizer_save_dictionary;
721 r->recognizer_free_dictionary = li_recognizer_free_dictionary;
722 r->recognizer_add_to_dictionary = li_recognizer_add_to_dictionary;
723 r->recognizer_delete_from_dictionary = li_recognizer_delete_from_dictionary;
724 r->recognizer_error = li_recognizer_error;
725 r->recognizer_translate = li_recognizer_translate;
726 r->recognizer_get_context = li_recognizer_get_context;
727 r->recognizer_set_context = li_recognizer_set_context;
728 r->recognizer_get_buffer = li_recognizer_get_buffer;
729 r->recognizer_set_buffer = li_recognizer_set_buffer;
730 r->recognizer_clear = li_recognizer_clear;
731 r->recognizer_get_extension_functions =
732 li_recognizer_get_extension_functions;
733 r->recognizer_get_gesture_names = li_recognizer_get_gesture_names;
734 r->recognizer_set_gesture_action =
735 li_recognizer_set_gesture_action;
736
737 /*Initialize LI Magic Number.*/
738
739 rec->li_magic = LI_MAGIC;
740
741 /*Initialize rClassifier.*/
742
743 rec->li_rc.file_name = nil;
744
745 for( i = 0; i < MAXSCLASSES; i++ ) {
746 rec->li_rc.ex[i] = nil;
747 rec->li_rc.cnames[i] = nil;
748 }
749
750 lialg_initialize(&rec->li_rc);
751
752 /*Get rid of error message. Not needed here.*/
753 li_err_msg = nil;
754
755 return(r);
756 }
757
758 /*free_rClassifier-Free the rClassifier.*/
759
760 static void
free_rClassifier(rClassifier * rc)761 free_rClassifier(rClassifier* rc)
762 {
763 int i;
764
765 if( rc->file_name != nil) {
766 free(rc->file_name);
767 }
768
769 for( i = 0; rc->ex[i] != nil; i++) {
770 delete_examples(rc->ex[i]);
771 free(rc->cnames[i]);
772 }
773
774 }
775
776 /*RECOGNIZER_FINALIZE-Deallocate the recognizer, finalize.*/
777
RECOGNIZER_FINALIZE(r)778 RECOGNIZER_FINALIZE(r)
779 {
780 li_recognizer* rec = (li_recognizer*)r->recognizer_specific;
781
782 /*Make sure this is a li_recognizer first*/
783 if( !CHECK_LI_MAGIC(rec) ) {
784 li_err_msg = "Not a LI recognizer";
785 return(-1);
786 }
787
788 /*Deallocate rClassifier resources.*/
789
790 free_rClassifier(&(rec->li_rc));
791
792 /*Deallocate the li_recognizer struct.*/
793 free(rec);
794
795 /*Deallocate the recognizer*/
796 delete_recognizer(r);
797
798 return(0);
799 }
800
801
802 /* **************************************************
803
804 Implementation of the Li/Yeung recognition algorithm
805
806 ************************************************** */
807
808 #define WORST_SCORE 0x7fffffff
809
810 /* Dynamic programming parameters */
811 #define DP_BAND 3
812 #define MIN_SIM 0
813 #define MAX_DIST 0x7fffffff
814 #define SIM_THLD 60 /* x 100 */
815 #define DIST_THLD 3200 /* x 100 */
816
817 /* Low-pass filter parameters -- empirically derived */
818 #define LP_FILTER_WIDTH 6
819 #define LP_FILTER_ITERS 8
820 #define LP_FILTER_THLD 250 /* x 100 */
821 #define LP_FILTER_MIN 5
822
823 /* Pseudo-extrema parameters -- empirically derived */
824 #define PE_AL_THLD 1500 /* x 100 */
825 #define PE_ATCR_THLD 135 /* x 100 */
826
827 /* Contour-angle derivation parameters */
828 #define T_ONE 1
829 #define T_TWO 20
830
831 /* Pre-processing and canonicalization parameters */
832 #define CANONICAL_X 108
833 #define CANONICAL_Y 128
834 #define DIST_SQ_THRESHOLD (3*3) /* copied from fv.h */
835 #define NCANONICAL 50
836
837 /* Tap-handling parameters */
838 #define TAP_CHAR "."
839 #define TAP_TIME_THLD 150 /* msec */
840 #define TAP_DIST_THLD 75 /* dx * dx + dy * dy */
841 #define TAP_PATHLEN 1000 /* x 100 */
842
843
844 /* region types */
845 #define RGN_CONVEX 0
846 #define RGN_CONCAVE 1
847 #define RGN_PLAIN 2
848 #define RGN_PSEUDO 3
849
850
851 typedef struct RegionList {
852 int start;
853 int end;
854 int type;
855 struct RegionList *next;
856 } region_list;
857
858
859 /* direction-code table; indexed by dx, dy */
860 static int lialg_dctbl[3][3] = {{1, 0, 7}, {2, 0x7FFFFFFF, 6}, {3, 4, 5}};
861
862 /* low-pass filter weights */
863 static int lialg_lpfwts[2 * LP_FILTER_WIDTH + 1];
864 static int lialg_lpfconst = -1;
865
866
867 static int lialg_preprocess_stroke(point_list *);
868 static point_list *lialg_compute_dominant_points(point_list *);
869 static point_list *lialg_interpolate_points(point_list *);
870 static void lialg_bresline(pen_point *, pen_point *, point_list *, int *);
871 static void lialg_compute_chain_code(point_list *);
872 static void lialg_compute_unit_chain_code(point_list *);
873 static region_list *lialg_compute_regions(point_list *);
874 static point_list *lialg_compute_dompts(point_list *, region_list *);
875 static int *lialg_compute_contour_angle_set(point_list *, region_list *);
876 static void lialg_score_stroke(point_list *, point_list *, int *, int *);
877 static int lialg_compute_similarity(point_list *, point_list *);
878 static int lialg_compute_distance(point_list *, point_list *);
879
880 static int lialg_read_classifier_digest(rClassifier *);
881
882 static int lialg_canonicalize_examples(rClassifier *);
883 static int lialg_canonicalize_example_stroke(point_list *);
884 static int lialg_compute_equipoints(point_list *);
885
886 static int lialg_compute_pathlen(point_list *);
887 static int lialg_compute_pathlen_subset(point_list *, int, int);
888 static int lialg_filter_points(point_list *);
889 static int lialg_translate_points(point_list *, int, int, int, int);
890 static void lialg_get_bounding_box(point_list *, int *, int *, int *, int *);
891 static void lialg_compute_lpf_parameters();
892 static int isqrt(int);
893 static int likeatan(int, int);
894 static int quadr(int);
895
896
897 /*************************************************************
898
899 Core routines for the Li/Yeung recognition algorithm
900
901 *************************************************************/
902
lialg_initialize(rClassifier * rec)903 static void lialg_initialize(rClassifier *rec) {
904 int i;
905
906 /* Initialize the dompts arrays. */
907 for (i = 0; i < MAXSCLASSES; i++) {
908 rec->dompts[i] = nil;
909 }
910 }
911
912
913 /*
914 * Main recognition routine -- called by HRE API.
915 */
lialg_recognize_stroke(rClassifier * rec,point_list * stroke)916 static char *lialg_recognize_stroke(rClassifier *rec, point_list *stroke) {
917 int i;
918 char *name = nil;
919 point_list *input_dompts = nil;
920 char *best_name = nil;
921 int best_score = WORST_SCORE;
922 char *curr_name;
923 point_list *curr_dompts;
924
925 /* (void)gettimeofday(&stv, nil);*/
926
927 if (stroke->npts < 1) goto done;
928
929 /* Check for tap. */
930
931 /* First thing is to filter out ``close points.'' */
932 if (lialg_filter_points(stroke) != 0) return(nil);
933
934 /* Unfortunately, we don't have the actual time that each point */
935 /* was recorded (i.e., dt is invalid). Hence, we have to use a */
936 /* heuristic based on total distance and the number of points. */
937 if (stroke->npts == 1 || lialg_compute_pathlen(stroke) < TAP_PATHLEN)
938 return(TAP_CHAR);
939
940 /* Pre-process input stroke. */
941 if (lialg_preprocess_stroke(stroke) != 0) goto done;
942
943 /* Compute its dominant points. */
944 input_dompts = lialg_compute_dominant_points(stroke);
945 if (input_dompts == nil) goto done;
946
947 /* Score input stroke against every class in classifier. */
948 for (i = 0, curr_name = rec->cnames[i], curr_dompts = rec->dompts[i];
949 i < MAXSCLASSES && curr_name != nil && curr_dompts != nil;
950 i++, curr_name = rec->cnames[i], curr_dompts = rec->dompts[i]) {
951 int sim;
952 int dist;
953 int curr_score;
954
955 lialg_score_stroke(input_dompts, curr_dompts, &sim, &dist);
956 curr_score = dist;
957
958 if (lidebug && curr_score < DIST_THLD)
959 fprint(2, "(%s, %d, %d) ", curr_name, sim, dist);
960
961 /* Is it the best so far? */
962 if (curr_score < best_score && curr_score <= DIST_THLD) {
963 best_score = curr_score;
964 best_name = curr_name;
965 }
966 }
967
968 if (lidebug)
969 fprint(2, "\n");
970
971 /* No errors. */
972 name = best_name;
973
974 done:
975 delete_examples(input_dompts);
976 return(name);
977 }
978
979
lialg_preprocess_stroke(point_list * points)980 static int lialg_preprocess_stroke(point_list *points) {
981 int minx, miny, maxx, maxy, xrange, yrange, scale, xoff, yoff;
982
983 /* Filter out points that are too close. */
984 /* We did this earlier, when we checked for a tap. */
985 /*
986 if (lialg_filter_points(points) != 0) return(-1);
987 */
988
989 /* assert(points->npts > 0);*/
990
991 /* Scale up to avoid conversion errors. */
992 lialg_get_bounding_box(points, &minx, &miny, &maxx, &maxy);
993 xrange = maxx - minx;
994 yrange = maxy - miny;
995 scale = ( ((100 * xrange + CANONICAL_X / 2) / CANONICAL_X) >
996 ((100 * yrange + CANONICAL_Y / 2) / CANONICAL_Y))
997 ? (100 * CANONICAL_X + xrange / 2) / xrange
998 : (100 * CANONICAL_Y + yrange / 2) / yrange;
999 if (lialg_translate_points(points, minx, miny, scale, scale) != 0)
1000 return(-1);
1001
1002 /* Center the stroke. */
1003 lialg_get_bounding_box(points, &minx, &miny, &maxx, &maxy);
1004 xrange = maxx - minx;
1005 yrange = maxy - miny;
1006 xoff = -((CANONICAL_X - xrange + 1) / 2);
1007 yoff = -((CANONICAL_Y - yrange + 1) / 2);
1008 if (lialg_translate_points(points, xoff, yoff, 100, 100) != 0) return(-1);
1009
1010 /* Store the x and y ranges in the point list. */
1011 xrange = maxx - minx;
1012 yrange = maxy - miny;
1013 points->xrange = xrange;
1014 points->yrange = yrange;
1015
1016 if (lidebug) {
1017 int i;
1018 fprint(2, "After pre-processing: %d %d %d %d\n",
1019 minx, miny, maxx, maxy);
1020 for (i = 0; i < points->npts; i++)
1021 fprint(2, " (%P)\n", points->pts[i].Point);
1022 fflush(stderr);
1023 }
1024
1025 return(0);
1026 }
1027
1028
lialg_compute_dominant_points(point_list * points)1029 static point_list *lialg_compute_dominant_points(point_list *points) {
1030 point_list *ipts;
1031 region_list *regions;
1032 point_list *dpts;
1033
1034 /* Interpolate points. */
1035 ipts = lialg_interpolate_points(points);
1036 if (ipts == nil) return(nil);
1037 if (lidebug) {
1038 int j;
1039 fprint(2, "After interpolation: %d ipts\n", ipts->npts);
1040 for (j = 0; j < ipts->npts; j++) {
1041 fprint(2, " (%P), %lud\n", ipts->pts[j].Point, ipts->pts[j].chaincode);
1042 }
1043 fflush(stderr);
1044 }
1045
1046 /* Compute regions. */
1047 regions = lialg_compute_regions(ipts);
1048 /* assert(regions != nil);*/
1049
1050 /* Compute dominant points. */
1051 dpts = lialg_compute_dompts(ipts, regions);
1052 if (lidebug) {
1053 int j;
1054 fprint(2, "Dominant points: ");
1055 for (j = 0; j < dpts->npts; j++) {
1056 fprint(2, "%P (%lud) ", dpts->pts[j].Point, dpts->pts[j].chaincode);
1057 }
1058 fprint(2, "\n");
1059 fflush(stderr);
1060 }
1061
1062 /* Delete region data structure. */
1063 {
1064 region_list *curr, *next;
1065 for (curr = regions; curr != nil; ) {
1066 next = curr->next;
1067 free(curr);
1068 curr = next;
1069 }
1070 }
1071 delete_examples(ipts);
1072 return(dpts);
1073 }
1074
1075 /* Input points are assumed to be integer-valued! */
lialg_interpolate_points(point_list * points)1076 static point_list *lialg_interpolate_points(point_list *points) {
1077 int i, j;
1078 int maxpts;
1079 point_list *newpts;
1080
1081 /* Compute an upper-bound on the number of interpolated points. */
1082 maxpts = 0;
1083 for (i = 0; i < (points->npts - 1); i++) {
1084 pen_point *pta = &(points->pts[i]);
1085 pen_point *ptb = &(points->pts[i+1]);
1086 maxpts += abs(pta->x - ptb->x) + abs(pta->y - ptb->y);
1087 }
1088
1089 /* Allocate an array of the requisite size. */
1090 maxpts += points->npts;
1091 /* newpts = (point_list *)safe_malloc(sizeof(point_list)); */
1092 newpts = malloc(sizeof(point_list));
1093 newpts->pts = mallocz(maxpts*sizeof(pen_point), 1);
1094 if (newpts->pts == nil) {
1095 free(newpts);
1096 return(nil);
1097 }
1098 newpts->npts = maxpts;
1099 newpts->next = nil;
1100
1101 /* Interpolate each of the segments. */
1102 j = 0;
1103 for (i = 0; i < (points->npts - 1); i++) {
1104 pen_point *startpt = &(points->pts[i]);
1105 pen_point *endpt = &(points->pts[i+1]);
1106
1107 lialg_bresline(startpt, endpt, newpts, &j);
1108
1109 j--; /* end point gets recorded as start point of next segment! */
1110 }
1111
1112 /* Add-in last point. */
1113 newpts->pts[j++] = points->pts[points->npts - 1];
1114 newpts->npts = j;
1115
1116 /* Compute the chain code for P (the list of points). */
1117 lialg_compute_unit_chain_code(newpts);
1118
1119 return(newpts);
1120 }
1121
1122
1123 /* This implementation is due to Kenny Hoff. */
lialg_bresline(pen_point * startpt,pen_point * endpt,point_list * newpts,int * j)1124 static void lialg_bresline(pen_point *startpt, pen_point *endpt,
1125 point_list *newpts, int *j) {
1126 int Ax, Ay, Bx, By, dX, dY, Xincr, Yincr;
1127
1128 Ax = startpt->x;
1129 Ay = startpt->y;
1130 Bx = endpt->x;
1131 By = endpt->y;
1132
1133 /* INITIALIZE THE COMPONENTS OF THE ALGORITHM THAT ARE NOT AFFECTED */
1134 /* BY THE SLOPE OR DIRECTION OF THE LINE */
1135 dX = abs(Bx-Ax); /* store the change in X and Y of the line endpoints */
1136 dY = abs(By-Ay);
1137
1138 /* DETERMINE "DIRECTIONS" TO INCREMENT X AND Y (REGARDLESS OF DECISION) */
1139 if (Ax > Bx) { Xincr=-1; } else { Xincr=1; } /* which direction in X? */
1140 if (Ay > By) { Yincr=-1; } else { Yincr=1; } /* which direction in Y? */
1141
1142 /* DETERMINE INDEPENDENT VARIABLE (ONE THAT ALWAYS INCREMENTS BY 1 (OR -1) ) */
1143 /* AND INITIATE APPROPRIATE LINE DRAWING ROUTINE (BASED ON FIRST OCTANT */
1144 /* ALWAYS). THE X AND Y'S MAY BE FLIPPED IF Y IS THE INDEPENDENT VARIABLE. */
1145 if (dX >= dY) { /* if X is the independent variable */
1146 int dPr = dY<<1; /* amount to increment decision if right is chosen (always) */
1147 int dPru = dPr - (dX<<1); /* amount to increment decision if up is chosen */
1148 int P = dPr - dX; /* decision variable start value */
1149
1150 /* process each point in the line one at a time (just use dX) */
1151 for (; dX>=0; dX--) {
1152 newpts->pts[*j].x = Ax;
1153 newpts->pts[*j].y = Ay;
1154 (*j)++;
1155
1156 if (P > 0) { /* is the pixel going right AND up? */
1157 Ax+=Xincr; /* increment independent variable */
1158 Ay+=Yincr; /* increment dependent variable */
1159 P+=dPru; /* increment decision (for up) */
1160 } else { /* is the pixel just going right? */
1161 Ax+=Xincr; /* increment independent variable */
1162 P+=dPr; /* increment decision (for right) */
1163 }
1164 }
1165 } else { /* if Y is the independent variable */
1166 int dPr = dX<<1; /* amount to increment decision if right is chosen (always) */
1167 int dPru = dPr - (dY<<1); /* amount to increment decision if up is chosen */
1168 int P = dPr - dY; /* decision variable start value */
1169
1170 /* process each point in the line one at a time (just use dY) */
1171 for (; dY>=0; dY--) {
1172 newpts->pts[*j].x = Ax;
1173 newpts->pts[*j].y = Ay;
1174 (*j)++;
1175
1176 if (P > 0) { /* is the pixel going up AND right? */
1177 Ax+=Xincr; /* increment dependent variable */
1178 Ay+=Yincr; /* increment independent variable */
1179 P+=dPru; /* increment decision (for up) */
1180 } else { /* is the pixel just going up? */
1181 Ay+=Yincr; /* increment independent variable */
1182 P+=dPr; /* increment decision (for right) */
1183 }
1184 }
1185 }
1186 }
1187
lialg_compute_chain_code(point_list * pts)1188 static void lialg_compute_chain_code(point_list *pts) {
1189 int i;
1190
1191 for (i = 0; i < (pts->npts - 1); i++) {
1192 pen_point *startpt = &(pts->pts[i]);
1193 pen_point *endpt = &(pts->pts[i+1]);
1194 int dx = endpt->x - startpt->x;
1195 int dy = endpt->y - startpt->y;
1196 int tmp = quadr(likeatan(dy, dx));
1197 int dircode = (12 - tmp) % 8;
1198
1199 startpt->chaincode = dircode;
1200 }
1201 }
1202
1203
lialg_compute_unit_chain_code(point_list * pts)1204 static void lialg_compute_unit_chain_code(point_list *pts) {
1205 int i;
1206
1207 for (i = 0; i < (pts->npts - 1); i++) {
1208 pen_point *startpt = &(pts->pts[i]);
1209 pen_point *endpt = &(pts->pts[i+1]);
1210 int dx = endpt->x - startpt->x;
1211 int dy = endpt->y - startpt->y;
1212 int dircode = lialg_dctbl[dx+1][dy+1];
1213
1214 startpt->chaincode = dircode;
1215 }
1216 }
1217
1218
lialg_compute_regions(point_list * pts)1219 static region_list *lialg_compute_regions(point_list *pts) {
1220 region_list *regions;
1221 region_list *curr_reg;
1222 int *R[2 + LP_FILTER_ITERS];
1223 int *junk;
1224 int *curr, *next;
1225 int i, j;
1226
1227 /* Initialize low-pass filter parameters if necessary. */
1228 if (lialg_lpfconst == -1)
1229 lialg_compute_lpf_parameters();
1230
1231 /* Allocate a 2 x pts->npts array for use in computing the (filtered) Angle set, A_n. */
1232 junk = malloc((2 + LP_FILTER_ITERS) * pts->npts*sizeof(int));
1233 for (i = 0; i < (2 + LP_FILTER_ITERS); i++)
1234 R[i] = junk + (i * pts->npts);
1235 curr = R[0];
1236
1237 /* Compute the Angle set, A, in the first element of array R. */
1238 /* Values in R are in degrees, x 100. */
1239 curr[0] = 18000; /* a_0 */
1240 for (i = 1; i < (pts->npts - 1); i++) {
1241 int d_i = pts->pts[i].chaincode;
1242 int d_iminusone = pts->pts[i-1].chaincode;
1243 int a_i;
1244
1245 if (d_iminusone < d_i)
1246 d_iminusone += 8;
1247
1248 a_i = (d_iminusone - d_i) % 8;
1249
1250 /* convert to degrees, x 100 */
1251 curr[i] = ((12 - a_i) % 8) * 45 * 100;
1252 }
1253 curr[pts->npts - 1] = 18000; /* a_L-1 */
1254
1255 /* Perform a number of filtering iterations. */
1256 next = R[1];
1257 for (j = 0; j < LP_FILTER_ITERS; j++, curr = R[j], next = R[j+1]) {
1258 for (i = 0; i < pts->npts; i++) {
1259 int k;
1260
1261 next[i] = 0;
1262
1263 for (k = i - LP_FILTER_WIDTH; k <= i + LP_FILTER_WIDTH; k++) {
1264 int oldval = (k < 0 || k >= pts->npts) ? 18000 : curr[k];
1265 next[i] += oldval * lialg_lpfwts[k - (i - LP_FILTER_WIDTH)]; /* overflow? */
1266 }
1267
1268 next[i] /= lialg_lpfconst;
1269 }
1270 }
1271
1272 /* Do final thresholding around PI. */
1273 /* curr and next are set-up correctly at end of previous loop! */
1274 for (i = 0; i < pts->npts; i++)
1275 next[i] = (abs(curr[i] - 18000) < LP_FILTER_THLD) ? 18000 : curr[i];
1276 curr = next;
1277
1278 /* Debugging. */
1279 if (lidebug > 1) {
1280 for (i = 0; i < pts->npts; i++) {
1281 fprint(2, "%3d: (%P) %lud ",
1282 i, pts->pts[i].Point, pts->pts[i].chaincode);
1283 for (j = 0; j < 2 + LP_FILTER_ITERS; j++)
1284 fprint(2, "%d ", R[j][i]);
1285 fprint(2, "\n");
1286 }
1287 }
1288
1289 /* Do the region segmentation. */
1290 {
1291 int start, end;
1292 int currtype;
1293
1294 #define RGN_TYPE(val) (((val)==18000)?RGN_PLAIN:((val)<18000?RGN_CONCAVE:RGN_CONVEX))
1295
1296 start = 0;
1297 currtype = RGN_TYPE(curr[0]);
1298 regions = malloc(sizeof(region_list));
1299 curr_reg = regions;
1300 curr_reg->start = start;
1301 curr_reg->end = 0;
1302 curr_reg->type = currtype;
1303 curr_reg->next = nil;
1304 for (i = 1; i < pts->npts; i++) {
1305 int nexttype = RGN_TYPE(curr[i]);
1306
1307 if (nexttype != currtype) {
1308 region_list *next_reg;
1309
1310 end = i - 1;
1311 curr_reg->end = end;
1312 if (lidebug > 1)
1313 fprint(2, " (%d, %d), %d\n", start, end, currtype);
1314
1315 start = i;
1316 currtype = nexttype;
1317 next_reg = malloc(sizeof(region_list));
1318 next_reg->start = start;
1319 next_reg->end = 0;
1320 next_reg->type = nexttype;
1321 next_reg->next = nil;
1322
1323 curr_reg->next = next_reg;
1324 curr_reg = next_reg;
1325 }
1326 }
1327 end = i - 1;
1328 curr_reg->end = end;
1329 if (lidebug > 1)
1330 fprint(2, " (%d, %d), %d\n", start, end, currtype);
1331
1332 /* Filter out convex/concave regions that are too short. */
1333 for (curr_reg = regions; curr_reg; curr_reg = curr_reg->next)
1334 if (curr_reg->type == RGN_PLAIN) {
1335 region_list *next_reg;
1336
1337 for (next_reg = curr_reg->next;
1338 next_reg != nil &&
1339 (next_reg->end - next_reg->start) < LP_FILTER_MIN;
1340 next_reg = curr_reg->next) {
1341 /* next_reg must not be plain, and it must be followed by a plain */
1342 /* assert(next_reg->type != RGN_PLAIN); */
1343 /* assert(next_reg->next != nil && (next_reg->next)->type == RGN_PLAIN); */
1344
1345 curr_reg->next = (next_reg->next)->next;
1346 curr_reg->end = (next_reg->next)->end;
1347
1348 free(next_reg->next);
1349 free(next_reg);
1350 }
1351 }
1352
1353 /* Add-in pseudo-extremes. */
1354 {
1355 region_list *tmp, *prev_reg;
1356
1357 tmp = regions;
1358 regions = nil;
1359 prev_reg = nil;
1360 for (curr_reg = tmp; curr_reg; curr_reg = curr_reg->next) {
1361 if (curr_reg->type == RGN_PLAIN) {
1362 int arclen = lialg_compute_pathlen_subset(pts,
1363 curr_reg->start,
1364 curr_reg->end);
1365 int dx = pts->pts[curr_reg->end].x -
1366 pts->pts[curr_reg->start].x;
1367 int dy = pts->pts[curr_reg->end].y -
1368 pts->pts[curr_reg->start].y;
1369 int chordlen = isqrt(10000 * (dx * dx + dy * dy));
1370 int atcr = (chordlen == 0) ? 0 : (100 * arclen + chordlen / 2) / chordlen;
1371
1372 if (lidebug)
1373 fprint(2, "%d, %d, %d\n", arclen, chordlen, atcr);
1374
1375 /* Split region if necessary. */
1376 if (arclen >= PE_AL_THLD && atcr >= PE_ATCR_THLD) {
1377 int mid = curr_reg->start + (curr_reg->end - curr_reg->start) / 2;
1378 int end = curr_reg->end;
1379 region_list *saved_next = curr_reg->next;
1380
1381 curr_reg->end = mid - 1;
1382 if (prev_reg == nil)
1383 regions = curr_reg;
1384 else
1385 prev_reg->next = curr_reg;
1386 prev_reg = curr_reg;
1387
1388 /* curr_reg = (region_list *)safe_malloc(sizeof(region_list));*/
1389 curr_reg = malloc(sizeof(region_list));
1390 curr_reg->start = mid;
1391 curr_reg->end = mid;
1392 curr_reg->type = RGN_PSEUDO;
1393 curr_reg->next = nil;
1394 prev_reg->next = curr_reg;
1395 prev_reg = curr_reg;
1396
1397 /* curr_reg = (region_list *)malloc(sizeof(region_list)); */
1398 curr_reg = malloc(sizeof(region_list));
1399 curr_reg->start = mid + 1;
1400 curr_reg->end = end;
1401 curr_reg->type = RGN_PLAIN;
1402 curr_reg->next = nil;
1403 prev_reg->next = curr_reg;
1404 prev_reg = curr_reg;
1405
1406 curr_reg->next = saved_next;
1407 continue;
1408 }
1409 }
1410
1411 if (prev_reg == nil)
1412 regions = curr_reg;
1413 else
1414 prev_reg->next = curr_reg;
1415 prev_reg = curr_reg;
1416 }
1417 }
1418 }
1419
1420 free(junk);
1421 return(regions);
1422 }
1423
1424
lialg_compute_dompts(point_list * pts,region_list * regions)1425 static point_list *lialg_compute_dompts(point_list *pts, region_list *regions) {
1426 point_list *dpts;
1427 int ndpts;
1428 int *cas;
1429 int nonplain;
1430 region_list *r;
1431 region_list *curr;
1432 int dp;
1433 int previx;
1434 int currix;
1435
1436 /* Compute contour angle set. */
1437 cas = lialg_compute_contour_angle_set(pts, regions);
1438
1439 /* Dominant points include: start_pt, end_pt, extrema_of_non_plain_regions, midpts of the preceding. */
1440 nonplain = 0;
1441 for (r = regions; r != nil; r = r->next)
1442 if (r->type != RGN_PLAIN)
1443 nonplain++;
1444 ndpts = 2 * (2 + nonplain) - 1;
1445 /* dpts = (point_list *)safe_malloc(sizeof(point_list)); */
1446 dpts = malloc(sizeof(point_list));
1447 dpts->pts = mallocz(ndpts*sizeof(pen_point), 1);
1448 if (dpts->pts == nil) {
1449 free(dpts);
1450 return(nil);
1451 }
1452 dpts->npts = ndpts;
1453 dpts->next = nil;
1454
1455 /* Pick out dominant points. */
1456
1457 /* Record start point. */
1458 dp = 0;
1459 previx = 0;
1460 dpts->pts[dp++] = pts->pts[previx];
1461
1462 for (curr = regions; curr != nil; curr = curr->next)
1463 if (curr->type != RGN_PLAIN) {
1464 int max_v = 0;
1465 int min_v = 0x7fffffff; /* maxint */
1466 int max_ix = -1;
1467 int min_ix = -1;
1468 int i;
1469
1470 for (i = curr->start; i <= curr->end; i++) {
1471 int v = cas[i];
1472 if (v > max_v) { max_v = v; max_ix = i; }
1473 if (v < min_v) { min_v = v; min_ix = i; }
1474 if (lidebug > 1)
1475 fprint(2, " %d\n", v);
1476 }
1477
1478 currix = (curr->type == RGN_CONVEX ? max_ix : min_ix);
1479
1480 /* Record midpoint. */
1481 dpts->pts[dp++] = pts->pts[previx + (currix - previx) / 2];
1482
1483 /* Record extreme point. */
1484 dpts->pts[dp++] = pts->pts[currix];
1485
1486 previx = currix;
1487 }
1488
1489 /* Record last mid-point and end point. */
1490 currix = pts->npts - 1;
1491 dpts->pts[dp++] = pts->pts[previx + (currix - previx) / 2];
1492 dpts->pts[dp] = pts->pts[currix];
1493
1494 /* Compute chain-code. */
1495 lialg_compute_chain_code(dpts);
1496
1497 free(cas);
1498 return(dpts);
1499 }
1500
1501
lialg_compute_contour_angle_set(point_list * pts,region_list * regions)1502 static int *lialg_compute_contour_angle_set(point_list *pts,
1503 region_list *regions) {
1504 int *V;
1505 region_list *curr_reg;
1506 int i;
1507
1508 V = malloc(pts->npts*sizeof(int));
1509
1510 V[0] = 18000;
1511 for (curr_reg = regions; curr_reg != nil; curr_reg = curr_reg->next) {
1512 for (i = curr_reg->start; i <= curr_reg->end; i++) {
1513 if (curr_reg->type == RGN_PLAIN) {
1514 V[i] = 18000;
1515 } else {
1516 /* For now, simply choose the mid-point. */
1517 int isMidPt = i == (curr_reg->start +
1518 (curr_reg->end - curr_reg->start) / 2);
1519 V[i] = (curr_reg->type == RGN_CONVEX)
1520 ? (isMidPt ? 18000 : 0)
1521 : (isMidPt ? 0 : 18000);
1522 }
1523 }
1524 }
1525 V[pts->npts - 1] = 18000;
1526
1527 return(V);
1528 }
1529
1530
1531 /*
1532 * First compute the similarity between the two strings.
1533 * If it's above a threshold, compute the distance between
1534 * the two and return it as the ``score.''
1535 * Otherwise, return the constant WORST_SCORE.
1536 *
1537 */
lialg_score_stroke(point_list * input_dompts,point_list * curr_dompts,int * sim,int * dist)1538 static void lialg_score_stroke(point_list *input_dompts, point_list *curr_dompts, int *sim, int *dist) {
1539 *sim = MIN_SIM;
1540 *dist = MAX_DIST;
1541
1542 *sim = lialg_compute_similarity(input_dompts, curr_dompts);
1543 if (*sim < SIM_THLD) goto done;
1544
1545 *dist = lialg_compute_distance(input_dompts, curr_dompts);
1546
1547 done:
1548 if (lidebug)
1549 fprint(2, "%d, %d\n", *sim, *dist);
1550 }
1551
1552
lialg_compute_similarity(point_list * input_dompts,point_list * curr_dompts)1553 static int lialg_compute_similarity(point_list *input_dompts, point_list *curr_dompts) {
1554 int sim;
1555 point_list *A, *B;
1556 int N, M;
1557 int **G;
1558 int *junk;
1559 int i, j;
1560
1561 /* A is the longer sequence, length N. */
1562 /* B is the shorter sequence, length M. */
1563 if (input_dompts->npts >= curr_dompts->npts) {
1564 A = input_dompts;
1565 N = input_dompts->npts;
1566 B = curr_dompts;
1567 M = curr_dompts->npts;
1568 } else {
1569 A = curr_dompts;
1570 N = curr_dompts->npts;
1571 B = input_dompts;
1572 M = input_dompts->npts;
1573 }
1574
1575 /* Allocate and initialize the Gain matrix, G. */
1576 /* The size of G is M x (N + 1). */
1577 /* Note that row 0 is unused. */
1578 /* Similarities are x 10. */
1579 G = malloc(M*sizeof(int *));
1580 junk = malloc(M * (N + 1) * sizeof(int));
1581 for (i = 0; i < M; i++)
1582 G[i] = junk + (i * (N + 1));
1583
1584 for (i = 1; i < M; i++) {
1585 int bval = B->pts[i-1].chaincode;
1586
1587 /* Source column. */
1588 G[i][0] = 0;
1589
1590 for (j = 1; j < N; j++) {
1591 int aval = A->pts[j-1].chaincode;
1592 int diff = abs(bval - aval);
1593 if (diff > 4) diff = 8 - diff;
1594
1595 G[i][j] = (diff == 0)
1596 ? 10
1597 : (diff == 1)
1598 ? 6
1599 : 0;
1600 }
1601
1602 /* Sink column. */
1603 G[i][N] = 0;
1604 }
1605
1606 /* Do the DP algorithm. */
1607 /* Proceed in column order, from highest column to the lowest. */
1608 /* Within each column, proceed from the highest row to the lowest. */
1609 /* Skip the highest column. */
1610 for (j = N - 1; j >= 0; j--)
1611 for (i = M - 1; i > 0; i--) {
1612 int max = G[i][j + 1];
1613
1614 if (i < (M - 1)) {
1615 int tmp = G[i + 1][j + 1];
1616 if (tmp > max) max = tmp;
1617 }
1618
1619 G[i][j] += max;
1620 }
1621
1622 sim = (10 * G[1][0] + (N - 1) / 2) / (N - 1);
1623
1624 if (G != nil)
1625 free(G);
1626 if (junk != nil)
1627 free(junk);
1628 return(sim);
1629 }
1630
1631
lialg_compute_distance(point_list * input_dompts,point_list * curr_dompts)1632 static int lialg_compute_distance(point_list *input_dompts,
1633 point_list *curr_dompts) {
1634 int dist;
1635 point_list *A, *B;
1636 int N, M;
1637 int **C;
1638 int *junk;
1639 int *BE;
1640 int *TE;
1641 int i, j;
1642
1643 /* A is the longer sequence, length N. */
1644 /* B is the shorter sequence, length M. */
1645 if (input_dompts->npts >= curr_dompts->npts) {
1646 A = input_dompts;
1647 N = input_dompts->npts;
1648 B = curr_dompts;
1649 M = curr_dompts->npts;
1650 }
1651 else {
1652 A = curr_dompts;
1653 N = curr_dompts->npts;
1654 B = input_dompts;
1655 M = input_dompts->npts;
1656 }
1657
1658 /* Construct the helper vectors, BE and TE, which say for each column */
1659 /* what are the ``bottom'' and ``top'' rows of interest. */
1660 BE = malloc((N + 1)*sizeof(int));
1661 TE = malloc((N + 1)*sizeof(int));
1662
1663 for (j = 1; j <= N; j++) {
1664 int bot, top;
1665
1666 bot = j + (M - DP_BAND);
1667 if (bot > M) bot = M;
1668 BE[j] = bot;
1669
1670 top = j - (N - DP_BAND);
1671 if (top < 1) top = 1;
1672 TE[j] = top;
1673 }
1674
1675 /* Allocate and initialize the Cost matrix, C. */
1676 /* The size of C is (M + 1) x (N + 1). */
1677 /* Note that row and column 0 are unused. */
1678 /* Costs are x 100. */
1679 /* C = (int **)safe_malloc((M + 1) * sizeof(int *)); */
1680 C = malloc((M + 1)*sizeof( int *));
1681 junk = malloc((M + 1) * (N + 1)*sizeof(int));
1682 for (i = 0; i <= M; i++)
1683 C[i] = junk + (i * (N + 1));
1684
1685 for (i = 1; i <= M; i++) {
1686 int bx = B->pts[i-1].x;
1687 int by = B->pts[i-1].y;
1688
1689 for (j = 1; j <= N; j++) {
1690 int ax = A->pts[j-1].x;
1691 int ay = A->pts[j-1].y;
1692 int dx = bx - ax;
1693 int dy = by - ay;
1694 int dist = isqrt(10000 * (dx * dx + dy * dy));
1695
1696 C[i][j] = dist;
1697 }
1698 }
1699
1700 /* Do the DP algorithm. */
1701 /* Proceed in column order, from highest column to the lowest. */
1702 /* Within each column, proceed from the highest row to the lowest. */
1703 for (j = N; j > 0; j--)
1704 for (i = M; i > 0; i--) {
1705 int min = MAX_DIST;
1706
1707 if (i > BE[j] || i < TE[j] || (j == N && i == M))
1708 continue;
1709
1710 if (j < N) {
1711 if (i >= TE[j+1]) {
1712 int tmp = C[i][j+1];
1713 if (tmp < min)
1714 min = tmp;
1715 }
1716
1717 if (i < M) {
1718 int tmp = C[i+1][j+1];
1719 if (tmp < min)
1720 min = tmp;
1721 }
1722 }
1723
1724 if (i < BE[j]) {
1725 int tmp = C[i+1][j];
1726 if (tmp < min) min = tmp;
1727 }
1728
1729 C[i][j] += min;
1730 }
1731
1732 dist = (C[1][1] + N / 2) / N;
1733
1734 if (C != nil) free(C);
1735 if (junk != nil) free(junk);
1736 if (BE != nil) free(BE);
1737 if (TE != nil) free(TE);
1738 return(dist);
1739 }
1740
1741
1742 /*************************************************************
1743
1744 Digest-processing routines
1745
1746 *************************************************************/
1747
lialg_read_classifier_digest(rClassifier * rec)1748 static int lialg_read_classifier_digest(rClassifier *rec) {
1749 int nclasses;
1750 FILE *fp;
1751
1752 /* Try to open the corresponding digest file. */
1753 {
1754 char *clx_path;
1755 char *dot;
1756
1757 /* Get a copy of the filename, with some room on the end. */
1758 /* clx_path = safe_malloc(strlen(rec->file_name) + 5); */
1759 clx_path = malloc((strlen(rec->file_name) + 5) *sizeof(char));
1760 strcpy(clx_path, rec->file_name);
1761
1762 /* Truncate the path after the last dot. */
1763 dot = strrchr(clx_path, '.');
1764 if (dot == nil) { free(clx_path); return(-1); }
1765 *(dot + 1) = 0;
1766
1767 /* Append the classifier-digest extension. */
1768 strcat(clx_path, "clx");
1769
1770 fp = fopen(clx_path, "r");
1771 if (fp == nil) {
1772 free(clx_path);
1773 return(-1);
1774 }
1775
1776 free(clx_path);
1777 }
1778
1779 /* Read-in the name and dominant points for each class. */
1780 for (nclasses = 0; !feof(fp); nclasses++) {
1781 point_list *dpts = nil;
1782 char class[BUFSIZ];
1783 int npts;
1784 int j;
1785
1786 if (fscanf(fp, "%s %d", class, &npts) != 2) {
1787 if (feof(fp)) break;
1788
1789 goto failed;
1790 }
1791 rec->cnames[nclasses] = strdup(class);
1792
1793 /* Allocate a dominant-points list. */
1794 /* dpts = (point_list *)safe_malloc(sizeof(point_list)); */
1795 dpts = malloc(sizeof(point_list));
1796 dpts->pts = mallocz(npts*sizeof(pen_point), 1);
1797 if (dpts->pts == nil) goto failed;
1798 dpts->npts = npts;
1799 dpts->next = nil;
1800
1801 /* Read in each point. */
1802 for (j = 0; j < npts; j++) {
1803 int x, y;
1804
1805 if (fscanf(fp, "%d %d", &x, &y) != 2) goto failed;
1806 dpts->pts[j].x = x;
1807 dpts->pts[j].y = y;
1808 }
1809
1810 /* Compute the chain-code. */
1811 lialg_compute_chain_code(dpts);
1812
1813 /* Store the list in the rec data structure. */
1814 rec->dompts[nclasses] = dpts;
1815
1816 continue;
1817
1818 failed:
1819 fprint(2, "read_classifier_digest failed...\n");
1820 for (; nclasses >= 0; nclasses--) {
1821 if (rec->cnames[nclasses] != nil) {
1822 free(rec->cnames[nclasses]);
1823 rec->cnames[nclasses] = nil;
1824 }
1825 if (rec->dompts[nclasses] != nil) {
1826 delete_examples(rec->dompts[nclasses]);
1827 rec->dompts[nclasses] = nil;
1828 }
1829 }
1830 if (dpts != nil)
1831 delete_examples(dpts);
1832 fclose(fp);
1833 return(-1);
1834 }
1835
1836 fclose(fp);
1837 return(0);
1838 }
1839
1840
1841 /*************************************************************
1842
1843 Canonicalization routines
1844
1845 *************************************************************/
1846
lialg_canonicalize_examples(rClassifier * rec)1847 static int lialg_canonicalize_examples(rClassifier *rec) {
1848 int i;
1849 int nclasses;
1850
1851 if (lidebug) {
1852 fprint(2, "lialg_canonicalize_examples working on %s\n",
1853 rec->file_name);
1854 }
1855 /* Initialize canonical-example arrays. */
1856 for (i = 0; i < MAXSCLASSES; i++) {
1857 rec->canonex[i] = nil;
1858 }
1859
1860 /* Figure out number of classes. */
1861 for (nclasses = 0;
1862 nclasses < MAXSCLASSES && rec->cnames[nclasses] != nil;
1863 nclasses++)
1864 ;
1865
1866 /* Canonicalize the examples for each class. */
1867 for (i = 0; i < nclasses; i++) {
1868 int j, k;
1869 int nex;
1870 point_list *pts, *tmp, *avg;
1871 int maxxrange, maxyrange;
1872 int minx, miny, maxx, maxy;
1873 int avgxrange, avgyrange, avgxoff, avgyoff, avgscale;
1874
1875
1876 if (lidebug) {
1877 fprint(2, "lialg_canonicalize_examples working on class %s\n",
1878 rec->cnames[i]);
1879 }
1880 /* Make a copy of the examples. */
1881 pts = nil;
1882 tmp = rec->ex[i];
1883 for (nex = 0; tmp != nil; nex++, tmp = tmp->next) {
1884 if ((pts = add_example(pts, tmp->npts, tmp->pts)) == nil) {
1885 delete_examples(pts);
1886 return(-1);
1887 }
1888 }
1889
1890 /* Canonicalize each example, and derive the max x and y ranges. */
1891 maxxrange = 0;
1892 maxyrange = 0;
1893 for (j = 0, tmp = pts; j < nex; j++, tmp = tmp->next) {
1894 if (lialg_canonicalize_example_stroke(tmp) != 0) {
1895 if (lidebug) {
1896 fprint(2, "lialg_canonicalize_example_stroke returned error\n");
1897 }
1898 return(-1);
1899 }
1900
1901 if (tmp->xrange > maxxrange) maxxrange = tmp->xrange;
1902 if (tmp->yrange > maxyrange) maxyrange = tmp->yrange;
1903 }
1904
1905 /* Normalize max ranges. */
1906 if (((100 * maxxrange + CANONICAL_X / 2) / CANONICAL_X) >
1907 ((100 * maxyrange + CANONICAL_Y / 2) / CANONICAL_Y)) {
1908 maxyrange = (maxyrange * CANONICAL_X + maxxrange / 2) / maxxrange;
1909 maxxrange = CANONICAL_X;
1910 }
1911 else {
1912 maxxrange = (maxxrange * CANONICAL_Y + maxyrange / 2) / maxyrange;
1913 maxyrange = CANONICAL_Y;
1914 }
1915
1916 /* Re-scale each example to max ranges. */
1917 for (j = 0, tmp = pts; j < nex; j++, tmp = tmp->next) {
1918 int scalex = (tmp->xrange == 0) ? 100 : (100 * maxxrange + tmp->xrange / 2) / tmp->xrange;
1919 int scaley = (tmp->yrange == 0) ? 100 : (100 * maxyrange + tmp->yrange / 2) / tmp->yrange;
1920 if (lialg_translate_points(tmp, 0, 0, scalex, scaley) != 0) {
1921 delete_examples(pts);
1922 return(-1);
1923 }
1924 }
1925
1926 /* Average the examples; leave average in first example. */
1927 avg = pts; /* careful aliasing!! */
1928 for (k = 0; k < NCANONICAL; k++) {
1929 int xsum = 0;
1930 int ysum = 0;
1931
1932 for (j = 0, tmp = pts; j < nex; j++, tmp = tmp->next) {
1933 xsum += tmp->pts[k].x;
1934 ysum += tmp->pts[k].y;
1935 }
1936
1937 avg->pts[k].x = (xsum + j / 2) / j;
1938 avg->pts[k].y = (ysum + j / 2) / j;
1939 }
1940
1941 /* Compute BB of averaged stroke and re-scale. */
1942 lialg_get_bounding_box(avg, &minx, &miny, &maxx, &maxy);
1943 avgxrange = maxx - minx;
1944 avgyrange = maxy - miny;
1945 avgscale = (((100 * avgxrange + CANONICAL_X / 2) / CANONICAL_X) >
1946 ((100 * avgyrange + CANONICAL_Y / 2) / CANONICAL_Y))
1947 ? (100 * CANONICAL_X + avgxrange / 2) / avgxrange
1948 : (100 * CANONICAL_Y + avgyrange / 2) / avgyrange;
1949 if (lialg_translate_points(avg, minx, miny, avgscale, avgscale) != 0) {
1950 delete_examples(pts);
1951 return(-1);
1952 }
1953
1954 /* Re-compute the x and y ranges and center the stroke. */
1955 lialg_get_bounding_box(avg, &minx, &miny, &maxx, &maxy);
1956 avgxrange = maxx - minx;
1957 avgyrange = maxy - miny;
1958 avgxoff = -((CANONICAL_X - avgxrange + 1) / 2);
1959 avgyoff = -((CANONICAL_Y - avgyrange + 1) / 2);
1960 if (lialg_translate_points(avg, avgxoff, avgyoff, 100, 100) != 0) {
1961 delete_examples(pts);
1962 return(-1);
1963 }
1964
1965 /* Create a point list to serve as the ``canonical representation. */
1966 if ((rec->canonex[i] = add_example(nil, avg->npts, avg->pts)) == nil) {
1967 delete_examples(pts);
1968 return(-1);
1969 }
1970 (rec->canonex[i])->xrange = maxx - minx;
1971 (rec->canonex[i])->yrange = maxy - miny;
1972
1973 if (lidebug) {
1974 fprint(2, "%s, avgpts = %d\n", rec->cnames[i], avg->npts);
1975 for (j = 0; j < avg->npts; j++) {
1976 fprint(2, " (%P)\n", avg->pts[j].Point);
1977 }
1978 }
1979
1980 /* Compute dominant points of canonical representation. */
1981 rec->dompts[i] = lialg_compute_dominant_points(avg);
1982
1983 /* Clean up. */
1984 delete_examples(pts);
1985 }
1986
1987 /* Sanity check. */
1988 for (i = 0; i < nclasses; i++) {
1989 char *best_name = lialg_recognize_stroke(rec, rec->canonex[i]);
1990
1991 if (best_name != rec->cnames[i])
1992 fprint(2, "%s, best = %s\n", rec->cnames[i], best_name);
1993 }
1994
1995 return(0);
1996 }
1997
1998
lialg_canonicalize_example_stroke(point_list * points)1999 static int lialg_canonicalize_example_stroke(point_list *points) {
2000 int minx, miny, maxx, maxy, xrange, yrange, scale;
2001
2002 /* Filter out points that are too close. */
2003 if (lialg_filter_points(points) != 0) return(-1);
2004
2005 /* Must be at least two points! */
2006 if (points->npts < 2) {
2007 if (lidebug) {
2008 fprint(2, "lialg_canonicalize_example_stroke: npts=%d\n",
2009 points->npts);
2010 }
2011 return(-1);
2012 }
2013
2014 /* Scale up to avoid conversion errors. */
2015 lialg_get_bounding_box(points, &minx, &miny, &maxx, &maxy);
2016 xrange = maxx - minx;
2017 yrange = maxy - miny;
2018 scale = (((100 * xrange + CANONICAL_X / 2) / CANONICAL_X) >
2019 ((100 * yrange + CANONICAL_Y / 2) / CANONICAL_Y))
2020 ? (100 * CANONICAL_X + xrange / 2) / xrange
2021 : (100 * CANONICAL_Y + yrange / 2) / yrange;
2022 if (lialg_translate_points(points, minx, miny, scale, scale) != 0) {
2023 if (lidebug) {
2024 fprint(2, "lialg_translate_points (minx=%d,miny=%d,scale=%d) returned error\n", minx, miny, scale);
2025 }
2026 return(-1);
2027 }
2028
2029 /* Compute an equivalent stroke with equi-distant points. */
2030 if (lialg_compute_equipoints(points) != 0) return(-1);
2031
2032 /* Re-translate the points to the origin. */
2033 lialg_get_bounding_box(points, &minx, &miny, &maxx, &maxy);
2034 if (lialg_translate_points(points, minx, miny, 100, 100) != 0) {
2035 if (lidebug) {
2036 fprint(2, "lialg_translate_points (minx=%d,miny=%d) returned error\n", minx, miny);
2037 }
2038 return(-1);
2039 }
2040
2041 /* Store the x and y ranges in the point list. */
2042 xrange = maxx - minx;
2043 yrange = maxy - miny;
2044 points->xrange = xrange;
2045 points->yrange = yrange;
2046
2047 if (lidebug) {
2048 int i;
2049 fprint(2, "Canonicalized: %d, %d, %d, %d\n", minx, miny, maxx, maxy);
2050 for (i = 0; i < points->npts; i++)
2051 fprint(2, " (%P)\n", points->pts[i].Point);
2052 fflush(stderr);
2053 }
2054
2055 return(0);
2056 }
2057
2058
lialg_compute_equipoints(point_list * points)2059 static int lialg_compute_equipoints(point_list *points) {
2060 pen_point *equipoints = mallocz(NCANONICAL*sizeof(pen_point), 1);
2061 int nequipoints = 0;
2062 int pathlen = lialg_compute_pathlen(points);
2063 int equidist = (pathlen + (NCANONICAL - 1) / 2) / (NCANONICAL - 1);
2064 int i;
2065 int dist_since_last_eqpt;
2066 int remaining_seglen;
2067 int dist_to_next_eqpt;
2068
2069 if (equipoints == nil) {
2070 fprint(2, "can't allocate memory in lialg_compute_equipoints");
2071 return(-1);
2072 }
2073
2074 if (lidebug) {
2075 fprint(2, "compute_equipoints: npts = %d, pathlen = %d, equidist = %d\n",
2076 points->npts, pathlen, equidist);
2077 fflush(stderr);
2078 }
2079
2080 /* First original point is an equipoint. */
2081 equipoints[0] = points->pts[0];
2082 nequipoints++;
2083 dist_since_last_eqpt = 0;
2084
2085 for (i = 1; i < points->npts; i++) {
2086 int dx1 = points->pts[i].x - points->pts[i-1].x;
2087 int dy1 = points->pts[i].y - points->pts[i-1].y;
2088 int endx = 100 * points->pts[i-1].x;
2089 int endy = 100 * points->pts[i-1].y;
2090 remaining_seglen = isqrt(10000 * (dx1 * dx1 + dy1 * dy1));
2091 dist_to_next_eqpt = equidist - dist_since_last_eqpt;
2092
2093 while (remaining_seglen >= dist_to_next_eqpt) {
2094 if (dx1 == 0) {
2095 /* x-coordinate stays the same */
2096 if (dy1 >= 0)
2097 endy += dist_to_next_eqpt;
2098 else
2099 endy -= dist_to_next_eqpt;
2100 }
2101 else {
2102 int slope = (100 * dy1 + dx1 / 2) / dx1;
2103 int tmp = isqrt(10000 + slope * slope);
2104 int dx = (100 * dist_to_next_eqpt + tmp / 2) / tmp;
2105 int dy = (slope * dx + 50) / 100;
2106
2107 if (dy < 0) dy = -dy;
2108 if (dx1 >= 0)
2109 endx += dx;
2110 else
2111 endx -= dx;
2112 if (dy1 >= 0)
2113 endy += dy;
2114 else
2115 endy -= dy;
2116 }
2117
2118 equipoints[nequipoints].x = (endx + 50) / 100;
2119 equipoints[nequipoints].y = (endy + 50) / 100;
2120 nequipoints++;
2121 /* assert(nequipoints <= NCANONICAL);*/
2122 dist_since_last_eqpt = 0;
2123 remaining_seglen -= dist_to_next_eqpt;
2124 dist_to_next_eqpt = equidist;
2125 }
2126
2127 dist_since_last_eqpt += remaining_seglen;
2128 }
2129
2130 /* Take care of last equipoint. */
2131 if (nequipoints == NCANONICAL) {
2132 /* Good. */
2133 } else if (nequipoints == (NCANONICAL - 1)) {
2134 /* Make last original point the last equipoint. */
2135 equipoints[nequipoints] = points->pts[points->npts - 1];
2136 } else {
2137 if (lidebug) {
2138 fprint(2,"lialg_compute_equipoints: nequipoints = %d\n",
2139 nequipoints);
2140 }
2141 /* assert(false);*/
2142 return(-1);
2143 }
2144
2145 points->npts = NCANONICAL;
2146 delete_pen_point_array(points->pts);
2147 points->pts = equipoints;
2148 return(0);
2149 }
2150
2151
2152 /*************************************************************
2153
2154 Utility routines
2155
2156 *************************************************************/
2157
2158 /* Result is x 100. */
lialg_compute_pathlen(point_list * points)2159 static int lialg_compute_pathlen(point_list *points) {
2160 return(lialg_compute_pathlen_subset(points, 0, points->npts - 1));
2161 }
2162
2163
2164 /* Result is x 100. */
lialg_compute_pathlen_subset(point_list * points,int start,int end)2165 static int lialg_compute_pathlen_subset(point_list *points,
2166 int start, int end) {
2167 int pathlen;
2168 int i;
2169
2170 pathlen = 0;
2171 for (i = start + 1; i <= end; i++) {
2172 int dx = points->pts[i].x - points->pts[i-1].x;
2173 int dy = points->pts[i].y - points->pts[i-1].y;
2174 int dist = isqrt(10000 * (dx * dx + dy * dy));
2175 pathlen += dist;
2176 }
2177
2178 return(pathlen);
2179 }
2180
2181
2182 /* Note that this does NOT update points->xrange and points->yrange! */
lialg_filter_points(point_list * points)2183 static int lialg_filter_points(point_list *points) {
2184 int filtered_npts;
2185 pen_point *filtered_pts = mallocz(points->npts*sizeof(pen_point), 1);
2186 int i;
2187
2188 if (filtered_pts == nil) {
2189 fprint(2, "can't allocate memory in lialg_filter_points");
2190 return(-1);
2191 }
2192
2193 filtered_pts[0] = points->pts[0];
2194 filtered_npts = 1;
2195 for (i = 1; i < points->npts; i++) {
2196 int j = filtered_npts - 1;
2197 int dx = points->pts[i].x - filtered_pts[j].x;
2198 int dy = points->pts[i].y - filtered_pts[j].y;
2199 int magsq = dx * dx + dy * dy;
2200
2201 if (magsq >= DIST_SQ_THRESHOLD) {
2202 filtered_pts[filtered_npts] = points->pts[i];
2203 filtered_npts++;
2204 }
2205 }
2206
2207 points->npts = filtered_npts;
2208 delete_pen_point_array(points->pts);
2209 points->pts = filtered_pts;
2210 return(0);
2211 }
2212
2213
2214 /* scalex and scaley are x 100. */
2215 /* Note that this does NOT update points->xrange and points->yrange! */
lialg_translate_points(point_list * points,int minx,int miny,int scalex,int scaley)2216 static int lialg_translate_points(point_list *points,
2217 int minx, int miny,
2218 int scalex, int scaley) {
2219 int i;
2220
2221 for (i = 0; i < points->npts; i++) {
2222 points->pts[i].x = ((points->pts[i].x - minx) * scalex + 50) / 100;
2223 points->pts[i].y = ((points->pts[i].y - miny) * scaley + 50) / 100;
2224 }
2225
2226 return(0);
2227 }
2228
2229
lialg_get_bounding_box(point_list * points,int * pminx,int * pminy,int * pmaxx,int * pmaxy)2230 static void lialg_get_bounding_box(point_list *points,
2231 int *pminx, int *pminy,
2232 int *pmaxx, int *pmaxy) {
2233 int minx, miny, maxx, maxy;
2234 int i;
2235
2236 minx = maxx = points->pts[0].x;
2237 miny = maxy = points->pts[0].y;
2238 for (i = 1; i < points->npts; i++) {
2239 pen_point *pt = &(points->pts[i]);
2240 if (pt->x < minx) minx = pt->x;
2241 else if (pt->x > maxx) maxx = pt->x;
2242 if (pt->y < miny) miny = pt->y;
2243 else if (pt->y > maxy) maxy = pt->y;
2244 }
2245
2246 *pminx = minx;
2247 *pminy = miny;
2248 *pmaxx = maxx;
2249 *pmaxy = maxy;
2250 }
2251
2252
2253 int wtvals[] = {100, 104, 117, 143, 189, 271, 422};
2254
lialg_compute_lpf_parameters(void)2255 static void lialg_compute_lpf_parameters(void) {
2256 int i;
2257
2258 for (i = LP_FILTER_WIDTH; i >= 0; i--) {
2259 // double x = 0.04 * (i * i);
2260 // double tmp = 100.0 * exp(x);
2261 // int wt = floor((double)tmp);
2262 int wt = wtvals[i];
2263 lialg_lpfwts[LP_FILTER_WIDTH - i] = wt;
2264 lialg_lpfwts[LP_FILTER_WIDTH + i] = wt;
2265 }
2266 lialg_lpfconst = 0;
2267 for (i = 0; i < (2 * LP_FILTER_WIDTH + 1); i++) {
2268 lialg_lpfconst += lialg_lpfwts[i];
2269 }
2270 }
2271
2272
2273 /* Code from Joseph Hall (jnhall@sat.mot.com). */
isqrt(int n)2274 static int isqrt(int n) {
2275 register int i;
2276 register long k0, k1, nn;
2277
2278 for (nn = i = n, k0 = 2; i > 0; i >>= 2, k0 <<= 1)
2279 ;
2280 nn <<= 2;
2281 for (;;) {
2282 k1 = (nn / k0 + k0) >> 1;
2283 if (((k0 ^ k1) & ~1) == 0)
2284 break;
2285 k0 = k1;
2286 }
2287 return (int) ((k1 + 1) >> 1);
2288 }
2289
2290
2291 /* Helper routines from Mark Hayter. */
likeatan(int tantop,int tanbot)2292 static int likeatan(int tantop, int tanbot) {
2293 int t;
2294 /* Use tan(theta)=top/bot --> order for t */
2295 /* t in range 0..0x40000 */
2296
2297 if ((tantop == 0) && (tanbot == 0))
2298 t = 0;
2299 else
2300 {
2301 t = (tantop << 16) / (abs(tantop) + abs(tanbot));
2302 if (tanbot < 0)
2303 t = 0x20000 - t;
2304 else
2305 if (tantop < 0) t = 0x40000 + t;
2306 }
2307 return t;
2308 }
2309
2310
quadr(int t)2311 static int quadr(int t) {
2312 return (8 - (((t + 0x4000) >> 15) & 7)) & 7;
2313 }
2314