li_recognizer.c

这是一个同样来自贝尔实验室的和UNIX有着渊源的操作系统, 其简洁的设计和实现易于我们学习和理解

		return(dist);
}


/*************************************************************

  Digest-processing routines

 *************************************************************/

static int lialg_read_classifier_digest(rClassifier *rec) {
	int nclasses;
	FILE *fp;

	/* Try to open the corresponding digest file. */
	{
				char *clx_path;
				char *dot;
		
				/* Get a copy of the filename, with some room on the end. */
				/*	clx_path = safe_malloc(strlen(rec->file_name) + 5); */
				clx_path = malloc((strlen(rec->file_name) + 5) *sizeof(char));
				strcpy(clx_path, rec->file_name);
		
				/* Truncate the path after the last dot. */
				dot = strrchr(clx_path, '.');
				if (dot == nil) { free(clx_path); return(-1); }
				*(dot + 1) = 0;
		
				/* Append the classifier-digest extension. */
				strcat(clx_path, "clx");
		
				fp = fopen(clx_path, "r");
				if (fp == nil) {
						free(clx_path);
						return(-1);
				}
		
				free(clx_path);
	}

	/* Read-in the name and dominant points for each class. */
	for (nclasses = 0; !feof(fp); nclasses++) {
		point_list *dpts = nil;
		char class[BUFSIZ];
		int npts;
		int j;

		if (fscanf(fp, "%s %d", class, &npts) != 2) {
			if (feof(fp)) break;

			goto failed;
		}
		rec->cnames[nclasses] = strdup(class);

		/* Allocate a dominant-points list. */
		/* dpts = (point_list *)safe_malloc(sizeof(point_list)); */
		dpts = malloc(sizeof(point_list));
		dpts->pts = mallocz(npts*sizeof(pen_point), 1);
		if (dpts->pts == nil) goto failed;
		dpts->npts = npts;
		dpts->next = nil;

		/* Read in each point. */
		for (j = 0; j < npts; j++) {
			int x, y;

			if (fscanf(fp, "%d %d", &x, &y) != 2) goto failed;
			dpts->pts[j].x = x;
			dpts->pts[j].y = y;
		}

		/* Compute the chain-code. */
		lialg_compute_chain_code(dpts);

		/* Store the list in the rec data structure. */
		rec->dompts[nclasses] = dpts;

		continue;

failed:
		fprint(2, "read_classifier_digest failed...\n");
		for (; nclasses >= 0; nclasses--) {
			if (rec->cnames[nclasses] != nil) {
				free(rec->cnames[nclasses]);
				rec->cnames[nclasses] = nil;
			}
			if (rec->dompts[nclasses] != nil) {
				delete_examples(rec->dompts[nclasses]);
				rec->dompts[nclasses] = nil;
			}
		}
		if (dpts != nil)
			delete_examples(dpts);
		fclose(fp);
		return(-1);
	}

	fclose(fp);
	return(0);
}


/*************************************************************

  Canonicalization routines

 *************************************************************/

