--- imach/src/imach.c	2003/05/03 01:18:24	1.75
+++ imach/src/imach.c	2003/06/13 21:44:43	1.84
@@ -1,4 +1,21 @@
-/* $Id: imach.c,v 1.75 2003/05/03 01:18:24 brouard Exp $
+/* $Id: imach.c,v 1.84 2003/06/13 21:44:43 brouard Exp $
+  $State: Exp $
+  $Log: imach.c,v $
+  Revision 1.84  2003/06/13 21:44:43  brouard
+  * imach.c (Repository): Replace "freqsummary" at a correct
+  place. It differs from routine "prevalence" which may be called
+  many times. Probs is memory consuming and must be used with
+  parcimony.
+  Version 0.95a2 (should output exactly the same maximization than 0.8a2)
+
+  Revision 1.83  2003/06/10 13:39:11  lievre
+  *** empty log message ***
+
+  Revision 1.82  2003/06/05 15:57:20  brouard
+  Add log in  imach.c and  fullversion number is now printed.
+
+*/
+/*
    Interpolated Markov Chain
 
   Short summary of the programme:
@@ -58,6 +75,7 @@
   read parameterfile
   read datafile
   concatwav
+  freqsummary
   if (mle >= 1)
     mlikeli
   print results files
@@ -119,7 +137,11 @@
 #define ODIRSEPARATOR '\\'
 #endif
 
-char version[80]="Imach version 0.95, May 2003, INED-EUROREVES ";
+/* $Id: imach.c,v 1.84 2003/06/13 21:44:43 brouard Exp $ */
+/* $State: Exp $ */
+
+char version[]="Imach version 0.95a2, June 2003, INED-EUROREVES ";
+char fullversion[]="$Revision: 1.84 $ $Date: 2003/06/13 21:44:43 $"; 
 int erreur; /* Error number */
 int nvar;
 int cptcovn=0, cptcovage=0, cptcoveff=0,cptcov;
@@ -141,7 +163,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 +349,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 +684,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); 
@@ -1144,7 +1185,42 @@ double func( double *x)
 	ll[s[mw[mi][i]][i]] += 2*weight[i]*lli;
       } /* end of wave */
     } /* end of individual */
-  }else{  /* ml=4 no inter-extrapolation */
+  }else if (mle==4){  /* ml=4 no inter-extrapolation */
+    for (i=1,ipmx=0, sw=0.; i<=imx; i++){
+      for (k=1; k<=cptcovn;k++) cov[2+k]=covar[Tvar[k]][i];
+      for(mi=1; mi<= wav[i]-1; mi++){
+	for (ii=1;ii<=nlstate+ndeath;ii++)
+	  for (j=1;j<=nlstate+ndeath;j++){
+	    oldm[ii][j]=(ii==j ? 1.0 : 0.0);
+	    savm[ii][j]=(ii==j ? 1.0 : 0.0);
+	  }
+	for(d=0; d<dh[mi][i]; d++){
+	  newm=savm;
+	  cov[2]=agev[mw[mi][i]][i]+d*stepm/YEARM;
+	  for (kk=1; kk<=cptcovage;kk++) {
+	    cov[Tage[kk]+2]=covar[Tvar[Tage[kk]]][i]*cov[2];
+	  }
+	
+	  out=matprod2(newm,oldm,1,nlstate+ndeath,1,nlstate+ndeath,
+		       1,nlstate+ndeath,pmij(pmmij,cov,ncovmodel,x,nlstate));
+	  savm=oldm;
+	  oldm=newm;
+	} /* end mult */
+      
+	s1=s[mw[mi][i]][i];
+	s2=s[mw[mi+1][i]][i];
+	if( s2 > nlstate){ 
+	  lli=log(out[s1][s2] - savm[s1][s2]);
+	}else{
+	  lli=log(out[s[mw[mi][i]][i]][s[mw[mi+1][i]][i]]); /* Original formula */
+	}
+	ipmx +=1;
+	sw += weight[i];
+	ll[s[mw[mi][i]][i]] += 2*weight[i]*lli;
+	/*printf("i=%6d s1=%1d s2=%1d mi=%1d mw=%1d dh=%3d prob=%10.6f w=%6.4f out=%10.6f sav=%10.6f\n",i,s1,s2,mi,mw[mi][i],dh[mi][i],exp(lli),weight[i],out[s1][s2],savm[s1][s2]);*/
+      } /* end of wave */
+    } /* end of individual */
+  }else{  /* ml=5 no inter-extrapolation no jackson =0.8a */
     for (i=1,ipmx=0, sw=0.; i<=imx; i++){
       for (k=1; k<=cptcovn;k++) cov[2+k]=covar[Tvar[k]][i];
       for(mi=1; mi<= wav[i]-1; mi++){
@@ -1166,16 +1242,20 @@ double func( double *x)
 	  oldm=newm;
 	} /* end mult */
       
+	s1=s[mw[mi][i]][i];
+	s2=s[mw[mi+1][i]][i];
 	lli=log(out[s[mw[mi][i]][i]][s[mw[mi+1][i]][i]]); /* Original formula */
 	ipmx +=1;
 	sw += weight[i];
 	ll[s[mw[mi][i]][i]] += 2*weight[i]*lli;
+	/*printf("i=%6d s1=%1d s2=%1d mi=%1d mw=%1d dh=%3d prob=%10.6f w=%6.4f out=%10.6f sav=%10.6f\n",i,s1,s2,mi,mw[mi][i],dh[mi][i],exp(lli),weight[i],out[s1][s2],savm[s1][s2]);*/
       } /* end of wave */
     } /* end of individual */
   } /* End of if */
   for(k=1,l=0.; k<=nlstate; k++) l += ll[k];
   /* printf("l1=%f l2=%f ",ll[1],ll[2]); */
   l= l*ipmx/sw; /* To get the same order of magnitude as if weight=1 for every body */
