/* Dissertation.c
 *
 * Copyright (C) 1992-2002 Paul Boersma
 *
 * This program is free software; you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation; either version 2 of the License, or (at
 * your option) any later version.
 *
 * This program is distributed in the hope that it will be useful, but
 * WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
 * See the GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this program; if not, write to the Free Software
 * Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
 */

/*
 * pb 1996/05/18
 * pb 2002/07/16 GPL
 */

#include "Art_Speaker.h"
#include "Sound.h"
#include "Polygon.h"

/* The procedures are free to change all graphics attributes, */
/* including line type and font, */
/* but a changed viewport should be set back to where it was. */
/* You can use the following graphics attributes; */
/* they represent the values upon entry of the procedures, */
/* i.e., they are implicit arguments to the procedures. */
static int oldLineType, oldFont, oldSize, oldStyle;
static double oldLineWidth;

#define WALL_GREY  0.8
#define WALL_HWIDTH  1

/********** FIGURE 4.1: MANY TUBES **********/

static void spring (Any graphics, double x, double bottom, double top) {
	const double smallRadius = 0.01, largeRadius = 0.03, dy = largeRadius - smallRadius;
	const double springLength = 6 * dy;
	double totalLength = top - bottom, lineLength = 0.5 * (totalLength - springLength);
	double springTop = top - lineLength, springBottom = bottom + lineLength, y = springTop;
	double xSmall = x + smallRadius, xLarge = xSmall;
	Graphics_line (graphics, x, top, x, y);
	Graphics_arc (graphics, xSmall, y, smallRadius, 180, 90);
	Graphics_arc (graphics, xLarge, y -= dy, largeRadius, 90, 270);
	Graphics_arc (graphics, xSmall, y -= dy, smallRadius, 270, 90);
	Graphics_arc (graphics, xLarge, y -= dy, largeRadius, 90, 270);
	Graphics_arc (graphics, xSmall, y -= dy, smallRadius, 270, 90);
	Graphics_arc (graphics, xLarge, y -= dy, largeRadius, 90, 270);
	Graphics_arc (graphics, xSmall, y -= dy, smallRadius, 270, 180);
	Graphics_line (graphics, x, bottom, x, y);
}
static void manyTubes (Graphics graphics)   {
	const double wallThickness = 0.02, yhalf = 0.35;
	static double x [] = { 0.0,   0.1,  0.2,   0.3,  0.4,  0.43, 0.51,   0.67,   0.77,  0.92,  1.0 };
	static double w [] = {    0.36,  0.3,  0.24, 0.05,  0.01,  0.14,   0.2,   0.12,   0.08,  0.16 };
	static float xcage [] = { 0.0, 0.0, 1.0, 1.0, -0.02, -0.02, 1.0, 1.0, 0.0 };
	static float ycage [] = { 0.0, 0.7, 0.7, 0.72, 0.72, -0.02, -0.02, 0.0, 0.0 };
	int itube = 0;
	Graphics_setWindow (graphics, 0, 1, 0, 1);
	Graphics_setLineWidth (graphics, 2.0);
	while (x [itube] < 1.0) {
		spring (graphics, 0.5 * (x [itube] + x [itube + 1]), 0, yhalf - 0.5 * w [itube] - wallThickness);
		spring (graphics, 0.5 * (x [itube] + x [itube + 1]), yhalf + 0.5 * w [itube] + wallThickness, 0.7);
		Graphics_setGrey (graphics, WALL_GREY);
		Graphics_fillRectangle (graphics, x [itube], x [itube + 1],
			yhalf - 0.5 * w [itube] - wallThickness, yhalf - 0.5 * w [itube]);
		Graphics_fillRectangle (graphics, x [itube], x [itube + 1],
			yhalf + 0.5 * w [itube], yhalf + 0.5 * w [itube] + wallThickness);
		Graphics_setGrey (graphics, 0.0);
		Graphics_rectangle (graphics, x [itube], x [itube + 1],
			yhalf - 0.5 * w [itube] - wallThickness, yhalf - 0.5 * w [itube]);
		Graphics_rectangle (graphics, x [itube], x [itube + 1],
			yhalf + 0.5 * w [itube], yhalf + 0.5 * w [itube] + wallThickness);
		if (itube > 0) {
			Graphics_line (graphics, x [itube], yhalf - 0.5 * w [itube - 1], x [itube], yhalf - 0.5 * w [itube]);
			Graphics_line (graphics, x [itube], yhalf + 0.5 * w [itube - 1], x [itube], yhalf + 0.5 * w [itube]);
		}
		itube ++;
	}
	Graphics_setGrey (graphics, WALL_GREY);
	Graphics_fillArea (graphics, 9, xcage, ycage);
	Graphics_setGrey (graphics, 0.0);
	Graphics_polyline (graphics, 9, xcage, ycage);
}

