if ((c[i] == '\0') && isnum && wfe)		return(1);	else		return(0);}void InitPStream(PStream *pst){	void		InitDWin(PStream *);	double		*dcalloc(int, int);	double		**ddcalloc(int, int, int, int);	double		***dddcalloc(int, int, int, int, int, int);	int		half, full;	register int	i, m;	InitDWin(pst);	half = pst->range * 2;	full = pst->range * 4 + 1;	pst->vSize	= (pst->order + 1) * pst->dw.num;	pst->sm.length = LENGTH;	while (pst->sm.length < pst->range + pst->dw.maxw[WRIGHT])		pst->sm.length *= 2;	pst->mean	= dcalloc(pst->vSize*2, 0);	pst->ivar	= pst->mean + pst->vSize;	pst->sm.mseq	= ddcalloc(pst->sm.length, pst->vSize, 0, 0);	pst->sm.ivseq	= ddcalloc(pst->sm.length, pst->vSize, 0, 0);	pst->sm.c	= ddcalloc(pst->sm.length, pst->order+1, 0, 0);	pst->sm.P	= dddcalloc(full, pst->sm.length, pst->order+1, half, 0, 0);	pst->sm.pi	= ddcalloc(pst->range+pst->dw.maxw[WRIGHT]+1, pst->order+1, pst->range, 0);	pst->sm.k	= ddcalloc(pst->range+pst->dw.maxw[WRIGHT]+1, pst->order+1, pst->range, 0);	for (i = 0; i < pst->sm.length; i++)		for (m = 0; m < pst->vSize; m++)			pst->sm.ivseq[i][m] = 0.0;	for (i = 0; i < pst->sm.length; i++)		for (m = 0; m <= pst->order; m++)			pst->sm.P[0][i][m] = INFTY;	pst->sm.t = pst->range - 1;	pst->sm.mask = pst->sm.length - 1;}void InitDWin(PStream *pst){	double		*dcalloc(int, int);	int		str2darray(char *, double **);	register int	i, j;	int		fsize, leng;	double		x, s4, s2, s0;	FILE		*fp;	/* memory allocation */	if ((pst->dw.width = (int **) calloc(pst->dw.num, sizeof(int *))) == NULL) {		fprintf(stderr, "Cannot Allocate Memory\n");		exit(1);	}	for (i = 0; i < pst->dw.num; i++)		if ((pst->dw.width[i] = (int *) calloc(2, sizeof(int))) == NULL) {			fprintf(stderr, "Cannot Allocate Memory\n");			exit(1);		}	if ((pst->dw.coef = (double **) calloc(pst->dw.num, sizeof(double *))) == NULL) {		fprintf(stderr, "Cannot Allocate Memory\n");		exit(1);	}	/* window for static parameter */	pst->dw.width[0][WLEFT] = pst->dw.width[0][WRIGHT] = 0;	pst->dw.coef[0] = dcalloc(1, 0);	pst->dw.coef[0][0] = 1;	/* set delta coefficients */	if (pst->dw.calccoef == 0) {		for (i = 1; i < pst->dw.num; i++) {			if (pst->dw.fn[i][0] == ' ') {				fsize = str2darray(pst->dw.fn[i], &(pst->dw.coef[i]));			}			else {						/* read from file */				if ((fp = fopen(pst->dw.fn[i], "r")) == NULL) {					fprintf(stderr, "file %s not found\n", pst->dw.fn[i]);					exit(1);				}				/* check the number of coefficients */				fseek(fp, 0L, 2);				fsize = ftell(fp) / sizeof(real);				fseek(fp, 0L, 0);				/* read coefficients */				pst->dw.coef[i] = dcalloc(fsize, 0);				freadf(pst->dw.coef[i], sizeof(**(pst->dw.coef)), fsize, fp);			}			/* set pointer */			leng = fsize / 2;			pst->dw.coef[i] += leng;			pst->dw.width[i][WLEFT] = -leng;			pst->dw.width[i][WRIGHT] = leng;			if (fsize % 2 == 0)				pst->dw.width[i][WRIGHT]--;		}	}	else if (pst->dw.calccoef == 1) {		for (i = 1; i < pst->dw.num; i++) {			leng = atoi(pst->dw.fn[i]);			if (leng < 1) {				fprintf(stderr, "Width for regression coefficient shuould be more than 1.\n");				exit(1);			}			pst->dw.width[i][WLEFT] = -leng;			pst->dw.width[i][WRIGHT] = leng;			pst->dw.coef[i] = dcalloc(leng*2 + 1, 0);			pst->dw.coef[i] += leng;		}		leng = atoi(pst->dw.fn[1]);		s2 = 1;		for (j = 2; j <= leng; j++) {			x = j * j;			s2 += x;		}		s2 += s2;		for (j = -leng; j <= leng; j++)			pst->dw.coef[1][j] = j / s2;		if (pst->dw.num > 2) {			leng = atoi(pst->dw.fn[2]);			s2 = s4 = 1;			for (j = 2; j <= leng; j++) {				x = j * j;				s2 += x;				s4 += x * x;			}			s2 += s2;			s4 += s4;			s0 = leng + leng + 1;			for (j = -leng; j <= leng; j++)				pst->dw.coef[2][j] = (s0*j*j - s2)/(s4*s0 - s2*s2)/2;		}	}	pst->dw.maxw[WLEFT] = pst->dw.maxw[WRIGHT] = 0;	for (i = 0; i < pst->dw.num; i++) {		if (pst->dw.maxw[WLEFT] > pst->dw.width[i][WLEFT])			pst->dw.maxw[WLEFT] = pst->dw.width[i][WLEFT];		if (pst->dw.maxw[WRIGHT] < pst->dw.width[i][WRIGHT])			pst->dw.maxw[WRIGHT] = pst->dw.width[i][WRIGHT];	}}double *dcalloc(int x, int xoff){	double		*ptr;	if ((ptr = (double *) calloc(x, sizeof(*ptr))) == NULL) {		fprintf(stderr, "Cannot Allocate Memory\n");		exit(1);	}	ptr += xoff;	return(ptr);}double **ddcalloc(int x, int y, int xoff, int yoff){	double		*dcalloc(int, int);	double		**ptr;	register int	i;	if ((ptr = (double **) calloc(x, sizeof(*ptr))) == NULL) {		fprintf(stderr, "Cannot Allocate Memory\n");		exit(1);	}	for (i = 0; i < x; i++)		ptr[i] = dcalloc(y, yoff);	ptr += xoff;	return(ptr);}double ***dddcalloc(int x, int y, int z, int xoff, int yoff, int zoff){	double		**ddcalloc(int, int, int, int);	double		***ptr;	register int	i;	if ((ptr = (double ***) calloc(x, sizeof(*ptr))) == NULL) {		fprintf(stderr, "Cannot Allocate Memory\n");		exit(1);	}	for (i = 0; i < x; i++)		ptr[i] = ddcalloc(y, z, yoff, zoff);	ptr += xoff;	return(ptr);}int str2darray(char *c, double **x){	int	i, size, sp;	char	*p, *buf;	while (isspace(*c))		c++;	if (*c == '\0') {		*x = NULL;		return(0);	}	size = 1;	sp = 0;	for (p = c; *p != '\0'; p++) {		if (!isspace(*p)) {			if (sp == 1) {				size++;				sp = 0;			}		}		else			sp = 1;	}	buf = getmem(strlen(c), sizeof(*buf));	*x = dgetmem(size);	for (i = 0; i < size; i++)		(*x)[i] = strtod(c, &c);	return(size);}/*--------------------------------------------------------------------*/double *mlpg(PStream *pst){	int		doupdate(PStream *, int);	void		calc_pi(PStream *, int);	void		calc_k(PStream *, int);	void		update_P(PStream *, int);	void		update_c(PStream *, int);	int		tcur, tmin, tmax;	register int	d, m, u;	pst->sm.t++;	tcur = pst->sm.t & pst->sm.mask;	tmin = (pst->sm.t - pst->range) & pst->sm.mask;	tmax = (pst->sm.t + pst->dw.maxw[WRIGHT]) & pst->sm.mask;	for (u = -pst->range*2; u <= pst->range*2; u++) {		for (m = 0; m <= pst->order; m++)			pst->sm.P[u][tmax][m] = 0.0;	}	for (m = 0; m < pst->vSize; m++) {		pst->sm.mseq[tmax][m] = pst->mean[m];		pst->sm.ivseq[tmax][m] = pst->ivar[m];	}	for (m = 0; m <= pst->order; m++) {		if (pst->iType != 2)			pst->sm.c[tmax][m] = pst->mean[m];		else			pst->sm.c[tmax][m] = pst->mean[m] * finv(pst->ivar[m]);		pst->sm.P[0][tmax][m] = finv(pst->ivar[m]);	}	for (d = 1; d < pst->dw.num; d++) {		if (doupdate(pst, d)) {			calc_pi(pst, d);			calc_k(pst, d);			update_P(pst, d);			update_c(pst, d);		}	}	pst->par = pst->sm.c[tmin];	return(pst->par);}int doupdate(PStream *pst, int d){	register int	j;	if (pst->sm.ivseq[pst->sm.t&pst->sm.mask][(pst->order+1)*d] == 0.0)		return(0);	for (j = pst->dw.width[d][WLEFT]; j <= pst->dw.width[d][WRIGHT]; j++)		if (pst->sm.P[0][(pst->sm.t+j)&pst->sm.mask][0] == INFTY)			return(0);	return(1);}void calc_pi(PStream *pst, int d){	register int	j, m, u;	for (m = 0; m <= pst->order; m++)		for (u = -pst->range; u <= pst->dw.maxw[WRIGHT]; u++) {			pst->sm.pi[u][m] = 0.0;			for (j = pst->dw.width[d][WLEFT]; j <= pst->dw.width[d][WRIGHT]; j++)				pst->sm.pi[u][m] += pst->sm.P[u-j][(pst->sm.t+j)&pst->sm.mask][m] * pst->dw.coef[d][j];		}}void calc_k(PStream *pst, int d){	register int	j, m, u;	double		*ivar, x;	ivar = pst->sm.ivseq[pst->sm.t&pst->sm.mask] + (pst->order+1)*d;	for (m = 0; m <= pst->order; m++) {		x = 0.0;		for (j = pst->dw.width[d][WLEFT]; j <= pst->dw.width[d][WRIGHT]; j++)			x += pst->dw.coef[d][j] * pst->sm.pi[j][m];		x = ivar[m] / (1.0 + ivar[m] * x);		for (u = -pst->range; u <= pst->dw.maxw[WRIGHT]; u++) {			pst->sm.k[u][m] = pst->sm.pi[u][m] * x;		}	}}void update_P(PStream *pst, int d){	register int	m, u, v;	for (m = 0; m <= pst->order; m++)		for (u = -pst->range; u <= pst->dw.maxw[WRIGHT]; u++)			for (v = u; v <= pst->dw.maxw[WRIGHT]; v++) {				pst->sm.P[v-u][(pst->sm.t+u)&pst->sm.mask][m] -= pst->sm.k[v][m] * pst->sm.pi[u][m];				if (v != u)					pst->sm.P[u-v][(pst->sm.t+v)&pst->sm.mask][m] = pst->sm.P[v-u][(pst->sm.t+u)&pst->sm.mask][m];			}}void update_c(PStream *pst, int d){	register int	j, m, u;	double		*mean, x;	mean = pst->sm.mseq[pst->sm.t&pst->sm.mask] + (pst->order+1)*d;	for (m = 0; m <= pst->order; m++) {		x = mean[m];		for (j = pst->dw.width[d][WLEFT]; j <= pst->dw.width[d][WRIGHT]; j++)			x -= pst->dw.coef[d][j] * pst->sm.c[(pst->sm.t+j)&pst->sm.mask][m];		for (u = -pst->range; u <= pst->dw.maxw[WRIGHT]; u++)			pst->sm.c[(pst->sm.t+u)&pst->sm.mask][m] += pst->sm.k[u][m] * x;	}}