p3dplot.c

来自「speech signal process tools」· C语言 代码 · 共 993 行 · 第 1/2 页

/*
 * This material contains unpublished, proprietary software of 
 * Entropic Research Laboratory, Inc. Any reproduction, distribution, 
 * or publication of this work must be authorized in writing by Entropic 
 * Research Laboratory, Inc., and must bear the notice: 
 *
 *    "Copyright (c) 1986-1990  Entropic Speech, Inc. 
 *    "Copyright (c) 1990-1993  Entropic Research Laboratory, Inc. 
 *                   All rights reserved"
 *
 * The copyright notice above does not evidence any actual or intended 
 * publication of this source code.     
 *
 * Written by:  Rodney W. Johnson, Entropic Speech, Inc.
 * Checked by:
 * Revised by:
 *
 * Brief description: Plotting with hidden-line removal.
 *
 */

static char *sccs_id = "@(#)p3dplot.c	1.12	6/24/93	ESI/ERL";

#include <stdio.h>
#include <math.h>
#include <esps/esps.h>
#include "plot3d.h"


extern void	splotline();
extern void	splotpoint();
extern void	plotline();
extern void	drawline();
extern void	aplotline();
extern void	pxplotline();
extern void	make_amat();
extern void	make_pmats();
extern int	base_trans();
extern point	ptrans();
extern void	atrans();
extern void	drawing_style();
extern int	pdirx(), pdiry();
extern void	string_extents();

extern int	debug_level;
extern int      force_monochrome_plot;

static point	interp();
static point	intersect();
void		pr_points();

static void	plot_scale();
static void	end_label();
static int	residue();



void
extrema(a, m, b, n, cp, pp)
    point   *a, *b, **cp;
    int	    m, n, *pp;
{
    point   *c;
    int	    i, j, k;
    point   a0, b0, a1, b1;

    c = (point *) malloc((unsigned) 2 * (m + n) * sizeof (point));

    k = 0;

    if (m == 0)
    {
	for (j = 0; j < n; j++)
	    c[k++] = b[j];
    }
    else
    if (n == 0)
    {
	for (i = 0; i < m; i++)
	    c[k++] = a[i];
    }
    else
    {
	for (i = 0; i < m && a[i].x < b[0].x; i++)
	{
	    c[k++] = a[i];
	}

	for (j = 0; j < n && b[j].x < a[0].x; j++)
	{
	    c[k++] = b[j];
	}

	if (i == 0 && j == 0)	/* a[0].x == b[0].x */
	{
	    a0 = a[0];
	    b0 = b[0];
	    c[k++] = (a0.z > b0.z) ? a0 : b0;
	    i++;
	    j++;
	}
	else
	if (i == 0 && j < n)	/* b[j-1].x < a[0].x <= b[j].x */
	{
	    a0 = a[0];
	    if (a0.x < b[j].x)
	    {
		b0 = interp(b[j-1], b[j], a0.x);
		if (a0.z >= b0.z)
		{
		    c[k++] = a0;
		}
	    }
	    else		/* a[0].x == b[j].x */
	    {
		b0 = b[j];
		c[k++] = (a0.z > b0.z) ? a0 : b0;
		j++;
	    }
	    i++;
	}
	else	/* (j == 0) */
	if (i < m)		/* a[i-1].x < b[0].x <= a[i].x */
	{
	    b0 = b[0];
	    if (b0.x < a[i].x)
	    {
		a0 = interp(a[i-1], a[i], b0.x);
		if (b0.z >= a0.z)
		{
		    c[k++] = b0;
		}
	    }
	    else		/* b[0].x == a[i].x */
	    {
		a0 = a[i];
		c[k++] = (b0.z > a0.z) ? b0 : a0;
		i++;
	    }
	    j++;
	}

	while (i < m && j < n)
	{
	    if (a[i].x < b[j].x)
	    {
		a1 = a[i];
		b1 = interp(b0, b[j], a1.x);
		if (a1.z < b1.z)
		{
		    if (a0.z > b0.z)
		    {
			c[k++] = intersect(a0, a1, b0, b1);
		    }
		}
		else
		if (a1.z > b1.z)
		{
		    if (a0.z < b0.z)
		    {
			c[k++] = intersect(a0, a1, b0, b1);
		    }
		    c[k++] = a1;
		}
		else	/* a1.z == b1.z */
		{
		    c[k++] = a1;
		}
		i++;
	    }
	    else
	    if (a[i].x > b[j].x)
	    {
		b1 = b[j];
		a1 = interp(a0, a[i], b1.x);
		if (b1.z < a1.z)
		{
		    if (b0.z > a0.z)
		    {
			c[k++] = intersect(b0, b1, a0, a1);
		    }
		}
		else
		if (b1.z > a1.z)
		{
		    if (b0.z < a0.z)
		    {
			c[k++] = intersect(b0, b1, a0, a1);
		    }
		    c[k++] = b1;
		}
		else	/* b1.z == a1.z */
		{
		    c[k++] = b1;
		}
		j++;
	    }
	    else	/* a[i].x == b[j].x */
	    {
		a1 = a[i];
		b1 = b[j];
		if (a1.z < b1.z)
		{
		    if (a0.z > b0.z)
		    {
			c[k++] = intersect(a0, a1, b0, b1);
		    }
		    c[k++] = b1;
		}
		else
		if (b1.z < a1.z)
		{
		    if (b0.z > a0.z)
		    {
			c[k++] = intersect(b0, b1, a0, a1);
		    }
		    c[k++] = a1;
		}
		else	/* a1.z == b1.z */
		{
		    c[k++] = a1;
		}
		i++;
		j++;
	    }
	    a0 = a1;
	    b0 = b1;
	}

	while (i < m)
	{
	    c[k++] = a[i];
	    i++;
	}

	while (j < n)
	{
	    c[k++] = b[j];
	    j++;
	}
    }
	
    *cp = c;
    *pp = k;
}

