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* 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 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* 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 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 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 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 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 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 294 static bool isa_li(recognizer r) 295 { return(CHECK_LI_MAGIC(r)); } 296 297 static int 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 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 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 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 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 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 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 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 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* 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 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 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* 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 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 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 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* 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** 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 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 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 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 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 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 */ 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 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 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! */ 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. */ 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 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 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 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 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 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 */ 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 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 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 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 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 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 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. */ 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. */ 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! */ 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! */ 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 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 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). */ 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. */ 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 2311 static int quadr(int t) { 2312 return (8 - (((t + 0x4000) >> 15) & 7)) & 7; 2313 } 2314