static int lialg_canonicalize_examples(rClassifier *rec) {
	int i;
	int nclasses;

	if (lidebug) {
		fprint(2, "lialg_canonicalize_examples working on %s\n",
				rec->file_name);
	}
	/* Initialize canonical-example arrays. */
	for (i = 0; i < MAXSCLASSES; i++) {
		rec->canonex[i] = nil;
	}

	/* Figure out number of classes. */
	for (nclasses = 0;
		  nclasses < MAXSCLASSES && rec->cnames[nclasses] != nil;
		  nclasses++)
		;

	/* Canonicalize the examples for each class. */
	for (i = 0; i < nclasses; i++) {
		int j, k;
		int nex;
		point_list *pts, *tmp, *avg;
		int maxxrange, maxyrange;
		int minx, miny, maxx, maxy;
		int avgxrange, avgyrange, avgxoff, avgyoff, avgscale;

		
		if (lidebug) {
			fprint(2, "lialg_canonicalize_examples working on class %s\n",
					rec->cnames[i]);
		}
		/* Make a copy of the examples. */
		pts = nil;
		tmp = rec->ex[i];
		for (nex = 0; tmp != nil; nex++, tmp = tmp->next) {
			if ((pts = add_example(pts, tmp->npts, tmp->pts)) == nil) {
				delete_examples(pts);
				return(-1);
			}
		}

		/* Canonicalize each example, and derive the max x and y ranges. */
		maxxrange = 0;
		maxyrange = 0;
		for (j = 0, tmp = pts; j < nex; j++, tmp = tmp->next) {
			if (lialg_canonicalize_example_stroke(tmp) != 0) {
  	        if (lidebug) {
					fprint(2, "lialg_canonicalize_example_stroke returned error\n");
				}
				return(-1);
			}

			if (tmp->xrange > maxxrange) maxxrange = tmp->xrange;
			if (tmp->yrange > maxyrange) maxyrange = tmp->yrange;
		}

		/* Normalize max ranges. */
		if (((100 * maxxrange + CANONICAL_X / 2) / CANONICAL_X) >
			((100 * maxyrange + CANONICAL_Y / 2) / CANONICAL_Y)) {
			maxyrange = (maxyrange * CANONICAL_X + maxxrange / 2) / maxxrange;
			maxxrange = CANONICAL_X;
		}
		else {
			maxxrange = (maxxrange * CANONICAL_Y + maxyrange / 2) / maxyrange;
			maxyrange = CANONICAL_Y;
		}

		/* Re-scale each example to max ranges. */
		for (j = 0, tmp = pts; j < nex; j++, tmp = tmp->next) {
			int scalex = (tmp->xrange == 0) ? 100 : (100 * maxxrange + tmp->xrange / 2) / tmp->xrange;
			int scaley = (tmp->yrange == 0) ? 100 : (100 * maxyrange + tmp->yrange / 2) / tmp->yrange;
			if (lialg_translate_points(tmp, 0, 0, scalex, scaley) != 0) {
				delete_examples(pts);
				return(-1);
			}
		}

		/* Average the examples; leave average in first example. */
		avg = pts;				/* careful aliasing!! */
		for (k = 0; k < NCANONICAL; k++) {
			int xsum = 0;
			int ysum = 0;

			for (j = 0, tmp = pts; j < nex; j++, tmp = tmp->next) {
						xsum += tmp->pts[k].x;
						ysum += tmp->pts[k].y;
			}

			avg->pts[k].x = (xsum + j / 2) / j;
			avg->pts[k].y = (ysum + j / 2) / j;
		}

		/* Compute BB of averaged stroke and re-scale. */
		lialg_get_bounding_box(avg, &minx, &miny, &maxx, &maxy);
		avgxrange = maxx - minx;
		avgyrange = maxy - miny;
		avgscale = (((100 * avgxrange + CANONICAL_X / 2) / CANONICAL_X) >
					((100 * avgyrange + CANONICAL_Y / 2) / CANONICAL_Y))
		  ? (100 * CANONICAL_X + avgxrange / 2) / avgxrange
		  : (100 * CANONICAL_Y + avgyrange / 2) / avgyrange;
		if (lialg_translate_points(avg, minx, miny, avgscale, avgscale) != 0) {
			delete_examples(pts);
			return(-1);
		}

		/* Re-compute the x and y ranges and center the stroke. */
		lialg_get_bounding_box(avg, &minx, &miny, &maxx, &maxy);
		avgxrange = maxx - minx;
		avgyrange = maxy - miny;
		avgxoff = -((CANONICAL_X - avgxrange + 1) / 2);
		avgyoff = -((CANONICAL_Y - avgyrange + 1) / 2);
		if (lialg_translate_points(avg, avgxoff, avgyoff, 100, 100) != 0) {
			delete_examples(pts);
			return(-1);
		}

		/* Create a point list to serve as the ``canonical representation. */
		if ((rec->canonex[i] = add_example(nil, avg->npts, avg->pts)) == nil) {
			delete_examples(pts);
			return(-1);
		}
		(rec->canonex[i])->xrange = maxx - minx;
		(rec->canonex[i])->yrange = maxy - miny;

		if (lidebug) {
			fprint(2, "%s, avgpts = %d\n", rec->cnames[i], avg->npts);
			for (j = 0; j < avg->npts; j++) {
						fprint(2, "  (%P)\n", avg->pts[j].Point);
			}
		}

		/* Compute dominant points of canonical representation. */
		rec->dompts[i] = lialg_compute_dominant_points(avg);

		/* Clean up. */
		delete_examples(pts);
	}

	/* Sanity check. */
	for (i = 0; i < nclasses; i++) {
		char *best_name = lialg_recognize_stroke(rec, rec->canonex[i]);

		if (best_name != rec->cnames[i])
			fprint(2, "%s, best = %s\n", rec->cnames[i], best_name);
	}

	return(0);
}