static point
interp(a, b, x)
    point   a, b;
    double  x;
{
    point   c;
    double  r;

    r = (x - a.x)/(b.x - a.x);
    c.x = x;
    c.z = (1.0 - r)*a.z + r*b.z;
    return c;
}

static point
intersect(a0, a1, b0, b1)
    point   a0, a1, b0, b1;
{
    point   c;
    double  adx, adz, bdx, bdz;
    double  r;

    adx = a1.x - a0.x;	adz = a1.z - a0.z;
    bdx = b1.x - b0.x;	bdz = b1.z - b0.z;
    r = ((b0.x - a0.x)*bdz - bdx*(b0.z - a0.z)) / (adx*bdz - bdx*adz);
    c.x = a0.x + r*(a1.x - a0.x);
    c.z = a0.z + r*(a1.z - a0.z);
    return c;
}


void
plot(data, x, y, nrecs, nitems)
    float   **data;
    double  *x, *y;
    int	    nrecs, nitems;
{
    int	    i, i0, i1, j, k;
    point   *line, *hiline, *new_hiline;
    int	    h0, h1, hilen;
    int	    dirymin, dirymax;

    drawing_style(BFN_SRC, GET_COL_TOP_FG, LNS_SOLID);

    make_pmats(1, 1);

    dirymin = pdiry(y[0]);
    dirymax = pdiry(y[nitems-1]);

    if (debug_level)
	(void) fprintf(stderr,
		       "plot: dirymin = %d, dirymax = %d\n", dirymin, dirymax);

    line = (point *) malloc((unsigned) (nitems + 2) * sizeof(point));
    for (i = 0; i < nrecs && pdirx(x[i]) < 0; i++)
    { }
    i0 = i-1;

    if (debug_level)
	(void) fprintf(stderr, "plot: i0 = %d\n", i0);

    for ( ; i < nrecs && pdirx(x[i]) == 0; i++)
    {
	for (j = 0; j+1 < nitems; j++)
	    splotline(x[i], y[j], data[i][j], x[i], y[j+1], data[i][j+1]);
    }
    i1 = i;

    if (debug_level)
	(void) fprintf(stderr, "plot: i1 = %d\n", i1);

/*!*//* May need edge segment connecting to previous slice. */
    hiline = (point *) malloc((unsigned) 1 * sizeof(point));
    hilen = 0;

    for ( ; i < nrecs; i++)
    {
	k = 0;

	h0 = 0;
	switch (dirymin)
	{
	case -1:
	    if (i > i1)
	    {
		line[k++] = ptrans(x[i-1], y[0], data[i-1][0]);
		for (h0 = 1; h0 < hilen && hiline[h0].x <= line[0].x; h0++)
		{ }
		h0--;
	    }
	    break;
	case 0:
	    if (i > i1)
	    {
		splotline(x[i-1], y[0], data[i-1][0], x[i], y[0], data[i][0]);
		for (h0 = 1; h0 < hilen && hiline[h0].x <= line[0].x; h0++)
		{ }
		h0--;
	    }
	    break;
	case 1:
	    if (i+1 < nrecs)
		line[k++] = ptrans(x[i+1], y[0], data[i+1][0]);
	    break;
	}

	for (j = 0; j < nitems; j++)
	    line[k++] =
		ptrans(x[i], y[j], data[i][j]);

	h1 = hilen;
	switch (dirymax)
	{
	case -1:
	    if (i+1 < nrecs)
		line[k++] = ptrans(x[i+1], y[nitems-1], data[i+1][nitems-1]);
	    break;
	case 0:
	    if (i > i1)
	    {
		splotline(x[i-1], y[nitems-1], data[i-1][nitems-1],
			  x[i], y[nitems-1], data[i][nitems-1]);
		for (h1 = hilen-2; h1 >= 0 && hiline[h1].x >= line[k-1].x; h1--)
		{ }
		h1 += 2;
	    }
	    break;
	case 1:
	    if (i > i1)
	    {
		line[k++] = ptrans(x[i-1], y[nitems-1], data[i-1][nitems-1]);
		for (h1 = hilen-2; h1 >= 0 && hiline[h1].x >= line[k-1].x; h1--)
		{ }
		h1 += 2;
	    }
	    break;
	}

        if (debug_level >= 3)
	{
	    fprintf(stderr, "plot: i = %d\n", i);
	    pr_points("line", k, line);
	    pr_points("hiline", h1 - h0, hiline + h0);
	}

	extrema(line, k, hiline + h0, h1 - h0, &new_hiline, &hilen);
	free((char *) hiline);
	hiline = new_hiline;

	if (debug_level >= 3)
	{
	    pr_points("new hiline", hilen, hiline);
	    for (j = 0; j+1< hilen && hiline[j+1].x <= line[k-1].x; j++)
	    { }
	    (void) fprintf(stderr, "plotting %d segments\n", j);
	}

	for (j = 0; j+1< hilen && hiline[j+1].x <= line[k-1].x; j++)
	    pxplotline(hiline[j].x, hiline[j].z, hiline[j+1].x, hiline[j+1].z);
    }

    make_pmats(-1, 1);

/*!*//* May need edge segment connecting to previous slice. */
    hiline = (point *) malloc((unsigned) 1 * sizeof(point));
    hilen = 0;

    for (i = i0; i >= 0; i--)
    {
	k = 0;

	h0 = 0;
	switch (dirymin)
	{
	case -1:
	    if (i < i0)
	    {
		line[k++] = ptrans(x[i+1], y[0], data[i+1][0]);
		for (h0 = 1; h0 < hilen && hiline[h0].x <= line[0].x; h0++)
		{ }
		h0--;
	    }
	    break;
	case 0:
	    if (i < i0)
	    {
		splotline(x[i+1], y[0], data[i+1][0], x[i], y[0], data[i][0]);
		for (h0 = 1; h0 < hilen && hiline[h0].x <= line[0].x; h0++)
		{ }
		h0--;
	    }
	    break;
	case 1:
	    if(i-1 >= 0)
		line[k++] = ptrans(x[i-1], y[0], data[i-1][0]);
	    break;
	}

	for (j = 0; j < nitems; j++)
	    line[k++] =
		ptrans(x[i], y[j], data[i][j]);

	h1 = hilen;
	switch (dirymax)
	{
	case -1:
	    if (i-1 >= 0)
		line[k++] = ptrans(x[i-1], y[nitems-1], data[i-1][nitems-1]);
	    break;
	case 0:
	    if (i < i0)
	    {
		splotline(x[i+1], y[nitems-1], data[i+1][nitems-1],
			  x[i], y[nitems-1], data[i][nitems-1]);
		for (h1 = hilen-2; h1 >= 0 && hiline[h1].x >= line[k-1].x; h1--)
		{ }
		h1 += 2;
	    }
	    break;
	case 1:
	    if (i < i0)
	    {
		line[k++] = ptrans(x[i+1], y[nitems-1], data[i+1][nitems-1]);
		for (h1 = hilen-2; h1 >= 0 && hiline[h1].x >= line[k-1].x; h1--)
		{ }
		h1 += 2;
	    }
	    break;
	}

        if (debug_level >= 3)
	{
	    fprintf(stderr, "plot: i = %d\n", i);
	    pr_points("line", k, line);
	    pr_points("hiline", h1 - h0, hiline + h0);
	}

	extrema(line, k, hiline + h0, h1 - h0, &new_hiline, &hilen);
	free((char *) hiline);
	hiline = new_hiline;

	if (debug_level >= 3)
	{
	    pr_points("new hiline", hilen, hiline);
	    for (j = 0; j+1< hilen; j++)
	    { }
	    (void) fprintf(stderr, "plotting %d segments\n", j);
	}

	for (j = 0; j+1< hilen; j++)
	    pxplotline(hiline[j].x, hiline[j].z, hiline[j+1].x, hiline[j+1].z);