From e3494cb4bc2ee37b34da757eacb26bdb3772e283 Mon Sep 17 00:00:00 2001
From: "N. Brouard" <brouard@ined.fr>
Date: Wed, 18 Nov 2015 17:41:20 +0000
Subject: [PATCH] Summary: Start working on projected prevalences

---
 src/imach.c | 32 ++++++++++++++++++++++++--------
 1 file changed, 24 insertions(+), 8 deletions(-)

diff --git a/src/imach.c b/src/imach.c
index 0445b3c..d464720 100644
--- a/src/imach.c
+++ b/src/imach.c
@@ -1,6 +1,16 @@
 /* $Id$
   $State$
   $Log$
+  Revision 1.209  2015/11/17 22:12:03  brouard
+  Summary: Adding ftolpl parameter
+  Author: N Brouard
+
+  We had difficulties to get smoothed confidence intervals. It was due
+  to the period prevalence which wasn't computed accurately. The inner
+  parameter ftolpl is now an outer parameter of the .imach parameter
+  file after estepm. If ftolpl is small 1.e-4 and estepm too,
+  computation are long.
+
   Revision 1.208  2015/11/17 14:31:57  brouard
   Summary: temporary
 
@@ -3211,7 +3221,7 @@ void  freqsummary(char fileres[], int iagemin, int iagemax, int **s, double **ag
       k2cpt=0;
       for (i=1; i<=imx; i++) {
 	bool=1;
-	if  (cptcovn>0) { /* Filter is here: Must be looked at for model=V1+V2+V3+V4 */
+	if (cptcovn>0) { /* Filter is here: Must be looked at for model=V1+V2+V3+V4 */
 	  for (z1=1; z1<=cptcoveff; z1++)       
             if (covar[Tvaraff[z1]][i]!= nbcode[Tvaraff[z1]][codtabm(j1,z1)]){
                 /* Tests if the value of each of the covariates of i is equal to filter j1 */
@@ -3221,7 +3231,7 @@ void  freqsummary(char fileres[], int iagemin, int iagemax, int **s, double **ag
                 j1,z1,nbcode[Tvaraff[z1]][codtabm(j1,z1)],j1);*/
               /* For j1=7 in V1+V2+V3+V4 = 0 1 1 0 and codtabm(7,3)=1 and nbcde[3][?]=1*/
             } 
-	}
+	} /* cptcovn > 0 */
  
 	if (bool==1){
 	  for(m=firstpass; m<=lastpass; m++){
@@ -3235,14 +3245,15 @@ void  freqsummary(char fileres[], int iagemin, int iagemax, int **s, double **ag
 		freq[s[m][i]][s[m+1][i]][iagemax+3] += weight[i];
 	      }
 	      
-	      if ((agev[m][i]>1) && (agev[m][i]< (iagemax+3))) {
+	      if ((agev[m][i]>1) && (agev[m][i]< (iagemax+3)) && (anint[m][i]!=9999) && (mint[m][i]!=99)) {
 		dateintsum=dateintsum+k2;
 		k2cpt++;
+		/* printf("i=%ld dateintmean = %lf dateintsum=%lf k2cpt=%lf k2=%lf\n",i, dateintsum/k2cpt, dateintsum,k2cpt, k2); */
 	      }
 	      /*}*/
-	  }
-	}
-      } /* end i */
+	  } /* end m */
+	} /* end bool */
+      } /* end i = 1 to imx */
        
       /*      fprintf(ficresp, "#Count between %.lf/%.lf/%.lf and %.lf/%.lf/%.lf\n",jprev1, mprev1,anprev1,jprev2, mprev2,anprev2);*/
       pstamp(ficresp);
@@ -3328,9 +3339,9 @@ void  freqsummary(char fileres[], int iagemin, int iagemax, int **s, double **ag
 	if(first==1)
 	  printf("Others in log...\n");
 	fprintf(ficlog,"\n");
-      }
+      } /* end loop i */
       /*}*/
-  }
+  } /* end j1 */
   dateintmean=dateintsum/k2cpt; 
  
   fclose(ficresp);
@@ -4963,6 +4974,11 @@ void printinggnuplot(char fileresu[], char optionfilefiname[], double ageminpar,
     /*#endif */
   m=pow(2,cptcoveff);
 
+  /* Projected Prevalences */
+/* plot "NAGI0w_V1V2_monthlyb2b-proj/F_NAGI0w_V1V2_monthlyb2b-proj.txt" u 6:((($1 == 1) && ($2==0) && ($3==2) &&($4==0))? $7/(1-$13):1/0) t 'p11' w line */
+/* replot ""  u 6:((($1 == 1) && ($2==0) && ($3==2) &&($4==0))? $8/(1-$14):1/0) t 'p21' w line */
+/* replot ""  u 6:((($1 == 1) && ($2==0) && ($3==2) &&($4==0)&&($9!=0))? $9/(1-$15):1/0) t 'p.1' w line */
+
   /* Contribution to likelihood */
   /* Plot the probability implied in the likelihood */
     fprintf(ficgp,"\n# Contributions to the Likelihood, mle >=1. For mle=4 no interpolation, pure matrix products.\n#\n");
-- 
2.43.0