static int lialg_canonicalize_example_stroke(point_list *points) {
	int minx, miny, maxx, maxy, xrange, yrange, scale;

	/* Filter out points that are too close. */
	if (lialg_filter_points(points) != 0) return(-1);

	/* Must be at least two points! */
	if (points->npts < 2) {
		if (lidebug) {
			fprint(2, "lialg_canonicalize_example_stroke: npts=%d\n",
					points->npts);
		}
		return(-1);
	}

	/* Scale up to avoid conversion errors. */
	lialg_get_bounding_box(points, &minx, &miny, &maxx, &maxy);
	xrange = maxx - minx;
	yrange = maxy - miny;
	scale = (((100 * xrange + CANONICAL_X / 2) / CANONICAL_X) >
			 ((100 * yrange + CANONICAL_Y / 2) / CANONICAL_Y))
	  ? (100 * CANONICAL_X + xrange / 2) / xrange
	  : (100 * CANONICAL_Y + yrange / 2) / yrange;
	if (lialg_translate_points(points, minx, miny, scale, scale) != 0) {
		if (lidebug) {
			fprint(2, "lialg_translate_points (minx=%d,miny=%d,scale=%d) returned error\n", minx, miny, scale);
		}
		return(-1);
	}

	/* Compute an equivalent stroke with equi-distant points. */
	if (lialg_compute_equipoints(points) != 0) return(-1);

	/* Re-translate the points to the origin. */
	lialg_get_bounding_box(points, &minx, &miny, &maxx, &maxy);
	if (lialg_translate_points(points, minx, miny, 100, 100) != 0) {
		if (lidebug) {
			fprint(2, "lialg_translate_points (minx=%d,miny=%d) returned error\n", minx, miny);
		}
		return(-1);
	}

	/* Store the x and y ranges in the point list. */
	xrange = maxx - minx;
	yrange = maxy - miny;
	points->xrange = xrange;
	points->yrange = yrange;

	if (lidebug) {
		int i;
		fprint(2, "Canonicalized:   %d, %d, %d, %d\n", minx, miny, maxx, maxy);
		for (i = 0; i < points->npts; i++)
			fprint(2, "      (%P)\n", points->pts[i].Point);
		fflush(stderr);
	}

	return(0);
}


static int lialg_compute_equipoints(point_list *points) {
	pen_point *equipoints = mallocz(NCANONICAL*sizeof(pen_point), 1);
	int nequipoints = 0;
	int pathlen = lialg_compute_pathlen(points);
	int equidist = (pathlen + (NCANONICAL - 1) / 2) / (NCANONICAL - 1);
	int i;
	int dist_since_last_eqpt;
	int remaining_seglen;
	int dist_to_next_eqpt;

	if (equipoints == nil) {
		fprint(2, "can't allocate memory in lialg_compute_equipoints");
		return(-1);
	}

	if (lidebug) {
		fprint(2, "compute_equipoints:  npts = %d, pathlen = %d, equidist = %d\n",
				points->npts, pathlen, equidist);
		fflush(stderr);
	}

	/* First original point is an equipoint. */
	equipoints[0] = points->pts[0];
	nequipoints++;
	dist_since_last_eqpt = 0;

	for (i = 1; i < points->npts; i++) {
		int dx1 = points->pts[i].x - points->pts[i-1].x;
		int dy1 = points->pts[i].y - points->pts[i-1].y;
		int endx = 100 * points->pts[i-1].x;
		int endy = 100 * points->pts[i-1].y;
		remaining_seglen = isqrt(10000 * (dx1 * dx1 + dy1 * dy1));
		dist_to_next_eqpt = equidist - dist_since_last_eqpt;

		while (remaining_seglen >= dist_to_next_eqpt) {
			if (dx1 == 0) {
				/* x-coordinate stays the same */
				if (dy1 >= 0)
					endy += dist_to_next_eqpt;
				else
					endy -= dist_to_next_eqpt;
			}
			else {
				int slope = (100 * dy1 + dx1 / 2) / dx1;
				int tmp = isqrt(10000 + slope * slope);
				int dx = (100 * dist_to_next_eqpt + tmp / 2) / tmp;
				int dy = (slope * dx + 50) / 100;

				if (dy < 0) dy = -dy;
				if (dx1 >= 0)
					endx += dx;
				else
					endx -= dx;
				if (dy1 >= 0)
					endy += dy;
				else
					endy -= dy;
			}

			equipoints[nequipoints].x = (endx + 50) / 100;
			equipoints[nequipoints].y = (endy + 50) / 100;
			nequipoints++;
/*	    assert(nequipoints <= NCANONICAL);*/
			dist_since_last_eqpt = 0;
			remaining_seglen -= dist_to_next_eqpt;
			dist_to_next_eqpt = equidist;
		}

		dist_since_last_eqpt += remaining_seglen;
	}

	/* Take care of last equipoint. */
	if (nequipoints == NCANONICAL) {
		/* Good. */
	} else if (nequipoints == (NCANONICAL - 1)) {
		/* Make last original point the last equipoint. */
		equipoints[nequipoints] = points->pts[points->npts - 1];
	} else {
	  if (lidebug) {
		fprint(2,"lialg_compute_equipoints: nequipoints = %d\n", 
				nequipoints);
	  }
/*	assert(false);*/
		return(-1);
	}

	points->npts = NCANONICAL;
	delete_pen_point_array(points->pts);
	points->pts = equipoints;
	return(0);
}