/********** FIGURE 4.4: TYPES OF TUBE BOUNDARIES **********/

static void rotate (float *x, float *y, double angle) {
	double help;
	if (angle == 0.0) return;
	angle *= NUMpi / 180;
	help = *x * cos (angle) - *y * sin (angle);
	*y = *x * sin (angle) + *y * cos (angle);
	*x = help;
} 
static void draw1Wall (Graphics graphics, double xc, double yc, double hlength, double angle) {
	double hwidth = WALL_HWIDTH, springLength = 5;
	float x [5], y [5];
	int i;
	Graphics_setGrey (graphics, WALL_GREY);
	x [4] = x [0] = x [3] = - hlength; x [1] = x [2] = hlength;
	y [4] = y [0] = y [1] = - hwidth; y [2] = y [3] = hwidth;
	for (i = 0; i <= 4; i ++) { rotate (& x [i], & y [i], angle); x [i] += xc; y [i] += yc; }
	Graphics_fillArea (graphics, 4, x, y);
	Graphics_setGrey (graphics, 0);
	Graphics_polyline (graphics, 5, x, y);
	x [4] = x [0] = x [3] = -0.3; x [1] = x [2] = 0.3;
	y [4] = y [0] = y [1] = hwidth; y [2] = y [3] = hwidth + springLength;
	for (i = 0; i <= 4; i ++) { rotate (& x [i], & y [i], angle); x [i] += xc; y [i] += yc; }
	Graphics_fillArea (graphics, 4, x, y);
}
static void draw2Walls (Graphics graphics, double xc, double yc, double hwidth, double hlength, double angle) {
	float x, y;
	x = 0; y = hwidth;
	rotate (& x, & y, angle);
	draw1Wall (graphics, xc + x, yc + y, hlength, angle);
	draw1Wall (graphics, xc - x, yc - y, hlength, angle + 180);
}
static void tubeConnections (Graphics graphics) {
	float centre = 83.333333;
	Graphics_setWindow (graphics, 0, 100, 0, 100);
	/* Tube closed at one end. */
	{
		float x = 7, d = 5, w = 6;
		draw2Walls (graphics, x, centre, w, d, 0);   /* The tube. */
		Graphics_line (graphics, x-d, centre - w, x-d, centre + w);   /* Back vertical wall. */
	}
	/* Two connected tubes. */
	{
		float x = 24, d1 = 5, d2 = 4, w1 = 3, w2 = 6;
		draw2Walls (graphics, x, centre, w1, d1, 0);   /* Left tube. */
		draw2Walls (graphics, x + d1 + d2, centre, w2, d2, 0);   /* Right tube. */
		x += d1;
		Graphics_line (graphics, x, centre - w1, x, centre - w2);   /* Lower connection. */
		Graphics_line (graphics, x, centre + w1, x, centre + w2);   /* Upper connection. */
	}
	/* Three connected tubes. */
	{
		float x1 = 48, x2 = x1 + 15, x3 = x1 + 15, y2 = centre + 10, y3 = centre - 10;
		float d = 4, w1 = 3.5, w2 = 2, w3 = 5;
		float xb, yb, xe, ye;
		draw2Walls (graphics, x1, centre, w1, d, 0);   /* Left tube. */
		draw2Walls (graphics, x2, y2, w2, d, 60);   /* Top right tube. */
		draw2Walls (graphics, x3, y3, w3, d, -60);   /* Bottom right tube. */
		xb = x1 + d; yb = centre + w1 - WALL_HWIDTH;
		xe = - d; ye = w2 - WALL_HWIDTH; rotate (& xe, & ye, 60);
		Graphics_line (graphics, xb, yb, x2 + xe, y2 + ye);   /* Connection from left to top. */
		xb = x1 + d; yb = centre - w1 + WALL_HWIDTH;
		xe = - d; ye = - w3 + WALL_HWIDTH; rotate (& xe, & ye, -60);
		Graphics_line (graphics, xb, yb, x3 + xe, y3 + ye);   /* Connection from left to bottom. */
		xb = - d; yb = - w2 + WALL_HWIDTH; rotate (& xb, & yb, 60);
		xe = - d; ye = w3 - WALL_HWIDTH; rotate (& xe, & ye, -60);
		Graphics_line (graphics, x2 + xb, y2 + yb, x3 + xe, y3 + ye);   /* From top to bottom. */
	}
	/* Open tube with radiation. */
	{
		float x = 85, d = 5, w = 6;
		draw2Walls (graphics, x, centre, w, d, 0);   /* The tube. */
		x += d - 1;   /* Radiation. */
		Graphics_line (graphics, x, centre, x + 7.071, centre);   /* Straight ahead. */
		Graphics_arc (graphics, x, centre - 12, 10, 45, 90);   /* Bending down. */
		Graphics_arc (graphics, x, centre + 12, 10, 270, 315);   /* Bending up. */
	}
}