+  /*exit(0); */
   return -l;
 }
 
@@ -1187,14 +1267,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));
 
@@ -1458,7 +1551,7 @@ void lubksb(double **a, int n, int *indx
 } 
 
 /************ Frequencies ********************/
-void  freqsummary(char fileres[], int iagemin, int iagemax, int **s, double **agev, int nlstate, int imx, int *Tvaraff, int **nbcode, int *ncodemax,double **mint,double **anint, double dateprev1,double dateprev2,double jprev1, double mprev1,double anprev1,double jprev2, double mprev2,double anprev2)
+void  freqsummary(char fileres[], int iagemin, int iagemax, int **s, double **agev, int nlstate, int imx, int *Tvaraff, int **nbcode, int *ncodemax,double **mint,double **anint)
 {  /* Some frequencies */
   
   int i, m, jk, k1,i1, j1, bool, z1,z2,j;
@@ -1512,7 +1605,7 @@ void  freqsummary(char fileres[], int ia
 	if (bool==1){
 	  for(m=firstpass; m<=lastpass; m++){
 	    k2=anint[m][i]+(mint[m][i]/12.);
-	    if ((k2>=dateprev1) && (k2<=dateprev2)) {
+	    /*if ((k2>=dateprev1) && (k2<=dateprev2)) {*/
 	      if(agev[m][i]==0) agev[m][i]=iagemax+1;
 	      if(agev[m][i]==1) agev[m][i]=iagemax+2;
 	      if (s[m][i]>0 && s[m][i]<=nlstate) prop[s[m][i]][(int)agev[m][i]] += weight[i];
@@ -1525,12 +1618,12 @@ void  freqsummary(char fileres[], int ia
 		dateintsum=dateintsum+k2;
 		k2cpt++;
 	      }
-	    }
+	      /*}*/
 	  }
 	}
       }
        
-      fprintf(ficresp, "#Count between %.lf/%.lf/%.lf and %.lf/%.lf/%.lf\n",jprev1, mprev1,anprev1,jprev2, mprev2,anprev2);
+      /*      fprintf(ficresp, "#Count between %.lf/%.lf/%.lf and %.lf/%.lf/%.lf\n",jprev1, mprev1,anprev1,jprev2, mprev2,anprev2);*/
 
       if  (cptcovn>0) {
 	fprintf(ficresp, "\n#********** Variable "); 
@@ -1591,7 +1684,7 @@ void  freqsummary(char fileres[], int ia
 	  if( i <= iagemax){
 	    if(pos>=1.e-5){
 	      fprintf(ficresp," %d %.5f %.0f %.0f",i,prop[jk][i]/posprop, prop[jk][i],posprop);
-	      probs[i][jk][j1]= pp[jk]/pos;
+	      /*probs[i][jk][j1]= pp[jk]/pos;*/
 	      /*printf("\ni=%d jk=%d j1=%d %.5f %.0f %.0f %f",i,jk,j1,pp[jk]/pos, pp[jk],pos,probs[i][jk][j1]);*/
 	    }
 	    else
@@ -1624,7 +1717,7 @@ void  freqsummary(char fileres[], int ia
 }
 
 /************ Prevalence ********************/
-void prevalence(double agemin, double agemax, int **s, double **agev, int nlstate, int imx, int *Tvar, int **nbcode, int *ncodemax,double **mint,double **anint, double dateprev1,double dateprev2, int firstpass, int lastpass)
+void prevalence(double ***probs, double agemin, double agemax, int **s, double **agev, int nlstate, int imx, int *Tvar, int **nbcode, int *ncodemax,double **mint,double **anint, double dateprev1,double dateprev2, int firstpass, int lastpass)
 {  
   /* Compute observed prevalence between dateprev1 and dateprev2 by counting the number of people
      in each health status at the date of interview (if between dateprev1 and dateprev2).
@@ -1743,21 +1836,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<wav[i];mi++){
       if (stepm <=0)
 	dh[mi][i]=1;
       else{
-	if (s[mw[mi+1][i]][i] > 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 +1858,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 +1871,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 +2073,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.;
 	  }
 	}
@@ -2356,7 +2451,7 @@ void varevsij(char optionfilefiname[], d
   fclose(ficresprobmorprev);
   fclose(ficgp);
   fclose(fichtm);
-}  
+}  /* end varevsij */
 
 /************ Variance of prevlim ******************/
 void varprevlim(char fileres[], double **varpl, double **matcov, double x[], double delti[], int nlstate, int stepm, double bage, double fage, double **oldm, double **savm, double **prlim, double ftolpl, int ij)
@@ -2792,7 +2887,7 @@ fprintf(fichtm," \n<ul><li><b>Graphs</b>
        fprintf(fichtm," ************\n<hr size=\"2\" color=\"#EC5E5E\">");
      }
      /* Pij */
-     fprintf(fichtm,"<br>- Pij or Conditional probabilities to be observed in state j being in state i %d (stepm) months before: pe%s%d1.png<br>
+     fprintf(fichtm,"<br>- Pij or Conditional probabilities to be observed in state j being in state i, %d (stepm) months before: pe%s%d1.png<br>
 <img src=\"pe%s%d1.png\">",stepm,strtok(optionfile, "."),jj1,strtok(optionfile, "."),jj1);     
      /* Quasi-incidences */
      fprintf(fichtm,"<br>- Pij or Conditional probabilities to be observed in state j being in state i %d (stepm) months before but expressed in per year i.e. quasi incidences if stepm is small and probabilities too: pe%s%d2.png<br>
@@ -2823,12 +2918,12 @@ health expectancies in states (1) and (2
  - Health expectancies with their variances (no covariance): <a href=\"t%s\">t%s</a> <br>\n
  - Standard deviation of stable prevalences: <a href=\"vpl%s\">vpl%s</a> <br>\n",rfileres,rfileres,fileres,fileres,fileres,fileres,fileres,fileres, estepm, fileres,fileres,fileres,fileres,fileres,fileres);
 
- if(popforecast==1) fprintf(fichtm,"\n
- - Prevalences forecasting: <a href=\"f%s\">f%s</a> <br>\n
- - Population forecasting (if popforecast=1): <a href=\"pop%s\">pop%s</a> <br>\n
-	<br>",fileres,fileres,fileres,fileres);
- else 
-   fprintf(fichtm,"\n No population forecast: popforecast = %d (instead of 1) or stepm = %d (instead of 1) or model=%s (instead of .)<br><br></li>\n",popforecast, stepm, model);
+/*  if(popforecast==1) fprintf(fichtm,"\n */
+/*  - Prevalences forecasting: <a href=\"f%s\">f%s</a> <br>\n */
+/*  - Population forecasting (if popforecast=1): <a href=\"pop%s\">pop%s</a> <br>\n */
+/* 	<br>",fileres,fileres,fileres,fileres); */
+/*  else  */
+/*    fprintf(fichtm,"\n No population forecast: popforecast = %d (instead of 1) or stepm = %d (instead of 1) or model=%s (instead of .)<br><br></li>\n",popforecast, stepm, model); */
 fprintf(fichtm," <ul><li><b>Graphs</b></li><p>");
 
  m=cptcoveff;
@@ -2845,7 +2940,7 @@ fprintf(fichtm," <ul><li><b>Graphs</b></
        fprintf(fichtm," ************\n<hr size=\"2\" color=\"#EC5E5E\">");
      }
      for(cpt=1; cpt<=nlstate;cpt++) {
-       fprintf(fichtm,"<br>- Observed and stationary prevalence (with confident
+       fprintf(fichtm,"<br>- Observed and period prevalence (with confident
 interval) in state (%d): v%s%d%d.png <br>
 <img src=\"v%s%d%d.png\">",cpt,strtok(optionfile, "."),cpt,jj1,strtok(optionfile, "."),cpt,jj1);  
      }
@@ -2946,9 +3041,9 @@ m=pow(2,cptcoveff);
     }
   }
   
-  /* CV preval stat */
+  /* CV preval stable (period) */
   for (k1=1; k1<= m ; k1 ++) { 
-    for (cpt=1; cpt<nlstate ; cpt ++) {
+    for (cpt=1; cpt<=nlstate ; cpt ++) {
       k=3;
       fprintf(ficgp,"\nset out \"p%s%d%d.png\" \n",strtok(optionfile, "."),cpt,k1);
       fprintf(ficgp,"set xlabel \"Age\" \nset ylabel \"Probability\" \nset ter png small\nset size 0.65,0.65\nplot [%.f:%.f] \"pij%s\" u ($1==%d ? ($3):1/0):($%d/($%d",ageminpar,agemaxpar,fileres,k1,k+cpt+1,k+1);
@@ -3092,7 +3187,7 @@ prevforecast(char fileres[], double anpr
   char fileresf[FILENAMELENGTH];
 
   agelim=AGESUP;
-  prevalence(ageminpar, agemax, s, agev, nlstate, imx, Tvar, nbcode, ncodemax, mint, anint, dateprev1, dateprev2, firstpass, lastpass);
+  prevalence(probs, ageminpar, agemax, s, agev, nlstate, imx, Tvar, nbcode, ncodemax, mint, anint, dateprev1, dateprev2, firstpass, lastpass);
  
   strcpy(fileresf,"f"); 
   strcat(fileresf,fileres);
@@ -3215,7 +3310,7 @@ populforecast(char fileres[], double anp
   agelim=AGESUP;
   calagedatem=(anpyram+mpyram/12.+jpyram/365.-dateintmean)*YEARM;
   
-  prevalence(ageminpar, agemax, s, agev, nlstate, imx, Tvar, nbcode, ncodemax, mint, anint, dateprev1, dateprev2, firstpass, lastpass);
+  prevalence(probs, ageminpar, agemax, s, agev, nlstate, imx, Tvar, nbcode, ncodemax, mint, anint, dateprev1, dateprev2, firstpass, lastpass);
   
   
   strcpy(filerespop,"pop"); 
@@ -3357,7 +3452,7 @@ populforecast(char fileres[], double anp
   free_ma3x(tabpop,1, AGESUP,1,NCOVMAX, 1,NCOVMAX);
   free_ma3x(tabpopprev,1, AGESUP,1,NCOVMAX, 1,NCOVMAX);
   fclose(ficrespop);
-}
+} /* End of popforecast */
 
 /***********************************************/
 /**************** Main Program *****************/
@@ -3420,7 +3515,7 @@ int main(int argc, char *argv[])
      gettimeofday(&start_time, (struct timezone*)0); */ /* at first time */
   getcwd(pathcd, size);
 
-  printf("\n%s",version);
+  printf("\n%s\n%s",version,fullversion);
   if(argc <=1){
     printf("\nEnter the parameter file name: ");
     scanf("%s",pathtot);
@@ -3786,11 +3881,20 @@ int main(int argc, char *argv[])
 
   for (i=1; i<=imx; i++) {
     for(m=2; (m<= maxwav); m++) {
-      if ((mint[m][i]== 99) && (s[m][i] <= nlstate)){
+      if (((int)mint[m][i]== 99) && (s[m][i] <= nlstate)){
 	anint[m][i]=9999;
 	s[m][i]=-1;
       }
-      if(moisdc[i]==99 && andc[i]==9999 & s[m][i]>nlstate) s[m][i]=-1;
+      if((int)moisdc[i]==99 && (int)andc[i]==9999 && s[m][i]>nlstate){
+	printf("Error! Date of death (month %2d and year %4d) of individual %d on line %d was unknown, you must set an arbitrary year of death or he/she is skipped and results are biased\n",(int)moisdc[i],(int)andc[i],num[i],i);
+	fprintf(ficlog,"Error! Date of death (month %2d and year %4d) of individual %d on line %d was unknown, you must set an arbitrary year of death or he/she is skipped and results are biased\n",(int)moisdc[i],(int)andc[i],num[i],i);
+	s[m][i]=-1;
+      }
+      if((int)moisdc[i]==99 && (int)andc[i]!=9999 && s[m][i]>nlstate){
+	printf("Error! Month of death of individual %d on line %d was unknown %2d, you should set it otherwise the information on the death is skipped and results are biased.\n",num[i],i,(int)moisdc[i]); 
+	fprintf(ficlog,"Error! Month of death of individual %d on line %d was unknown %f, you should set it otherwise the information on the death is skipped and results are biased.\n",num[i],i,moisdc[i]); 
+	s[m][i]=-1; /* We prefer to skip it (and to skip it in version 0.8a1 too */
+      }
     }
   }
 
@@ -3800,11 +3904,11 @@ int main(int argc, char *argv[])
       if(s[m][i] >0){
 	if (s[m][i] >= nlstate+1) {
 	  if(agedc[i]>0)
-	    if(moisdc[i]!=99 && andc[i]!=9999)
+	    if((int)moisdc[i]!=99 && (int)andc[i]!=9999)
 	      agev[m][i]=agedc[i];
 	  /*if(moisdc[i]==99 && andc[i]==9999) s[m][i]=-1;*/
 	    else {
-	      if (andc[i]!=9999){
+	      if ((int)andc[i]!=9999){
 		printf("Warning negative age at death: %d line:%d\n",num[i],i);
 		fprintf(ficlog,"Warning negative age at death: %d line:%d\n",num[i],i);
 		agev[m][i]=-1;
@@ -3815,7 +3919,7 @@ int main(int argc, char *argv[])
 				 years but with the precision of a
 				 month */
 	  agev[m][i]=(mint[m][i]/12.+1./24.+anint[m][i])-(moisnais[i]/12.+1./24.+annais[i]);
-	  if(mint[m][i]==99 || anint[m][i]==9999)
+	  if((int)mint[m][i]==99 || (int)anint[m][i]==9999)
 	    agev[m][i]=1;
 	  else if(agev[m][i] <agemin){ 
 	    agemin=agev[m][i];
@@ -3911,6 +4015,7 @@ int main(int argc, char *argv[])
     
   /* Calculates basic frequencies. Computes observed prevalence at single age
      and prints on file fileres'p'. */
+  freqsummary(fileres, agemin, agemax, s, agev, nlstate, imx,Tvaraff,nbcode, ncodemax,mint,anint);
 
     pmmij= matrix(1,nlstate+ndeath,1,nlstate+ndeath); /* creation */
     oldms= matrix(1,nlstate+ndeath,1,nlstate+ndeath); /* creation */
@@ -3954,7 +4059,7 @@ int main(int argc, char *argv[])
 	}
     }
   }
-  if(mle==1){
+  if(mle!=0){
     /* Computing hessian and covariance matrix */
     ftolhess=ftol; /* Usually correct */
     hesscov(matcov, p, npar, delti, ftolhess, func);
@@ -4085,8 +4190,8 @@ int main(int argc, char *argv[])
   fprintf(ficparo,"popforecast=%d popfile=%s popfiledate=%.lf/%.lf/%.lf last-popfiledate=%.lf/%.lf/%.lf\n",popforecast,popfile,jpyram,mpyram,anpyram,jpyram1,mpyram1,anpyram1);
   fprintf(ficres,"popforecast=%d popfile=%s popfiledate=%.lf/%.lf/%.lf last-popfiledate=%.lf/%.lf/%.lf\n",popforecast,popfile,jpyram,mpyram,anpyram,jpyram1,mpyram1,anpyram1);
 
-  probs= ma3x(1,AGESUP,1,NCOVMAX, 1,NCOVMAX);
-  freqsummary(fileres, agemin, agemax, s, agev, nlstate, imx,Tvaraff,nbcode, ncodemax,mint,anint,dateprev1,dateprev2,jprev1, mprev1,anprev1,jprev2, mprev2,anprev2);
+  freqsummary(fileres, agemin, agemax, s, agev, nlstate, imx,Tvaraff,nbcode, ncodemax,mint,anint);
+  /*,dateprev1,dateprev2,jprev1, mprev1,anprev1,jprev2, mprev2,anprev2);*/
 
   /*------------ gnuplot -------------*/
   strcpy(optionfilegnuplot,optionfilefiname);
@@ -4113,12 +4218,13 @@ int main(int argc, char *argv[])
 Title=%s <br>Datafile=%s Firstpass=%d Lastpass=%d Stepm=%d Weight=%d Model=%s<br>\n
 \n
 Total number of observations=%d <br>\n
+Youngest age at first (selected) pass %.2f, oldest age %.2f<br>\n
 Interval (in months) between two waves: Min=%d Max=%d Mean=%.2lf<br>\n
 <hr  size=\"2\" color=\"#EC5E5E\">
  <ul><li><h4>Parameter files</h4>\n
  - Copy of the parameter file: <a href=\"o%s\">o%s</a><br>\n
  - Log file of the run: <a href=\"%s\">%s</a><br>\n
- - Gnuplot file name: <a href=\"%s\">%s</a></ul>\n",version,title,datafile,firstpass,lastpass,stepm, weightopt,model,imx,jmin,jmax,jmean,fileres,fileres,filelog,filelog,optionfilegnuplot,optionfilegnuplot);
+ - Gnuplot file name: <a href=\"%s\">%s</a></ul>\n",version,title,datafile,firstpass,lastpass,stepm, weightopt,model,imx,agemin,agemax,jmin,jmax,jmean,fileres,fileres,filelog,filelog,optionfilegnuplot,optionfilegnuplot);
    fclose(fichtm);
 
   printinghtml(fileres,title,datafile, firstpass, lastpass, stepm, weightopt,model,imx,jmin,jmax,jmean,rfileres,popforecast,estepm,jprev1,mprev1,anprev1,jprev2,mprev2,anprev2);
@@ -4249,6 +4355,7 @@ Interval (in months) between two waves:
 
   fclose(ficrespij);
 
+  probs= ma3x(1,AGESUP,1,NCOVMAX, 1,NCOVMAX);
 
   /*---------- Forecasting ------------------*/
   /*if((stepm == 1) && (strcmp(model,".")==0)){*/
@@ -4296,7 +4403,7 @@ Interval (in months) between two waves:
   fprintf(ficlog,"Computing Variance-covariance of DFLEs: file '%s' \n", fileresv);
 
   /* Computes prevalence between agemin (i.e minimal age computed) and no more ageminpar */
-  prevalence(agemin, agemax, s, agev, nlstate, imx, Tvar, nbcode, ncodemax, mint, anint, dateprev1, dateprev2, firstpass, lastpass);
+  prevalence(probs, agemin, agemax, s, agev, nlstate, imx, Tvar, nbcode, ncodemax, mint, anint, dateprev1, dateprev2, firstpass, lastpass);
   /*  printf("ageminpar=%f, agemax=%f, s[lastpass][imx]=%d, agev[lastpass][imx]=%f, nlstate=%d, imx=%d,  mint[lastpass][imx]=%f, anint[lastpass][imx]=%f,dateprev1=%f, dateprev2=%f, firstpass=%d, lastpass=%d\n",\
 ageminpar, agemax, s[lastpass][imx], agev[lastpass][imx], nlstate, imx, mint[lastpass][imx],anint[lastpass][imx], dateprev1, dateprev2, firstpass, lastpass);
   */