--- imach/src/imach.c 2017/07/17 08:53:49 1.277
+++ imach/src/imach.c 2018/04/19 14:49:16 1.283
@@ -1,6 +1,26 @@
-/* $Id: imach.c,v 1.277 2017/07/17 08:53:49 brouard Exp $
+/* $Id: imach.c,v 1.283 2018/04/19 14:49:16 brouard Exp $
$State: Exp $
$Log: imach.c,v $
+ Revision 1.283 2018/04/19 14:49:16 brouard
+ Summary: Some minor bugs fixed
+
+ Revision 1.282 2018/02/27 22:50:02 brouard
+ *** empty log message ***
+
+ Revision 1.281 2018/02/27 19:25:23 brouard
+ Summary: Adding second argument for quitting
+
+ Revision 1.280 2018/02/21 07:58:13 brouard
+ Summary: 0.99r15
+
+ New Makefile with recent VirtualBox 5.26. Bug in sqrt negatve in imach.c
+
+ Revision 1.279 2017/07/20 13:35:01 brouard
+ Summary: temporary working
+
+ Revision 1.278 2017/07/19 14:09:02 brouard
+ Summary: Bug for mobil_average=0 and prevforecast fixed(?)
+
Revision 1.277 2017/07/17 08:53:49 brouard
Summary: BOM files can be read now
@@ -1034,12 +1054,12 @@ typedef struct {
#define ODIRSEPARATOR '\\'
#endif
-/* $Id: imach.c,v 1.277 2017/07/17 08:53:49 brouard Exp $ */
+/* $Id: imach.c,v 1.283 2018/04/19 14:49:16 brouard Exp $ */
/* $State: Exp $ */
#include "version.h"
char version[]=__IMACH_VERSION__;
-char copyright[]="February 2016,INED-EUROREVES-Institut de longevite-Japan Society for the Promotion of Science (Grant-in-Aid for Scientific Research 25293121), Intel Software 2015-2018";
-char fullversion[]="$Revision: 1.277 $ $Date: 2017/07/17 08:53:49 $";
+char copyright[]="April 2018,INED-EUROREVES-Institut de longevite-Japan Society for the Promotion of Science (Grant-in-Aid for Scientific Research 25293121), Intel Software 2015-2018";
+char fullversion[]="$Revision: 1.283 $ $Date: 2018/04/19 14:49:16 $";
char strstart[80];
char optionfilext[10], optionfilefiname[FILENAMELENGTH];
int erreur=0, nberr=0, nbwarn=0; /* Error number, number of errors number of warnings */
@@ -2513,15 +2533,18 @@ void powell(double p[], double **xi, int
double **prevalim(double **prlim, int nlstate, double x[], double age, double **oldm, double **savm, double ftolpl, int *ncvyear, int ij, int nres)
{
- /* Computes the prevalence limit in each live state at age x and for covariate combination ij
- (and selected quantitative values in nres)
- by left multiplying the unit
- matrix by transitions matrix until convergence is reached with precision ftolpl */
- /* Wx= Wx-1 Px-1= Wx-2 Px-2 Px-1 = Wx-n Px-n ... Px-2 Px-1 I */
- /* Wx is row vector: population in state 1, population in state 2, population dead */
- /* or prevalence in state 1, prevalence in state 2, 0 */
- /* newm is the matrix after multiplications, its rows are identical at a factor */
- /* Initial matrix pimij */
+ /**< Computes the prevalence limit in each live state at age x and for covariate combination ij
+ * (and selected quantitative values in nres)
+ * by left multiplying the unit
+ * matrix by transitions matrix until convergence is reached with precision ftolpl
+ * Wx= Wx-1 Px-1= Wx-2 Px-2 Px-1 = Wx-n Px-n ... Px-2 Px-1 I
+ * Wx is row vector: population in state 1, population in state 2, population dead
+ * or prevalence in state 1, prevalence in state 2, 0
+ * newm is the matrix after multiplications, its rows are identical at a factor.
+ * Inputs are the parameter, age, a tolerance for the prevalence limit ftolpl.
+ * Output is prlim.
+ * Initial matrix pimij
+ */
/* {0.85204250825084937, 0.13044499163996345, 0.017512500109187184, */
/* 0.090851990222114765, 0.88271245433047185, 0.026435555447413338, */
/* 0, 0 , 1} */
@@ -4348,7 +4371,7 @@ void freqsummary(char fileres[], double
double ***freq; /* Frequencies */
double *x, *y, a=0.,b=0.,r=1., sa=0., sb=0.; /* for regression, y=b+m*x and r is the correlation coefficient */
int no=0, linreg(int ifi, int ila, int *no, const double x[], const double y[], double* a, double* b, double* r, double* sa, double * sb);
- double *meanq;
+ double *meanq, *idq;
double **meanqt;
double *pp, **prop, *posprop, *pospropt;
double pos=0., posproptt=0., pospropta=0., k2, dateintsum=0,k2cpt=0;
@@ -4361,6 +4384,7 @@ void freqsummary(char fileres[], double
pospropt=vector(1,nlstate); /* Counting the number of transition starting from a live state */
/* prop=matrix(1,nlstate,iagemin,iagemax+3); */
meanq=vector(1,nqfveff); /* Number of Quantitative Fixed Variables Effective */
+ idq=vector(1,nqfveff); /* Number of Quantitative Fixed Variables Effective */
meanqt=matrix(1,lastpass,1,nqtveff);
strcpy(fileresp,"P_");
strcat(fileresp,fileresu);
@@ -4469,13 +4493,15 @@ Title=%s
Datafile=%s Firstpass=%d La
posprop[i]=0;
pospropt[i]=0;
}
- /* for (z1=1; z1<= nqfveff; z1++) { */
- /* meanq[z1]+=0.; */
+ for (z1=1; z1<= nqfveff; z1++) { /* zeroing for each combination j1 as well as for the total */
+ idq[z1]+=0.;
+ meanq[z1]+=0.;
+ }
+ /* for (z1=1; z1<= nqtveff; z1++) { */
/* for(m=1;m<=lastpass;m++){ */
- /* meanqt[m][z1]=0.; */
- /* } */
- /* } */
-
+ /* meanqt[m][z1]=0.; */
+ /* } */
+ /* } */
/* dateintsum=0; */
/* k2cpt=0; */
@@ -4485,7 +4511,7 @@ Title=%s
Datafile=%s Firstpass=%d La
if(j !=0){
if(anyvaryingduminmodel==0){ /* If All fixed covariates */
if (cptcoveff >0) { /* Filter is here: Must be looked at for model=V1+V2+V3+V4 */
- /* for (z1=1; z1<= nqfveff; z1++) { */
+ /* for (z1=1; z1<= nqfveff; z1++) { */
/* meanq[z1]+=coqvar[Tvar[z1]][iind]; /\* Computes mean of quantitative with selected filter *\/ */
/* } */
for (z1=1; z1<=cptcoveff; z1++) { /* loops on covariates in the model */
@@ -4551,6 +4577,11 @@ Title=%s
Datafile=%s Firstpass=%d La
/* freq[s[m][iind]][s[m+1][iind]][(int)((agebegin+ageend)/2.)] += weight[iind]; */
freq[s[m][iind]][s[m+1][iind]][iagemax+3] += weight[iind]; /* Total is in iagemax+3 *//* At age of beginning of transition, where status is known */
}
+ for (z1=1; z1<= nqfveff; z1++) {
+ idq[z1]++;
+ meanq[z1]+=covar[ncovcol+z1][iind]; /* *weight[iind];*/ /* Computes mean of quantitative with selected filter */
+ /* meanq[z1]+=coqvar[Tvar[z1]][iind]; /\* Computes mean of quantitative with selected filter *\/ */
+ }
} /* end if between passes */
if ((agev[m][iind]>1) && (agev[m][iind]< (iagemax+3)) && (anint[m][iind]!=9999) && (mint[m][iind]!=99) && (j==0)) {
dateintsum=dateintsum+k2; /* on all covariates ?*/
@@ -4598,6 +4629,19 @@ Title=%s
Datafile=%s Firstpass=%d La
fprintf(ficresphtmfr, "**********\n");
fprintf(ficlog, "**********\n");
}
+ /*
+ Printing means of quantitative variables if any
+ */
+ for (z1=1; z1<= nqfveff; z1++) {
+ fprintf(ficresphtmfr,"V quantitative id %d, number of idividuals= %f, sum=%f", z1, idq[z1], meanq[z1]);
+ fprintf(ficresphtmfr,", mean=%f
\n",meanq[z1]/idq[z1]);
+ }
+ /* for (z1=1; z1<= nqtveff; z1++) { */
+ /* for(m=1;m<=lastpass;m++){ */
+ /* fprintf(ficresphtmfr,"V quantitative id %d, pass id=%d, mean=%f
\n", z1, m, meanqt[m][z1]); */
+ /* } */
+ /* } */
+
fprintf(ficresphtm,"
");
if((cptcoveff==0 && nj==1)|| nj==2 ) /* no covariate and first pass */
fprintf(ficresp, " Age");
@@ -4878,6 +4922,7 @@ Title=%s
Datafile=%s Firstpass=%d La
fclose(ficresp);
fclose(ficresphtm);
fclose(ficresphtmfr);
+ free_vector(idq,1,nqfveff);
free_vector(meanq,1,nqfveff);
free_matrix(meanqt,1,lastpass,1,nqtveff);
free_vector(x, iagemin-AGEMARGE, iagemax+4+AGEMARGE);
@@ -5769,10 +5814,11 @@ void concatwav(int wav[], int **dh, int
/************ Variance ******************/
void varevsij(char optionfilefiname[], double ***vareij, double **matcov, double x[], double delti[], int nlstate, int stepm, double bage, double fage, double **oldm, double **savm, double **prlim, double ftolpl, int *ncvyearp, int ij, int estepm, int cptcov, int cptcod, int popbased, int mobilav, char strstart[], int nres)
{
- /* Variance of health expectancies */
- /* double **prevalim(double **prlim, int nlstate, double *xp, double age, double **oldm, double ** savm,double ftolpl);*/
- /* double **newm;*/
- /* int movingaverage(double ***probs, double bage,double fage, double ***mobaverage, int mobilav)*/
+ /** Variance of health expectancies
+ * double **prevalim(double **prlim, int nlstate, double *xp, double age, double **oldm, double ** savm,double ftolpl);
+ * double **newm;
+ * int movingaverage(double ***probs, double bage,double fage, double ***mobaverage, int mobilav)
+ */
/* int movingaverage(); */
double **dnewm,**doldm;
@@ -5780,11 +5826,11 @@ void concatwav(int wav[], int **dh, int
int i, j, nhstepm, hstepm, h, nstepm ;
int k;
double *xp;
- double **gp, **gm; /* for var eij */
- double ***gradg, ***trgradg; /*for var eij */
- double **gradgp, **trgradgp; /* for var p point j */
- double *gpp, *gmp; /* for var p point j */
- double **varppt; /* for var p point j nlstate to nlstate+ndeath */
+ double **gp, **gm; /**< for var eij */
+ double ***gradg, ***trgradg; /**< for var eij */
+ double **gradgp, **trgradgp; /**< for var p point j */
+ double *gpp, *gmp; /**< for var p point j */
+ double **varppt; /**< for var p point j nlstate to nlstate+ndeath */
double ***p3mat;
double age,agelim, hf;
/* double ***mobaverage; */
@@ -5845,7 +5891,7 @@ void concatwav(int wav[], int **dh, int
/* fprintf(fichtm, "#Local time at start: %s", strstart);*/
fprintf(fichtm,"\n Computing probabilities of dying over estepm months as a weighted average (i.e global mortality independent of initial healh state)
\n");
fprintf(fichtm,"\n
%s
\n",digitp);
- /* } */
+
varppt = matrix(nlstate+1,nlstate+ndeath,nlstate+1,nlstate+ndeath);
pstamp(ficresvij);
fprintf(ficresvij,"# Variance and covariance of health expectancies e.j \n# (weighted average of eij where weights are ");
@@ -5900,9 +5946,12 @@ void concatwav(int wav[], int **dh, int
for(i=1; i<=npar; i++){ /* Computes gradient x + delta*/
xp[i] = x[i] + (i==theta ?delti[theta]:0);
}
-
+ /**< Computes the prevalence limit with parameter theta shifted of delta up to ftolpl precision and
+ * returns into prlim .
+ */
prevalim(prlim,nlstate,xp,age,oldm,savm,ftolpl,ncvyearp,ij, nres);
-
+
+ /* If popbased = 1 we use crossection prevalences. Previous step is useless but prlim is created */
if (popbased==1) {
if(mobilav ==0){
for(i=1; i<=nlstate;i++)
@@ -5912,23 +5961,28 @@ void concatwav(int wav[], int **dh, int
prlim[i][i]=mobaverage[(int)age][i][ij];
}
}
-
- hpxij(p3mat,nhstepm,age,hstepm,xp,nlstate,stepm,oldm,savm, ij,nres); /* Returns p3mat[i][j][h] for h=1 to nhstepm */
+ /**< Computes the shifted transition matrix \f$ {}{h}_p^{ij}_x\f$ at horizon h.
+ */
+ hpxij(p3mat,nhstepm,age,hstepm,xp,nlstate,stepm,oldm,savm, ij,nres); /* Returns p3mat[i][j][h] for h=0 to nhstepm */
+ /**< And for each alive state j, sums over i \f$ w^i_x {}{h}_p^{ij}_x\f$, which are the probability
+ * at horizon h in state j including mortality.
+ */
for(j=1; j<= nlstate; j++){
for(h=0; h<=nhstepm; h++){
for(i=1, gp[h][j]=0.;i<=nlstate;i++)
gp[h][j] += prlim[i][i]*p3mat[i][j][h];
}
}
- /* Next for computing probability of death (h=1 means
+ /* Next for computing shifted+ probability of death (h=1 means
computed over hstepm matrices product = hstepm*stepm months)
- as a weighted average of prlim.
+ as a weighted average of prlim(i) * p(i,j) p.3=w1*p13 + w2*p23 .
*/
for(j=nlstate+1;j<=nlstate+ndeath;j++){
for(i=1,gpp[j]=0.; i<= nlstate; i++)
gpp[j] += prlim[i][i]*p3mat[i][j][1];
- }
- /* end probability of death */
+ }
+
+ /* Again with minus shift */
for(i=1; i<=npar; i++) /* Computes gradient x - delta */
xp[i] = x[i] - (i==theta ?delti[theta]:0);
@@ -5961,19 +6015,23 @@ void concatwav(int wav[], int **dh, int
for(i=1,gmp[j]=0.; i<= nlstate; i++)
gmp[j] += prlim[i][i]*p3mat[i][j][1];
}
- /* end probability of death */
-
+ /* end shifting computations */
+
+ /**< Computing gradient matrix at horizon h
+ */
for(j=1; j<= nlstate; j++) /* vareij */
for(h=0; h<=nhstepm; h++){
gradg[h][theta][j]= (gp[h][j]-gm[h][j])/2./delti[theta];
}
-
- for(j=nlstate+1; j<= nlstate+ndeath; j++){ /* var mu */
+ /**< Gradient of overall mortality p.3 (or p.j)
+ */
+ for(j=nlstate+1; j<= nlstate+ndeath; j++){ /* var mu mortality from j */
gradgp[theta][j]= (gpp[j]-gmp[j])/2./delti[theta];
}
} /* End theta */
-
+
+ /* We got the gradient matrix for each theta and state j */
trgradg =ma3x(0,nhstepm,1,nlstate,1,npar); /* veij */
for(h=0; h<=nhstepm; h++) /* veij */
@@ -5984,13 +6042,19 @@ void concatwav(int wav[], int **dh, int
for(j=nlstate+1; j<=nlstate+ndeath;j++) /* mu */
for(theta=1; theta <=npar; theta++)
trgradgp[j][theta]=gradgp[theta][j];
-
+ /**< as well as its transposed matrix
+ */
hf=hstepm*stepm/YEARM; /* Duration of hstepm expressed in year unit. */
for(i=1;i<=nlstate;i++)
for(j=1;j<=nlstate;j++)
vareij[i][j][(int)age] =0.;
-
+
+ /* Computing trgradg by matcov by gradg at age and summing over h
+ * and k (nhstepm) formula 15 of article
+ * Lievre-Brouard-Heathcote
+ */
+
for(h=0;h<=nhstepm;h++){
for(k=0;k<=nhstepm;k++){
matprod2(dnewm,trgradg[h],1,nlstate,1,npar,1,npar,matcov);
@@ -6001,7 +6065,11 @@ void concatwav(int wav[], int **dh, int
}
}
- /* pptj */
+ /* pptj is p.3 or p.j = trgradgp by cov by gradgp, variance of
+ * p.j overall mortality formula 49 but computed directly because
+ * we compute the grad (wix pijx) instead of grad (pijx),even if
+ * wix is independent of theta.
+ */
matprod2(dnewmp,trgradgp,nlstate+1,nlstate+ndeath,1,npar,1,npar,matcov);
matprod2(doldmp,dnewmp,nlstate+1,nlstate+ndeath,1,npar,nlstate+1,nlstate+ndeath,gradgp);
for(j=nlstate+1;j<=nlstate+ndeath;j++)
@@ -6610,7 +6678,12 @@ To be simple, these graphs help to under
}
/* Eigen vectors */
- v11=(1./sqrt(1+(v1-lc1)*(v1-lc1)/cv12/cv12));
+ if(1+(v1-lc1)*(v1-lc1)/cv12/cv12 <1.e-5){
+ printf(" Error sqrt of a negative number: %lf\n",1+(v1-lc1)*(v1-lc1)/cv12/cv12);
+ fprintf(ficlog," Error sqrt of a negative number: %lf\n",1+(v1-lc1)*(v1-lc1)/cv12/cv12);
+ v11=(1./sqrt(fabs(1+(v1-lc1)*(v1-lc1)/cv12/cv12)));
+ }else
+ v11=(1./sqrt(1+(v1-lc1)*(v1-lc1)/cv12/cv12));
/*v21=sqrt(1.-v11*v11); *//* error */
v21=(lc1-v1)/cv12*v11;
v12=-v21;
@@ -6641,8 +6714,8 @@ To be simple, these graphs help to under
fprintf(ficgp,"\nset label \"%d\" at %11.3e,%11.3e center",(int) age, mu1,mu2);
fprintf(ficgp,"\n# Age %d, p%1d%1d - p%1d%1d",(int) age, k1,l1,k2,l2);
fprintf(ficgp,"\nplot [-pi:pi] %11.3e+ %.3f*(%11.3e*%11.3e*cos(t)+%11.3e*%11.3e*sin(t)), %11.3e +%.3f*(%11.3e*%11.3e*cos(t)+%11.3e*%11.3e*sin(t)) not", \
- mu1,std,v11,sqrt(lc1),v12,sqrt(fabs(lc2)), \
- mu2,std,v21,sqrt(lc1),v22,sqrt(fabs(lc2))); /* For gnuplot only */
+ mu1,std,v11,sqrt(fabs(lc1)),v12,sqrt(fabs(lc2)), \
+ mu2,std,v21,sqrt(fabs(lc1)),v22,sqrt(fabs(lc2))); /* For gnuplot only */
}else{
first=0;
fprintf(fichtmcov," %d (%.3f),",(int) age, c12);
@@ -6817,7 +6890,7 @@ divided by h: hPij
for(cpt=1; cpt<=nlstate;cpt++){
fprintf(fichtm,"
\n- Survival functions from state %d in each live state and total.\
Or probability to survive in various states (1 to %d) being in state %d at different ages. \
- %s_%d%d-%d.svg
", cpt, nlstate, cpt, subdirf2(optionfilefiname,"LIJT_"),cpt,k1,nres,subdirf2(optionfilefiname,"LIJT_"),cpt,k1,nres,subdirf2(optionfilefiname,"LIJT_"),cpt,k1,nres);
+ %s_%d-%d-%d.svg
", cpt, nlstate, cpt, subdirf2(optionfilefiname,"LIJT_"),cpt,k1,nres,subdirf2(optionfilefiname,"LIJT_"),cpt,k1,nres,subdirf2(optionfilefiname,"LIJT_"),cpt,k1,nres);
}
/* Period (stable) prevalence in each health state */
for(cpt=1; cpt<=nlstate;cpt++){
@@ -7945,8 +8018,8 @@ set ter svg size 640, 480\nunset log y\n
double *agemingoodr, *agemaxgoodr;
- /* modcovmax=2*cptcoveff;/\* Max number of modalities. We suppose */
- /* a covariate has 2 modalities, should be equal to ncovcombmax *\/ */
+ /* modcovmax=2*cptcoveff; Max number of modalities. We suppose */
+ /* a covariate has 2 modalities, should be equal to ncovcombmax */
sumnewp = vector(1,ncovcombmax);
sumnewm = vector(1,ncovcombmax);
@@ -8270,11 +8343,11 @@ set ter svg size 640, 480\nunset log y\n
for(j=1; j<=nlstate+ndeath;j++) {
ppij=0.;
for(i=1; i<=nlstate;i++) {
- /* if (mobilav>=1) */
- ppij=ppij+p3mat[i][j][h]*prev[(int)agec][i][k];
- /* else { */ /* even if mobilav==-1 we use mobaverage */
- /* ppij=ppij+p3mat[i][j][h]*probs[(int)(agec)][i][k]; */
- /* } */
+ if (mobilav>=1)
+ ppij=ppij+p3mat[i][j][h]*prev[(int)agec][i][k];
+ else { /* even if mobilav==-1 we use mobaverage, probs may not sums to 1 */
+ ppij=ppij+p3mat[i][j][h]*probs[(int)(agec)][i][k];
+ }
fprintf(ficresf," %.3f", p3mat[i][j][h]);
} /* end i */
fprintf(ficresf," %.3f", ppij);
@@ -10753,8 +10826,13 @@ int main(int argc, char *argv[])
if(pathr[0] == '\0') break; /* Dirty */
}
}
+ else if (argc<=2){
+ strcpy(pathtot,argv[1]);
+ }
else{
strcpy(pathtot,argv[1]);
+ strcpy(z,argv[2]);
+ printf("\nargv[2]=%s z=%c\n",argv[2],z[0]);
}
/*if(getcwd(pathcd, MAXLINE)!= NULL)printf ("Error pathcd\n");*/
/*cygwin_split_path(pathtot,path,optionfile);
@@ -10832,8 +10910,6 @@ int main(int argc, char *argv[])
exit(70);
}
-
-
strcpy(filereso,"o");
strcat(filereso,fileresu);
if((ficparo=fopen(filereso,"w"))==NULL) { /* opened on subdirectory */
@@ -10842,7 +10918,20 @@ int main(int argc, char *argv[])
fflush(ficlog);
goto end;
}
+ /*-------- Rewriting parameter file ----------*/
+ strcpy(rfileres,"r"); /* "Rparameterfile */
+ strcat(rfileres,optionfilefiname); /* Parameter file first name */
+ strcat(rfileres,"."); /* */
+ strcat(rfileres,optionfilext); /* Other files have txt extension */
+ if((ficres =fopen(rfileres,"w"))==NULL) {
+ printf("Problem writing new parameter file: %s\n", rfileres);goto end;
+ fprintf(ficlog,"Problem writing new parameter file: %s\n", rfileres);goto end;
+ fflush(ficlog);
+ goto end;
+ }
+ fprintf(ficres,"#IMaCh %s\n",version);
+
/* Reads comments: lines beginning with '#' */
numlinepar=0;
/* Is it a BOM UTF-8 Windows file? */
@@ -10874,6 +10963,7 @@ int main(int argc, char *argv[])
numlinepar++;
fputs(line,stdout);
fputs(line,ficparo);
+ fputs(line,ficres);
fputs(line,ficlog);
continue;
}else
@@ -10883,18 +10973,24 @@ int main(int argc, char *argv[])
title, datafile, &lastobs, &firstpass,&lastpass)) !=EOF){
if (num_filled != 5) {
printf("Should be 5 parameters\n");
+ fprintf(ficlog,"Should be 5 parameters\n");
}
numlinepar++;
printf("title=%s datafile=%s lastobs=%d firstpass=%d lastpass=%d\n", title, datafile, lastobs, firstpass,lastpass);
+ fprintf(ficparo,"title=%s datafile=%s lastobs=%d firstpass=%d lastpass=%d\n", title, datafile, lastobs, firstpass,lastpass);
+ fprintf(ficres,"title=%s datafile=%s lastobs=%d firstpass=%d lastpass=%d\n", title, datafile, lastobs, firstpass,lastpass);
+ fprintf(ficlog,"title=%s datafile=%s lastobs=%d firstpass=%d lastpass=%d\n", title, datafile, lastobs, firstpass,lastpass);
}
/* Second parameter line */
while(fgets(line, MAXLINE, ficpar)) {
- /* If line starts with a # it is a comment */
+ /* while(fscanf(ficpar,"%[^\n]", line)) { */
+ /* If line starts with a # it is a comment. Strangely fgets reads the EOL and fputs doesn't */
if (line[0] == '#') {
numlinepar++;
- fputs(line,stdout);
- fputs(line,ficparo);
- fputs(line,ficlog);
+ printf("%s",line);
+ fprintf(ficres,"%s",line);
+ fprintf(ficparo,"%s",line);
+ fprintf(ficlog,"%s",line);
continue;
}else
break;
@@ -10904,8 +11000,13 @@ int main(int argc, char *argv[])
if (num_filled != 11) {
printf("Not 11 parameters, for example:ftol=1.e-8 stepm=12 ncovcol=2 nqv=1 ntv=2 nqtv=1 nlstate=2 ndeath=1 maxwav=3 mle=1 weight=1\n");
printf("but line=%s\n",line);
+ fprintf(ficlog,"Not 11 parameters, for example:ftol=1.e-8 stepm=12 ncovcol=2 nqv=1 ntv=2 nqtv=1 nlstate=2 ndeath=1 maxwav=3 mle=1 weight=1\n");
+ fprintf(ficlog,"but line=%s\n",line);
}
printf("ftol=%e stepm=%d ncovcol=%d nqv=%d ntv=%d nqtv=%d nlstate=%d ndeath=%d maxwav=%d mle=%d weight=%d\n",ftol, stepm, ncovcol, nqv, ntv, nqtv, nlstate, ndeath, maxwav, mle, weightopt);
+ fprintf(ficparo,"ftol=%e stepm=%d ncovcol=%d nqv=%d ntv=%d nqtv=%d nlstate=%d ndeath=%d maxwav=%d mle=%d weight=%d\n",ftol, stepm, ncovcol, nqv, ntv, nqtv, nlstate, ndeath, maxwav, mle, weightopt);
+ fprintf(ficres,"ftol=%e stepm=%d ncovcol=%d nqv=%d ntv=%d nqtv=%d nlstate=%d ndeath=%d maxwav=%d mle=%d weight=%d\n",ftol, stepm, ncovcol, nqv, ntv, nqtv, nlstate, ndeath, maxwav, mle, weightopt);
+ fprintf(ficlog,"ftol=%e stepm=%d ncovcol=%d nqv=%d ntv=%d nqtv=%d nlstate=%d ndeath=%d maxwav=%d mle=%d weight=%d\n",ftol, stepm, ncovcol, nqv, ntv, nqtv, nlstate, ndeath, maxwav, mle, weightopt);
}
/* ftolpl=6*ftol*1.e5; /\* 6.e-3 make convergences in less than 80 loops for the prevalence limit *\/ */
/*ftolpl=6.e-4; *//* 6.e-3 make convergences in less than 80 loops for the prevalence limit */
@@ -10914,22 +11015,18 @@ int main(int argc, char *argv[])
/* If line starts with a # it is a comment */
if (line[0] == '#') {
numlinepar++;
- fputs(line,stdout);
- fputs(line,ficparo);
- fputs(line,ficlog);
+ printf("%s",line);
+ fprintf(ficres,"%s",line);
+ fprintf(ficparo,"%s",line);
+ fprintf(ficlog,"%s",line);
continue;
}else
break;
}
if((num_filled=sscanf(line,"model=1+age%[^.\n]", model)) !=EOF){
- if (num_filled == 0){
- printf("ERROR %d: Model should be at minimum 'model=1+age.' WITHOUT space:'%s'\n",num_filled, line);
- fprintf(ficlog,"ERROR %d: Model should be at minimum 'model=1+age.' WITHOUT space:'%s'\n",num_filled, line);
- model[0]='\0';
- goto end;
- } else if (num_filled != 1){
- printf("ERROR %d: Model should be at minimum 'model=1+age.' %s\n",num_filled, line);
- fprintf(ficlog,"ERROR %d: Model should be at minimum 'model=1+age.' %s\n",num_filled, line);
+ if (num_filled != 1){
+ printf("ERROR %d: Model should be at minimum 'model=1+age' %s\n",num_filled, line);
+ fprintf(ficlog,"ERROR %d: Model should be at minimum 'model=1+age' %s\n",num_filled, line);
model[0]='\0';
goto end;
}
@@ -10942,20 +11039,23 @@ int main(int argc, char *argv[])
}
/* printf(" model=1+age%s modeltemp= %s, model=%s\n",model, modeltemp, model);fflush(stdout); */
printf("model=1+age+%s\n",model);fflush(stdout);
+ fprintf(ficparo,"model=1+age+%s\n",model);fflush(stdout);
+ fprintf(ficres,"model=1+age+%s\n",model);fflush(stdout);
+ fprintf(ficlog,"model=1+age+%s\n",model);fflush(stdout);
}
/* fscanf(ficpar,"title=%s datafile=%s lastobs=%d firstpass=%d lastpass=%d\nftol=%lf stepm=%d ncovcol=%d nlstate=%d ndeath=%d maxwav=%d mle=%d weight=%d model=1+age+%s\n",title, datafile, &lastobs, &firstpass,&lastpass,&ftol, &stepm, &ncovcol, &nlstate,&ndeath, &maxwav, &mle, &weightopt,model); */
/* numlinepar=numlinepar+3; /\* In general *\/ */
/* printf("title=%s datafile=%s lastobs=%d firstpass=%d lastpass=%d\nftol=%e stepm=%d ncovcol=%d nlstate=%d ndeath=%d maxwav=%d mle=%d weight=%d\nmodel=1+age+%s\n", title, datafile, lastobs, firstpass,lastpass,ftol, stepm, ncovcol, nlstate,ndeath, maxwav, mle, weightopt,model); */
- fprintf(ficparo,"title=%s datafile=%s lastobs=%d firstpass=%d lastpass=%d\nftol=%e stepm=%d ncovcol=%d nqv=%d ntv=%d nqtv=%d nlstate=%d ndeath=%d maxwav=%d mle=%d weight=%d\nmodel=1+age+%s.\n", title, datafile, lastobs, firstpass,lastpass,ftol,stepm,ncovcol, nqv, ntv, nqtv, nlstate,ndeath,maxwav, mle, weightopt,model);
- fprintf(ficlog,"title=%s datafile=%s lastobs=%d firstpass=%d lastpass=%d\nftol=%e stepm=%d ncovcol=%d nqv=%d ntv=%d nqtv=%d nlstate=%d ndeath=%d maxwav=%d mle=%d weight=%d\nmodel=1+age+%s.\n", title, datafile, lastobs, firstpass,lastpass,ftol,stepm,ncovcol, nqv, ntv, nqtv, nlstate,ndeath,maxwav, mle, weightopt,model);
+ /* fprintf(ficparo,"title=%s datafile=%s lastobs=%d firstpass=%d lastpass=%d\nftol=%e stepm=%d ncovcol=%d nqv=%d ntv=%d nqtv=%d nlstate=%d ndeath=%d maxwav=%d mle=%d weight=%d\nmodel=1+age+%s.\n", title, datafile, lastobs, firstpass,lastpass,ftol,stepm,ncovcol, nqv, ntv, nqtv, nlstate,ndeath,maxwav, mle, weightopt,model); */
+ /* fprintf(ficlog,"title=%s datafile=%s lastobs=%d firstpass=%d lastpass=%d\nftol=%e stepm=%d ncovcol=%d nqv=%d ntv=%d nqtv=%d nlstate=%d ndeath=%d maxwav=%d mle=%d weight=%d\nmodel=1+age+%s.\n", title, datafile, lastobs, firstpass,lastpass,ftol,stepm,ncovcol, nqv, ntv, nqtv, nlstate,ndeath,maxwav, mle, weightopt,model); */
fflush(ficlog);
/* if(model[0]=='#'|| model[0]== '\0'){ */
if(model[0]=='#'){
- printf("Error in 'model' line: model should start with 'model=1+age+' and end with '.' \n \
- 'model=1+age+.' or 'model=1+age+V1.' or 'model=1+age+age*age+V1+V1*age.' or \n \
- 'model=1+age+V1+V2.' or 'model=1+age+V1+V2+V1*V2.' etc. \n"); \
+ printf("Error in 'model' line: model should start with 'model=1+age+' and end without space \n \
+ 'model=1+age+' or 'model=1+age+V1.' or 'model=1+age+age*age+V1+V1*age' or \n \
+ 'model=1+age+V1+V2' or 'model=1+age+V1+V2+V1*V2' etc. \n"); \
if(mle != -1){
- printf("Fix the model line and run imach with mle=-1 to get a correct template of the parameter file.\n");
+ printf("Fix the model line and run imach with mle=-1 to get a correct template of the parameter vectors and subdiagonal covariance matrix.\n");
exit(1);
}
}
@@ -11177,16 +11277,6 @@ Please run with mle=-1 to get a correct
fflush(ficlog);
- /*-------- Rewriting parameter file ----------*/
- strcpy(rfileres,"r"); /* "Rparameterfile */
- strcat(rfileres,optionfilefiname); /* Parameter file first name*/
- strcat(rfileres,"."); /* */
- strcat(rfileres,optionfilext); /* Other files have txt extension */
- if((ficres =fopen(rfileres,"w"))==NULL) {
- printf("Problem writing new parameter file: %s\n", rfileres);goto end;
- fprintf(ficlog,"Problem writing new parameter file: %s\n", rfileres);goto end;
- }
- fprintf(ficres,"#%s\n",version);
} /* End of mle != -3 */
/* Main data
@@ -11820,7 +11910,7 @@ Please run with mle=-1 to get a correct
printf("\n");
/*--------- results files --------------*/
- fprintf(ficres,"title=%s datafile=%s lastobs=%d firstpass=%d lastpass=%d\nftol=%e stepm=%d ncovcol=%d nqv=%d ntv=%d nqtv=%d nlstate=%d ndeath=%d maxwav=%d mle= 0 weight=%d\nmodel=1+age+%s.\n", title, datafile, lastobs, firstpass,lastpass,ftol, stepm, ncovcol, nqv, ntv, nqtv, nlstate, ndeath, maxwav, weightopt,model);
+ /* fprintf(ficres,"title=%s datafile=%s lastobs=%d firstpass=%d lastpass=%d\nftol=%e stepm=%d ncovcol=%d nqv=%d ntv=%d nqtv=%d nlstate=%d ndeath=%d maxwav=%d mle= 0 weight=%d\nmodel=1+age+%s.\n", title, datafile, lastobs, firstpass,lastpass,ftol, stepm, ncovcol, nqv, ntv, nqtv, nlstate, ndeath, maxwav, weightopt,model); */
fprintf(ficres,"# Parameters nlstate*nlstate*ncov a12*1 + b12 * age + ...\n");
@@ -12618,6 +12708,8 @@ Please run with mle=-1 to get a correct
fclose(ficlog);
/*------ End -----------*/
+
+/* Executes gnuplot */
printf("Before Current directory %s!\n",pathcd);
#ifdef WIN32
@@ -12686,4 +12778,6 @@ end:
printf("\nType q for exiting: "); fflush(stdout);
scanf("%s",z);
}
+ printf("End\n");
+ exit(0);
}