/********** FIGURE 4.7 **********/

static void vocalTractConstruction (Graphics graphics) {
	Art art = Art_create ();
	Speaker speaker = Speaker_create ("Female", 2);
	double f = speaker -> relativeSize * 1e-3;
	double intX [1 + 16], intY [1 + 16], extX [1 + 11], extY [1 + 11];
	double bodyX, bodyY;
	Graphics_Viewport oldViewport;

	Graphics_setLineType (graphics, Graphics_DOTTED);
	Graphics_setLineWidth (graphics, 2.0);
	Art_Speaker_toVocalTract (art, speaker, intX, intY, extX, extY, & bodyX, & bodyY);
	oldViewport = Graphics_insetViewport (graphics, 0.1, 0.9, 0.1, 0.9);
	Graphics_setWindow (graphics, -0.05, 0.05, -0.05, 0.05);
	/* Lines from joint. */
	{
		double jaw_x = -75 * f, jaw_y = 53 * f;
		double teeth_x = intX [11], teeth_y = intY [11];   /* Lower teeth. */
		double jaw_a = atan2 (teeth_y - jaw_y, teeth_x - jaw_x);
		double body_a = atan2 (bodyY - jaw_y, bodyX - jaw_x);
		Graphics_setTextAlignment (graphics, Graphics_RIGHT, Graphics_HALF);
		Graphics_fillCircle_mm (graphics, jaw_x, jaw_y, 1.0);
		Graphics_text (graphics, jaw_x - f, jaw_y, "#J");
		Graphics_setTextAlignment (graphics, Graphics_CENTRE, Graphics_BOTTOM);
		Graphics_text (graphics, jaw_x, jaw_y + 2*f, "%%joint");
		Graphics_line (graphics, jaw_x, jaw_y, 0.05, jaw_y);
		Graphics_line (graphics, jaw_x, jaw_y, teeth_x, teeth_y);
		Graphics_arcArrow (graphics, jaw_x, jaw_y, 0.04, jaw_a * 180 / NUMpi, 0, 1, 0);
		Graphics_setTextAlignment (graphics, Graphics_LEFT, Graphics_BOTTOM);
		Graphics_text (graphics, -74*f + 0.04 * cos (jaw_a), 53*f + 0.04 * sin (jaw_a), "-0.3");
		Graphics_line (graphics, jaw_x, jaw_y, bodyX, bodyY);
		Graphics_arcArrow (graphics, jaw_x, jaw_y, 0.03, body_a * 180 / NUMpi, 0, 1, 0);
		Graphics_setTextAlignment (graphics, Graphics_LEFT, Graphics_BOTTOM);
		Graphics_text (graphics, -73*f + 0.03 * cos (body_a), 53*f + 0.03 * sin (body_a), "-0.6");
	}
	Graphics_circle (graphics, bodyX, bodyY, 20 * f);
	Graphics_line (graphics, bodyX, bodyY, intX [6], intY [6]);
	Graphics_line (graphics, bodyX, bodyY, intX [7], intY [7]);
	Graphics_setTextAlignment (graphics, Graphics_LEFT, Graphics_HALF);
	Graphics_fillCircle_mm (graphics, intX [4], intY [4], 1.0);
	Graphics_text (graphics, intX [4], intY [4], "#H");
	Graphics_fillCircle_mm (graphics, bodyX, bodyY, 1.0);
	Graphics_text (graphics, bodyX, bodyY, "#B");
	Graphics_text (graphics, (bodyX+intX[7])/2, (bodyY+intY[7])/2, "%%R__body");
	Graphics_fillCircle_mm (graphics, 0, 0, 1.0);
	Graphics_text (graphics, 0, 0, "#O");
	{
		double xd = bodyX - 20*f * cos (23*NUMpi/180), yd = bodyY - 20*f * sin (23*NUMpi/180);
		Graphics_fillCircle_mm (graphics, xd, yd, 1.0);
		Graphics_text (graphics, xd, yd, "#D");
		Graphics_line (graphics, bodyX, bodyY, xd, yd);
		Graphics_line (graphics, xd, yd, intX [4], intY [4]);
	}
	{
		Graphics_setTextAlignment (graphics, Graphics_LEFT, Graphics_BOTTOM);
		Graphics_text (graphics, intX[2], intY[2], "%%hyoid bone");
		Graphics_setTextAlignment (graphics, Graphics_CENTRE, Graphics_HALF);
		Graphics_setTextRotation (graphics, -70);
		Graphics_text (graphics, intX[2], (extY[1]+extY[2])/2, "%%larynx");
		Graphics_setTextRotation (graphics, 90);
		Graphics_setTextAlignment (graphics, Graphics_CENTRE, Graphics_BOTTOM);
		Graphics_text (graphics, extX[4], extY[4], "%%rear pharyngeal wall");
		Graphics_setTextAlignment (graphics, Graphics_LEFT, Graphics_TOP);
		Graphics_text (graphics, 0, 5*f, "%%R__palate");
		Graphics_setTextRotation (graphics, 40);
		Graphics_text (graphics, extX[6]+5*f, extY[6]+6*f, "%%velum");
		Graphics_setTextRotation (graphics, 0);
		{
			double palateRadius = f * sqrt (23*23 + 31*31);
			Graphics_text (graphics, 0, palateRadius, "%%palate");
			Graphics_line (graphics, 0, 0, 0, palateRadius);
		}
		Graphics_setTextAlignment (graphics, Graphics_CENTRE, Graphics_HALF);
		Graphics_text (graphics, -15*f, 20*f, "%%tongue");
		Graphics_text (graphics, -15*f, 16*f, "%%body");
		Graphics_text (graphics, intX[5]+10*f, intY[5]+2*f, "%%tongue");
		Graphics_text (graphics, intX[5]+10*f, intY[5]-2*f, "%%root");
		Graphics_text (graphics, 18*f, 16*f, "%%tongue");
		Graphics_text (graphics, 18*f, 12*f, "%%tip");
		Graphics_text (graphics, extX[11], extY[11]+11*f, "%%upper");
		Graphics_text (graphics, extX[11], extY[11]+7*f, "%%lip");
		Graphics_text (graphics, intX[12], intY[12]-22*f, "%%chin");
	}
	/* Points on outer contour. */
	{
		struct { int h, v; } align [] = { 0,0, Graphics_RIGHT, Graphics_HALF,
			Graphics_RIGHT, Graphics_TOP, Graphics_RIGHT, Graphics_HALF,
			Graphics_LEFT, Graphics_HALF, Graphics_RIGHT, Graphics_HALF,
			Graphics_CENTRE, Graphics_BOTTOM, Graphics_CENTRE, Graphics_BOTTOM,
			Graphics_LEFT, Graphics_HALF, Graphics_RIGHT, Graphics_BOTTOM,
			Graphics_CENTRE, Graphics_BOTTOM, Graphics_CENTRE, Graphics_BOTTOM };
		int i;
		for (i = 1; i <= 11; i ++) {
			Graphics_setTextAlignment (graphics, align [i]. h, align [i]. v);
			Graphics_printf (graphics, extX [i], extY [i], "%d", i);
		}
	}
	Graphics_fillCircle_mm (graphics, extX [4], extY [4], 1.0);
	/* Points on inner contour. */
	{
		struct { int h, v; } align [] = { 0,0, Graphics_RIGHT, Graphics_HALF,
			Graphics_RIGHT, Graphics_HALF, Graphics_RIGHT, Graphics_HALF,
			Graphics_CENTRE, Graphics_TOP, Graphics_RIGHT, Graphics_HALF,
			Graphics_RIGHT, Graphics_HALF, Graphics_CENTRE, Graphics_BOTTOM,
			Graphics_LEFT, Graphics_HALF, Graphics_RIGHT, Graphics_HALF,
			Graphics_CENTRE, Graphics_TOP, Graphics_LEFT, Graphics_HALF,
			Graphics_CENTRE, Graphics_TOP, Graphics_LEFT, Graphics_BOTTOM,
			Graphics_LEFT, Graphics_HALF };
		int i;
		for (i = 1; i <= 14; i ++) {
			Graphics_setTextAlignment (graphics, align [i]. h, align [i]. v);
			Graphics_printf (graphics, intX [i], intY [i], "%d", i);
		}
	}
	Graphics_resetViewport (graphics, oldViewport);
	Graphics_setLineType (graphics, Graphics_NORMAL);
	Graphics_setLineWidth (graphics, 4.0);
	Art_Speaker_draw (art, speaker, graphics);

	forget (art);
	forget (speaker);
}

