--- imach/src/imach.c 2003/05/02 18:51:41 1.74 +++ imach/src/imach.c 2003/06/05 15:41:51 1.81 @@ -1,4 +1,4 @@ -/* $Id: imach.c,v 1.74 2003/05/02 18:51:41 brouard Exp $ +/* $Id: imach.c,v 1.81 2003/06/05 15:41:51 brouard Exp $ Interpolated Markov Chain Short summary of the programme: @@ -119,7 +119,23 @@ #define ODIRSEPARATOR '\\' #endif -char version[80]="Imach version 0.95, February 2003, INED-EUROREVES "; +/* $Id: imach.c,v 1.81 2003/06/05 15:41:51 brouard Exp $ */ +/* $Log: imach.c,v $ + * Revision 1.81 2003/06/05 15:41:51 brouard + * *** empty log message *** + * +/* Revision 1.80 2003/06/05 15:34:14 brouard +/* Trying to add the true revision in the program and log +/* +/* Revision 1.79 2003/06/05 15:17:23 brouard +/* *** empty log message *** +/* */ +/* $Revision: 1.81 $ */ +/* $Date: 2003/06/05 15:41:51 $ */ +/* $State: Exp $ */ + +char version[80]="Imach version 0.95a1, June 2003, INED-EUROREVES "; +char fullversion[]="$Revision: 1.81 $ $Date: 2003/06/05 15:41:51 $"; int erreur; /* Error number */ int nvar; int cptcovn=0, cptcovage=0, cptcoveff=0,cptcov; @@ -141,7 +157,7 @@ double jmean; /* Mean space between 2 wa double **oldm, **newm, **savm; /* Working pointers to matrices */ double **oldms, **newms, **savms; /* Fixed working pointers to matrices */ FILE *fic,*ficpar, *ficparo,*ficres, *ficrespl, *ficrespij, *ficrest,*ficresf,*ficrespop; -FILE *ficlog; +FILE *ficlog, *ficrespow; FILE *ficgp,*ficresprob,*ficpop, *ficresprobcov, *ficresprobcor; FILE *ficresprobmorprev; FILE *fichtm; /* Html File */ @@ -327,6 +343,21 @@ void free_vector(double*v, int nl, int n } /************************ivector *******************************/ +char *cvector(long nl,long nh) +{ + char *v; + v=(char *) malloc((size_t)((nh-nl+1+NR_END)*sizeof(char))); + if (!v) nrerror("allocation failure in cvector"); + return v-nl+NR_END; +} + +/******************free ivector **************************/ +void free_cvector(char *v, long nl, long nh) +{ + free((FREE_ARG)(v+nl-NR_END)); +} + +/************************ivector *******************************/ int *ivector(long nl,long nh) { int *v; @@ -647,11 +678,15 @@ void powell(double p[], double **xi, int del=0.0; printf("\nPowell iter=%d -2*LL=%.12f",*iter,*fret); fprintf(ficlog,"\nPowell iter=%d -2*LL=%.12f",*iter,*fret); - for (i=1;i<=n;i++) + fprintf(ficrespow,"%d %.12f",*iter,*fret); + for (i=1;i<=n;i++) { printf(" %d %.12f",i, p[i]); - fprintf(ficlog," %d %.12f",i, p[i]); + fprintf(ficlog," %d %.12lf",i, p[i]); + fprintf(ficrespow," %.12lf", p[i]); + } printf("\n"); fprintf(ficlog,"\n"); + fprintf(ficrespow,"\n"); for (i=1;i<=n;i++) { for (j=1;j<=n;j++) xit[j]=xi[j][i]; fptt=(*fret); @@ -1187,14 +1222,27 @@ void mlikeli(FILE *ficres,double p[], in int i,j, iter; double **xi; double fret; + char filerespow[FILENAMELENGTH]; xi=matrix(1,npar,1,npar); for (i=1;i<=npar;i++) for (j=1;j<=npar;j++) xi[i][j]=(i==j ? 1.0 : 0.0); printf("Powell\n"); fprintf(ficlog,"Powell\n"); + strcpy(filerespow,"pow"); + strcat(filerespow,fileres); + if((ficrespow=fopen(filerespow,"w"))==NULL) { + printf("Problem with resultfile: %s\n", filerespow); + fprintf(ficlog,"Problem with resultfile: %s\n", filerespow); + } + fprintf(ficrespow,"# Powell\n# iter -2*LL"); + for (i=1;i<=nlstate;i++) + for(j=1;j<=nlstate+ndeath;j++) + if(j!=i)fprintf(ficrespow," p%1d%1d",i,j); + fprintf(ficrespow,"\n"); powell(p,xi,npar,ftol,&iter,&fret,func); - printf("\n#Number of iterations = %d, -2 Log likelihood = %.12f\n",iter,func(p)); + fclose(ficrespow); + printf("\n#Number of iterations = %d, -2 Log likelihood = %.12f\n",iter,func(p)); fprintf(ficlog,"\n#Number of iterations = %d, -2 Log likelihood = %.12f \n",iter,func(p)); fprintf(ficres,"#Number of iterations = %d, -2 Log likelihood = %.12f \n",iter,func(p)); @@ -1743,21 +1791,21 @@ void concatwav(int wav[], int **dh, int wav[i]=mi; if(mi==0){ if(first==0){ - printf("Warning, no any valid information for:%d line=%d and may be others, see log file\n",num[i],i); + printf("Warning! None valid information for:%d line=%d (skipped) and may be others, see log file\n",num[i],i); first=1; } if(first==1){ - fprintf(ficlog,"Warning, no any valid information for:%d line=%d\n",num[i],i); + fprintf(ficlog,"Warning! None valid information for:%d line=%d (skipped)\n",num[i],i); } } /* end mi==0 */ - } + } /* End individuals */ for(i=1; i<=imx; i++){ for(mi=1; mi nlstate) { + if (s[mw[mi+1][i]][i] > nlstate) { /* A death */ if (agedc[i] < 2*AGESUP) { j= rint(agedc[i]*12-agev[mw[mi][i]][i]*12); if(j==0) j=1; /* Survives at least one month after exam */ @@ -1765,9 +1813,9 @@ void concatwav(int wav[], int **dh, int if (j >= jmax) jmax=j; if (j <= jmin) jmin=j; sum=sum+j; - /*if (j<0) printf("j=%d num=%d \n",j,i); */ + /*if (j<0) printf("j=%d num=%d \n",j,i);*/ /* printf("%d %d %d %d\n", s[mw[mi][i]][i] ,s[mw[mi+1][i]][i],j,i);*/ - /*printf("%d %lf %d %d %d\n", i,agev[mw[mi][i]][i],j,s[mw[mi][i]][i] ,s[mw[mi+1][i]][i]);*/ + if(j<0)printf("Error! Negative delay (%d to death) between waves %d and %d of individual %d at line %d who is aged %.1f with statuses from %d to %d\n ",j,mw[mi][i],mw[mi+1][i],num[i], i,agev[mw[mi][i]][i],s[mw[mi][i]][i] ,s[mw[mi+1][i]][i]); } } else{ @@ -1778,6 +1826,7 @@ void concatwav(int wav[], int **dh, int else if (j <= jmin)jmin=j; /* if (j<10) printf("j=%d jmin=%d num=%d ",j,jmin,i); */ /*printf("%d %lf %d %d %d\n", i,agev[mw[mi][i]][i],j,s[mw[mi][i]][i] ,s[mw[mi+1][i]][i]);*/ + if(j<0)printf("Error! Negative delay (%d) between waves %d and %d of individual %d at line %d who is aged %.1f with statuses from %d to %d\n ",j,mw[mi][i],mw[mi+1][i],num[i], i,agev[mw[mi][i]][i],s[mw[mi][i]][i] ,s[mw[mi+1][i]][i]); sum=sum+j; } jk= j/stepm; @@ -1979,6 +2028,7 @@ void evsij(char fileres[], double ***eij for(i=1;i<=nlstate;i++){ cptj=cptj+1; for(h=0, gm[h][cptj]=0.; h<=nhstepm-1; h++){ + gm[h][cptj] = (p3mat[i][j][h]+p3mat[i][j][h+1])/2.; } } @@ -2792,7 +2842,7 @@ fprintf(fichtm," \n