/*************************************************************

  Utility routines

 *************************************************************/

/* Result is x 100. */
static int lialg_compute_pathlen(point_list *points) {
	return(lialg_compute_pathlen_subset(points, 0, points->npts - 1));
}


/* Result is x 100. */
static int lialg_compute_pathlen_subset(point_list *points,
										   int start, int end) {
	int pathlen;
	int i;

	pathlen = 0;
	for (i = start + 1; i <= end; i++) {
		int dx = points->pts[i].x - points->pts[i-1].x;
		int dy = points->pts[i].y - points->pts[i-1].y;
		int dist = isqrt(10000 * (dx * dx + dy * dy));
		pathlen += dist;
	}

	return(pathlen);
}


/* Note that this does NOT update points->xrange and points->yrange! */
static int lialg_filter_points(point_list *points) {
	int filtered_npts;
	pen_point *filtered_pts = mallocz(points->npts*sizeof(pen_point), 1);
	int i;

	if (filtered_pts == nil) {
		fprint(2, "can't allocate memory in lialg_filter_points");
		return(-1);
	}

	filtered_pts[0] = points->pts[0];
	filtered_npts = 1;
	for (i = 1; i < points->npts; i++) {
		int j = filtered_npts - 1;
		int dx = points->pts[i].x - filtered_pts[j].x;
		int dy = points->pts[i].y - filtered_pts[j].y;
		int magsq = dx * dx + dy * dy;

		if (magsq >= DIST_SQ_THRESHOLD) {
			filtered_pts[filtered_npts] = points->pts[i];
			filtered_npts++;
		}
	}

	points->npts = filtered_npts;
	delete_pen_point_array(points->pts);
	points->pts = filtered_pts;
	return(0);
}


/* scalex and scaley are x 100. */
/* Note that this does NOT update points->xrange and points->yrange! */
static int lialg_translate_points(point_list *points,
								   int minx, int miny,
								   int scalex, int scaley) {
	int i;

	for (i = 0; i < points->npts; i++) {
				points->pts[i].x = ((points->pts[i].x - minx) * scalex + 50) / 100;
				points->pts[i].y = ((points->pts[i].y - miny) * scaley + 50) / 100;
	}

	return(0);
}


static void lialg_get_bounding_box(point_list *points,
									int *pminx, int *pminy,
									int *pmaxx, int *pmaxy) {
	int minx, miny, maxx, maxy;
	int i;

	minx = maxx = points->pts[0].x;
	miny = maxy = points->pts[0].y;
	for (i = 1; i < points->npts; i++) {
				pen_point *pt = &(points->pts[i]);
				if (pt->x < minx) minx = pt->x;
				else if (pt->x > maxx) maxx = pt->x;
				if (pt->y < miny) miny = pt->y;
				else if (pt->y > maxy) maxy = pt->y;
	}

	*pminx = minx;
	*pminy = miny;
	*pmaxx = maxx;
	*pmaxy = maxy;
}


int wtvals[] = {100, 104, 117, 143, 189, 271, 422};

static void lialg_compute_lpf_parameters(void) {
	int i;

		for (i = LP_FILTER_WIDTH; i >= 0; i--) {
//		double x = 0.04 * (i * i);
//		double tmp = 100.0 * exp(x);
//		int wt = floor((double)tmp);
				int wt = wtvals[i];
				lialg_lpfwts[LP_FILTER_WIDTH - i] = wt;
				lialg_lpfwts[LP_FILTER_WIDTH + i] = wt;
		}
		lialg_lpfconst = 0;
		for (i = 0; i < (2 * LP_FILTER_WIDTH + 1); i++) {
				lialg_lpfconst += lialg_lpfwts[i];
		}
}


/* Code from Joseph Hall ([e-mail]). */
static int isqrt(int n) {
		register int i;
		register long k0, k1, nn;

		for (nn = i = n, k0 = 2; i > 0; i >>= 2, k0 <<= 1)
				;
		nn <<= 2;
	for (;;) {
				k1 = (nn / k0 + k0) >> 1;
				if (((k0 ^ k1) & ~1) == 0)
				break;
				k0 = k1;
		}
		return (int) ((k1 + 1) >> 1);
}


/* Helper routines from Mark Hayter. */
static int likeatan(int tantop, int tanbot) { 
		int t;
		/* Use tan(theta)=top/bot --> order for t */
		/* t in range 0..0x40000 */

		if ((tantop == 0) && (tanbot == 0)) 
				t = 0;
	else
		{
				t = (tantop << 16) / (abs(tantop) + abs(tanbot));
				if (tanbot < 0) 
						t = 0x20000 - t;
				else 
						if (tantop < 0) t = 0x40000 + t;
		}
		return t;
}


static int quadr(int t) {
	return (8 - (((t + 0x4000) >> 15) & 7)) & 7;
}