/********** FIGURE 5.2: PUMPING AND SUCKING **********/

static void pumpingAndSucking (Graphics graphics) {
	float centre = 83.333333, w1 = 6, w2 = 5, w3 = 6, d = 5, x = 2 + d;
	Graphics_setWindow (graphics, 0, 100, 0, 100);
	/* Pump. */
	draw2Walls (graphics, x, centre, w1, d, 0);
	x += d;
	Graphics_line (graphics, x, centre - w1, x, centre - w2);
	Graphics_line (graphics, x, centre + w1, x, centre + w2);
	Graphics_arrow (graphics, x + 2, centre, x - 2, centre);
	x += d;
	draw2Walls (graphics, x, centre, w2, d, 0);
	Graphics_arrow (graphics, x, centre + w2, x, centre + 1);
	Graphics_arrow (graphics, x, centre - w2, x, centre - 1);
	x += d;
	Graphics_line (graphics, x, centre - w2, x, centre - w3);
	Graphics_line (graphics, x, centre + w2, x, centre + w3);
	Graphics_arrow (graphics, x - 2, centre, x + 2, centre);
	x += d;
	draw2Walls (graphics, x, centre, w3, d, 0);
	/* Suck. */
	x += 20;
	draw2Walls (graphics, x, centre, w1, d, 0);
	x += d;
	Graphics_line (graphics, x, centre - w1, x, centre - w2);
	Graphics_line (graphics, x, centre + w1, x, centre + w2);
	Graphics_arrow (graphics, x - 2, centre, x + 2, centre);
	x += d;
	draw2Walls (graphics, x, centre, w2, d, 0);
	Graphics_arrow (graphics, x+2, centre + w2, x+2, centre + w2 + 4);
	Graphics_arrow (graphics, x+2, centre - w2, x+2, centre - w2 - 4);
	x += d;
	Graphics_line (graphics, x, centre - w2, x, centre - w3);
	Graphics_line (graphics, x, centre + w2, x, centre + w3);
	Graphics_arrow (graphics, x + 2, centre, x - 2, centre);
	x += d;
	draw2Walls (graphics, x, centre, w3, d, 0);
}
static void tubeLabels (Graphics graphics, double x, double y, double angle, double d, int tube) {
	float dx, dy;
	Graphics_setTextRotation (graphics, angle);
	Graphics_setFontStyle (graphics, Graphics_ITALIC);
	Graphics_setTextAlignment (graphics, Graphics_LEFT, Graphics_HALF);
	dx = -d; dy = 0.0; rotate (& dx, & dy, angle);
	Graphics_text (graphics, x+dx, y+dy, "l");
	Graphics_setTextAlignment (graphics, Graphics_RIGHT, Graphics_HALF);
	dx = +d; dy = 0.0; rotate (& dx, & dy, angle);
	Graphics_text (graphics, x+dx, y+dy, "r");
	Graphics_setFontStyle (graphics, Graphics_NORMAL);
	Graphics_setTextAlignment (graphics, Graphics_CENTRE, Graphics_HALF);
	Graphics_printf (graphics, x, y, "%d", tube);
	Graphics_setTextRotation (graphics, 0.0);
}

/********** FIGURE 5.8-5.11: BOUNDARY TYPES **********/

static void closedBoundary (Graphics graphics) {
	float centre = 83.333333, x = 7, d = 5, w = 6;
	Graphics_setWindow (graphics, 0, 100, 0, 100);
	draw2Walls (graphics, x, centre, w, d, 0);   /* The tube. */
	Graphics_line (graphics, x-d, centre - w, x-d, centre + w);   /* Back vertical wall. */
	tubeLabels (graphics, x, centre, 0.0, d, 1);
}
static void openBoundary (Graphics graphics) {   /***** FIGURE 4.9 *****/
	float centre = 83.333333, x = 7, d = 5, w = 6;
	Graphics_setWindow (graphics, 0, 100, 0, 100);
	draw2Walls (graphics, x, centre, w, d, 0);   /* The tube. */
	Graphics_line (graphics, x+d, centre, x+d + 7.071, centre);   /* Straight ahead. */
	Graphics_arc (graphics, x+d, centre - 12, 10, 45, 90);   /* Bending down. */
	Graphics_arc (graphics, x+d, centre + 12, 10, 270, 315);   /* Bending up. */
	tubeLabels (graphics, x, centre, 0.0, d, 1);
}
static void twoWayBoundary (Graphics graphics) {   /***** FIGURE 4.10 *****/
	float centre = 83.333333, x = 7, d1 = 5, d2 = 4, w1 = 3, w2 = 6;
	Graphics_setWindow (graphics, 0, 100, 0, 100);
	draw2Walls (graphics, x, centre, w1, d1, 0);   /* Left tube. */
	tubeLabels (graphics, x, centre, 0.0, d1, 1);
	x += d1;
	Graphics_line (graphics, x, centre - w1, x, centre - w2);   /* Lower connection. */
	Graphics_line (graphics, x, centre + w1, x, centre + w2);   /* Upper connection. */
	x += d2;
	draw2Walls (graphics, x, centre, w2, d2, 0);   /* Right tube. */
	tubeLabels (graphics, x, centre, 0.0, d2, 2);
}
static void threeWayBoundary (Graphics graphics) {   /***** FIGURE 4.11 *****/
	float centre = 83.333333, x1 = 7, x2 = x1 + 15, x3 = x1 + 15, y2 = centre + 10, y3 = centre - 10;
	float d = 4, w1 = 4, w2 = 3, w3 = 5;
	float xb, yb, xe, ye;
	Graphics_setWindow (graphics, 0, 100, 0, 100);
	draw2Walls (graphics, x1, centre, w1, d, 0);   /* Left tube. */
	tubeLabels (graphics, x1, centre, 0.0, d, 1);
	draw2Walls (graphics, x2, y2, w2, d, 60);   /* Top right tube. */
	tubeLabels (graphics, x2, y2, 60.0, d, 2);
	draw2Walls (graphics, x3, y3, w3, d, -60);   /* Bottom right tube. */
	tubeLabels (graphics, x3, y3, -60.0, d, 3);
	xb = x1 + d; yb = centre + w1 - WALL_HWIDTH;
	xe = - d; ye = w2 - WALL_HWIDTH; rotate (& xe, & ye, 60);
	Graphics_line (graphics, xb, yb, x2 + xe, y2 + ye);   /* Connection from left to top. */
	xb = x1 + d; yb = centre - w1 + WALL_HWIDTH;
	xe = - d; ye = - w3 + WALL_HWIDTH; rotate (& xe, & ye, -60);
	Graphics_line (graphics, xb, yb, x3 + xe, y3 + ye);   /* Connection from left to bottom. */
	xb = - d; yb = - w2 + WALL_HWIDTH; rotate (& xb, & yb, 60);
	xe = - d; ye = w3 - WALL_HWIDTH; rotate (& xe, & ye, -60);
	Graphics_line (graphics, x2 + xb, y2 + yb, x3 + xe, y3 + ye);   /* From top to bottom. */
}
/********** CHAPTER 6 **********/
static double twopi = 2 * NUMpi;
#define pow2(x)  exp (x * NUMln2)
#define pi2  (NUMpi * NUMpi)

/* Windows centred around 0 with duration 1. */
static double wRectangular (double x)
	{ return fabs (x) < 0.5 ? 1 : 0; }
static double wTriangular (double x)   /* or Parzen */
	{ return fabs (x) < 0.5 ? 1 - fabs (2 * x) : 0; }
static double wParabolic (double x)   /* or Welch */
	{ return fabs (x) < 0.5 ? 1 - 4 * x * x : 0; }
static double wHanning (double x)
	{ return fabs (x) < 0.5 ? 0.5 + 0.5 * cos (twopi * x) : 0; }
static double wHamming (double x)
	{ return fabs (x) < 0.5 ? 0.54 + 0.46 * cos (twopi * x) : 0; }
static double wGauss (double x)
	{ return exp (- pi2 * x * x); }
static double wKaiser1 (double x)
	{ return fabs (x) < 0.5 ? NUMbesselI (0, (0.5*pi2+0.5) * sqrt (1 - 4*x*x)) / NUMbesselI (0, 0.5*pi2+0.5) : 0; }
static double wKaiser2 (double x)
	{ return fabs (x) < 1 ? NUMbesselI (0, (2*pi2+0.5) * sqrt (1 - x*x)) / NUMbesselI (0, 2*pi2+0.5) : 0; }
static void smoothingError (Any graphics, double (*window) (double), double windowArea) {
	int i;
	Polygon polygon = Polygon_create (4001);
	for (i = 1; i <= polygon -> numberOfPoints; i ++) {
		double octFT = (i - 1) / 1000.0, FT = pow2 (octFT);
		double dx = 1 / FT, x, sum, error, largestError = 0;
		/* Try phase 0. */
		sum = 0;
		for (x = 0; x > -10.0; x -= dx) sum += window (x);   /* Add negative values. */
		sum += sum - window (0);   /* Add positive values. */
		sum *= dx / windowArea;   /* Divide by window area (= windowArea / dx). */ 
		error = fabs (sum - 1);
		if (error > largestError) largestError = error;
		/* Try phase 0.5. */
		sum = 0;
		for (x = - 0.5 * dx; x > -10.0; x -= dx) sum += window (x);   /* Add negative values. */
		sum += sum;   /* Add positive values. */
		sum *= dx / windowArea;
		error = fabs (sum - 1);
		if (error > largestError) largestError = error;
		/* Try discontinuity at edge. */
		if (window == wHamming || window == wKaiser1 || window == wKaiser2) {
			double leftEdge = window == wKaiser2 ? -1 : -0.5, rightEdge = - leftEdge, save;
			/* Lower limit: edge not included. */
			sum = 0;
			for (x = leftEdge + dx; x < rightEdge; x += dx) sum += window (x);
			save = sum;
			sum *= dx / windowArea;
			error = fabs (sum - 1);
			if (error > largestError) largestError = error;
			/* Upper limit: edge included. */
			sum = save + window (leftEdge);
			sum *= dx / windowArea;
			error = fabs (sum - 1);
			if (error > largestError) largestError = error;
		}
		/* Convert to logarithmic scale. */
		largestError = - log10 (largestError);
		if (largestError > 10) largestError = 10;
		polygon -> x [i] = octFT;
		polygon -> y [i] = largestError;
	}
	Polygon_draw (polygon, graphics, 0.0, 4.0, 0.0, 10.0);
	forget (polygon);
}
static double erfc (double x)
{
	double t = 1 / (1 + 0.5 * fabs (x));
	double result = t * exp (- x * x - 1.26551223 + t * (1.00002368 + t * (0.37409196 +
		t * (0.09678418 + t * (-0.18628806 + t * (0.27886807 + t * (-1.13520398 +
		t * (1.48851587 + t * (-0.82215223 + t * 0.17087277)))))))));
	return x >= 0.0 ? result : 2 - result;
}
static void smoothingErrors (Graphics graphics) {   /***** FIGURE 6.3 *****/
	int i;
	smoothingError (graphics, wHanning, 0.5);
	smoothingError (graphics, wHamming, 0.54);
	smoothingError (graphics, wRectangular, 1);
	/*smoothingError (graphics, wParabolic, 2.0 / 3);*/
	smoothingError (graphics, wGauss, sqrt (1 / NUMpi));
	/*smoothingError (graphics, wGauss2,
		(sqrt (1 / NUMpi) * (1 - erfc (sqrt (NUMpi*NUMpi))) - 2 * exp (-NUMpi*NUMpi)) /
		(1 - exp (-NUMpi*NUMpi)));*/
	/*smoothingError (graphics, wGauss3,
		(sqrt (1 / NUMpi) * (1 - erfc (sqrt (2.25*NUMpi*NUMpi))) - 3 * exp (-2.25*NUMpi*NUMpi)) /
		(1 - exp (-2.25*NUMpi*NUMpi)));*/
	/*smoothingError (graphics, wKaiser1, 0.92866 * sqrt (1 / NUMpi));*/
	smoothingError (graphics, wKaiser2, 0.9813318062 * sqrt (1 / NUMpi));
	Graphics_drawInnerBox (graphics);
	Graphics_markBottom (graphics, 0.0, 0, 1, 0, "1");
	Graphics_markBottom (graphics, 1.0, 0, 1, 1, "2");
	Graphics_markBottom (graphics, 2.0, 0, 1, 1, "4");
	Graphics_markBottom (graphics, 3.0, 0, 1, 1, "8");
	Graphics_markBottom (graphics, 4.0, 0, 1, 0, "16");
	Graphics_textBottom (graphics, 1, "Number of periods per window");
	Graphics_markLeft (graphics, 0.0, 0, 1, 0, "1");
	for (i = 1; i <= 10; i ++) {
		char mark [10];
		sprintf (mark, "10^^-%d", i);
		Graphics_markLeft (graphics, i, 0, 1, 1, mark);
	}
	Graphics_textLeft (graphics, 1, "Maximum relative error");
	Graphics_textTop (graphics, 0, "Maximum error in smoothing a periodic signal");
}
static double rwHanning (double x) {
	if (x < -1 || x > 1) return 0;
	return (1-fabs (x)) * (2.0/3 + 1.0/3 * cos (twopi*x)) + sin (twopi*fabs(x)) / twopi;
}
static double rwHamming (double x) {
	if (x < -1 || x > 1) return 0;
	return ((1-fabs (x)) * (0.2916 + 0.1058 * cos (twopi*x)) + 0.3910 * sin (twopi*fabs(x)) / twopi) / 0.3974;
}
static double rwRectangular (double x) {
	if (x < -1 || x > 1) return 0;
	return 1 - fabs (x);
}
static double rwWelch (double x) {
	if (x < -1 || x > 1) return 0;
	return (1 - fabs (x)) * cos (NUMpi * x) + sin (NUMpi * fabs (x)) / NUMpi;
}
static double rwGauss (double x) {
	return exp (- NUMpi * NUMpi / 2 * x * x);
}
static void hnr_windowRipple (Any graphics, double (*w) (double), double (*rw) (double))
{
	int i;
	Polygon polygon = Polygon_create (4001);
	for (i = 1; i <= polygon -> numberOfPoints; i ++) {
		double octFT = (i - 1) * 0.001, FT = 3 * pow2 (octFT);
		double dx = 1 / FT, x, ra, ra0, rx, hnr, rx1, rx2;
		ra = ra0 = 0;
		for (x = 0.00; x > -10.0; x -= dx)
			{ ra += w (x) * w (x - dx); ra0 += w (x) * w (x); }
		ra *= 2, ra0 *= 2;
		ra0 -= w (0.00) * w (0.00);   /* This was counted twice. */
		rx1 = ra / ra0 / rw (dx);
		if (rx1 > 1) rx1 = 1 / rx1;
		ra = ra0 = 0;
		for (x = 0.00 - 0.5 * dx; x > -10.0; x -= dx)
			{ ra += w (x) * w (x - dx); ra0 += w (x) * w (x); }
		ra *= 2;
		ra += w (0.00 - 0.5 * dx) * w (0.00 + 0.5 * dx);
		ra0 *= 2;
		rx2 = ra / ra0 / rw (dx);
		if (rx2 > 1) rx2 = 1 / rx2;
		rx = rx1 < rx2 ? rx1 : rx2;
		hnr = rx >= 1 ? 200 : 10 * log10 (rx / (1 - rx));
		if (hnr > 150) { polygon -> numberOfPoints = i - 1; break; }
		polygon -> x [i] = octFT;
		polygon -> y [i] = hnr;
	}
	Polygon_draw (polygon, graphics, 0.0, 4.0, 0.0, 150.0);
	forget (polygon);
}
static void hnr_windowRipples (Any graphics)   /***** FIGURE 6.8 *****/
{
	hnr_windowRipple (graphics, wHanning, rwHanning);
	hnr_windowRipple (graphics, wHamming, rwHamming);
	hnr_windowRipple (graphics, wRectangular, rwRectangular);
	hnr_windowRipple (graphics, wGauss, rwGauss);
	Graphics_drawInnerBox (graphics);
	Graphics_markBottom (graphics, 0.0, 0, 1, 0, "3.0");
	Graphics_markBottom (graphics, 1.0, 0, 1, 1, "6.0");
	Graphics_markBottom (graphics, 2.0, 0, 1, 1, "12.0");
	Graphics_markBottom (graphics, 3.0, 0, 1, 1, "24.0");
	Graphics_markBottom (graphics, 4.0, 0, 1, 0, "48.0");
	Graphics_textBottom (graphics, 1, "Number of periods per window");
	Graphics_marksLeft (graphics, 6, 1, 1, 1);
	Graphics_textLeft (graphics, 1, "Minimum measured harmonics-to-noise ratio (dB)");
	Graphics_textTop (graphics, 1, "Worst-case resolution of HNR determination");
	Graphics_markRight (graphics, 36, 0, 0, 0, "rectangular");
	Graphics_markRight (graphics, 48, 0, 0, 0, "Hamming");
	Graphics_markRight (graphics, 78, 0, 0, 0, "Hanning");
	Graphics_markTop (graphics, 1.61, 0, 0, 0, "Gauss");
}

extern int Dissertation_draw (Any graphics, int chapter, int figure);
int Dissertation_draw (Any graphics, int chapter, int figure)
{
	oldLineType = Graphics_inqLineType (graphics);
	oldLineWidth = Graphics_inqLineWidth (graphics);
	oldFont = Graphics_inqFont (graphics);
	oldSize = Graphics_inqFontSize (graphics);
	oldStyle = Graphics_inqFontStyle (graphics);
	switch (chapter) {
		case 4: switch (figure) {
			case 1: manyTubes (graphics); break;
			case 4: tubeConnections (graphics); break;
			case 7: vocalTractConstruction (graphics); break;
			default: return 0;
		} break;
		case 5: switch (figure) {
			case 2: pumpingAndSucking (graphics); break;
			case 8: closedBoundary (graphics); break;
			case 9: openBoundary (graphics); break;
			case 10: twoWayBoundary (graphics); break;
			case 11: threeWayBoundary (graphics); break;
			default: return 0;
		} break;
		case 6: switch (figure) {
			case 3: smoothingErrors (graphics); break;
			case 8: hnr_windowRipples (graphics); break;
			default: return 0;
		} break;
		default: return 0;
	}
	Graphics_setLineType (graphics, oldLineType);
	Graphics_setLineWidth (graphics, oldLineWidth);
	Graphics_setFont (graphics, oldFont);
	Graphics_setFontSize (graphics, oldSize);
	Graphics_setFontStyle (graphics, oldStyle);
	return 1;
}

/* End of file Dissertation.c */