Diff for /imach/src/imach.c between versions 1.237 and 1.238

version 1.237, 2016/08/26 09:20:19 version 1.238, 2016/08/26 14:23:35
Line 1 Line 1
 /* $Id$  /* $Id$
   $State$    $State$
   $Log$    $Log$
     Revision 1.238  2016/08/26 14:23:35  brouard
     Summary: Starting tests of 0.99
   
   Revision 1.237  2016/08/26 09:20:19  brouard    Revision 1.237  2016/08/26 09:20:19  brouard
   Summary: to valgrind    Summary: to valgrind
   
Line 1145  double *Tvalsel; /**< Selected modality Line 1148  double *Tvalsel; /**< Selected modality
 int *Typevar; /**< 0 for simple covariate (dummy, quantitative, fixed or varying), 1 for age product, 2 for  product */  int *Typevar; /**< 0 for simple covariate (dummy, quantitative, fixed or varying), 1 for age product, 2 for  product */
 int *Fixed; /** Fixed[k] 0=fixed, 1 varying, 2 fixed with age product, 3 varying with age product */   int *Fixed; /** Fixed[k] 0=fixed, 1 varying, 2 fixed with age product, 3 varying with age product */ 
 int *Dummy; /** Dummy[k] 0=dummy (0 1), 1 quantitative (single or product without age), 2 dummy with age product, 3 quant with age product */   int *Dummy; /** Dummy[k] 0=dummy (0 1), 1 quantitative (single or product without age), 2 dummy with age product, 3 quant with age product */ 
   int *DummyV; /** Dummy[v] 0=dummy (0 1), 1 quantitative */
   int *FixedV; /** FixedV[v] 0 fixed, 1 varying */
 int *Tage;  int *Tage;
 int anyvaryingduminmodel=0; /**< Any varying dummy in Model=1 yes, 0 no, to avoid a loop on waves in freq */   int anyvaryingduminmodel=0; /**< Any varying dummy in Model=1 yes, 0 no, to avoid a loop on waves in freq */ 
 int *Tmodelind; /** Tmodelind[Tvaraff[3]]=9 for V1 position,Tvaraff[1]@9={4, 3, 1, 0, 0, 0, 0, 0, 0}, model=V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1*/  int *Tmodelind; /** Tmodelind[Tvaraff[3]]=9 for V1 position,Tvaraff[1]@9={4, 3, 1, 0, 0, 0, 0, 0, 0}, model=V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1*/
Line 1154  int *Ndum; /** Freq of modality (tricode Line 1159  int *Ndum; /** Freq of modality (tricode
 /* int **codtab;*/ /**< codtab=imatrix(1,100,1,10); */  /* int **codtab;*/ /**< codtab=imatrix(1,100,1,10); */
 int **Tvard;  int **Tvard;
 int *Tprod;/**< Gives the k position of the k1 product */  int *Tprod;/**< Gives the k position of the k1 product */
   /* Tprod[k1=1]=3(=V1*V4) for V2+V1+V1*V4+age*V3  */
 int *Tposprod; /**< Gives the k1 product from the k position */  int *Tposprod; /**< Gives the k1 product from the k position */
 /* Tprod[k1=1]=3(=V1*V4) for V2+V1+V1*V4+age*V3     /* if  V2+V1+V1*V4+age*V3+V3*V2   TProd[k1=2]=5 (V3*V2) */
    if  V2+V1+V1*V4+age*V3+V3*V2   TProd[k1=2]=5 (V3*V2)     /* Tposprod[k]=k1 , Tposprod[3]=1, Tposprod[5(V3*V2)]=2 (2nd product without age) */
    Tposprod[k]=k1 , Tposprod[3]=1, Tposprod[5]=2   
 */  
 int cptcovprod, *Tvaraff, *invalidvarcomb;  int cptcovprod, *Tvaraff, *invalidvarcomb;
 double *lsurv, *lpop, *tpop;  double *lsurv, *lpop, *tpop;
   
Line 4918  void  concatwav(int wav[], int **dh, int Line 4922  void  concatwav(int wav[], int **dh, int
   double ***p3mat;    double ***p3mat;
   double eip;    double eip;
   
   pstamp(ficreseij);    /* pstamp(ficreseij); */
   fprintf(ficreseij,"# (a) Life expectancies by health status at initial age and (b) health expectancies by health status at initial age\n");    fprintf(ficreseij,"# (a) Life expectancies by health status at initial age and (b) health expectancies by health status at initial age\n");
   fprintf(ficreseij,"# Age");    fprintf(ficreseij,"# Age");
   for(i=1; i<=nlstate;i++){    for(i=1; i<=nlstate;i++){
Line 5283  void  concatwav(int wav[], int **dh, int Line 5287  void  concatwav(int wav[], int **dh, int
    fprintf(ficlog,"Computing total mortality p.j=w1*p1j+w2*p2j+..: result on file '%s' \n",fileresprobmorprev);     fprintf(ficlog,"Computing total mortality p.j=w1*p1j+w2*p2j+..: result on file '%s' \n",fileresprobmorprev);
    pstamp(ficresprobmorprev);     pstamp(ficresprobmorprev);
    fprintf(ficresprobmorprev,"# probabilities of dying before estepm=%d months for people of exact age and weighted probabilities w1*p1j+w2*p2j+... stand dev in()\n",estepm);     fprintf(ficresprobmorprev,"# probabilities of dying before estepm=%d months for people of exact age and weighted probabilities w1*p1j+w2*p2j+... stand dev in()\n",estepm);
      fprintf(ficresprobmorprev,"# Selected quantitative variables and dummies");
      for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */
        fprintf(ficresprobmorprev," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
      }
      for(j=1;j<=cptcoveff;j++) 
        fprintf(ficresprobmorprev,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(ij,j)]);
      fprintf(ficresprobmorprev,"\n");
   
    fprintf(ficresprobmorprev,"# Age cov=%-d",ij);     fprintf(ficresprobmorprev,"# Age cov=%-d",ij);
    for(j=nlstate+1; j<=(nlstate+ndeath);j++){     for(j=nlstate+1; j<=(nlstate+ndeath);j++){
      fprintf(ficresprobmorprev," p.%-d SE",j);       fprintf(ficresprobmorprev," p.%-d SE",j);
Line 6270  void printinggnuplot(char fileresu[], ch Line 6282  void printinggnuplot(char fileresu[], ch
   strcpy(dirfileres,optionfilefiname);    strcpy(dirfileres,optionfilefiname);
   strcpy(optfileres,"vpl");    strcpy(optfileres,"vpl");
   /* 1eme*/    /* 1eme*/
   for (cpt=1; cpt<= nlstate ; cpt ++) { /* For each live state */    for (cpt=1; cpt<= nlstate ; cpt ++){ /* For each live state */
     for (k1=1; k1<= m ; k1 ++) /* For each valid combination of covariate */      for (k1=1; k1<= m ; k1 ++){ /* For each valid combination of covariate */
       for(nres=1; nres <= nresult; nres++){ /* For each resultline */        for(nres=1; nres <= nresult; nres++){ /* For each resultline */
     /* plot [100000000000000000000:-100000000000000000000] "mysbiaspar/vplrmysbiaspar.txt to check */          /* plot [100000000000000000000:-100000000000000000000] "mysbiaspar/vplrmysbiaspar.txt to check */
       if(TKresult[nres]!= k1)          if(TKresult[nres]!= k1)
         continue;            continue;
       /* We are interested in selected combination by the resultline */          /* We are interested in selected combination by the resultline */
       printf("\n# 1st: Period (stable) prevalence with CI: 'VPL_' files and live state =%d ", cpt);          printf("\n# 1st: Period (stable) prevalence with CI: 'VPL_' files and live state =%d ", cpt);
       fprintf(ficgp,"\n# 1st: Period (stable) prevalence with CI: 'VPL_' files  and live state =%d ", cpt);          fprintf(ficgp,"\n# 1st: Period (stable) prevalence with CI: 'VPL_' files  and live state =%d ", cpt);
       for (k=1; k<=cptcoveff; k++){    /* For each covariate k get corresponding value lv for combination k1 */          for (k=1; k<=cptcoveff; k++){    /* For each covariate k get corresponding value lv for combination k1 */
         lv= decodtabm(k1,k,cptcoveff); /* Should be the value of the covariate corresponding to k1 combination */            lv= decodtabm(k1,k,cptcoveff); /* Should be the value of the covariate corresponding to k1 combination */
         /* decodtabm(1,1,4) = 1 because h=1  k= (1) 1  1  1 */            /* decodtabm(1,1,4) = 1 because h=1  k= (1) 1  1  1 */
         /* decodtabm(1,2,4) = 1 because h=1  k=  1 (1) 1  1 */            /* decodtabm(1,2,4) = 1 because h=1  k=  1 (1) 1  1 */
         /* decodtabm(13,3,4)= 2 because h=13 k=  1  1 (2) 2 */            /* decodtabm(13,3,4)= 2 because h=13 k=  1  1 (2) 2 */
         vlv= nbcode[Tvaraff[k]][lv]; /* vlv is the value of the covariate lv, 0 or 1 */            vlv= nbcode[Tvaraff[k]][lv]; /* vlv is the value of the covariate lv, 0 or 1 */
         /* For each combination of covariate k1 (V1=1, V3=0), we printed the current covariate k and its value vlv */            /* For each combination of covariate k1 (V1=1, V3=0), we printed the current covariate k and its value vlv */
         printf(" V%d=%d ",Tvaraff[k],vlv);            printf(" V%d=%d ",Tvaraff[k],vlv);
         fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);            fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);
       }          }
       for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */          for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
         printf(" V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);            printf(" V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
         fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);            fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
       }           }       
       printf("\n#\n");          printf("\n#\n");
       fprintf(ficgp,"\n#\n");          fprintf(ficgp,"\n#\n");
       if(invalidvarcomb[k1]){          if(invalidvarcomb[k1]){
         fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);             fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1); 
         continue;            continue;
       }          }
               
       fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" \n",subdirf2(optionfilefiname,"V_"),cpt,k1);          fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" \n",subdirf2(optionfilefiname,"V_"),cpt,k1);
       fprintf(ficgp,"\n#set out \"V_%s_%d-%d.svg\" \n",optionfilefiname,cpt,k1);          fprintf(ficgp,"\n#set out \"V_%s_%d-%d.svg\" \n",optionfilefiname,cpt,k1);
       fprintf(ficgp,"set xlabel \"Age\" \n\          fprintf(ficgp,"set xlabel \"Age\" \n\
 set ylabel \"Probability\" \n             \  set ylabel \"Probability\" \n             \
 set ter svg size 640, 480\n                                             \  set ter svg size 640, 480\n                                             \
 plot [%.f:%.f] \"%s\" every :::%d::%d u 1:2 \"%%lf",ageminpar,fage,subdirf2(fileresu,"VPL_"),k1-1,k1-1);  plot [%.f:%.f] \"%s\" every :::%d::%d u 1:2 \"%%lf",ageminpar,fage,subdirf2(fileresu,"VPL_"),k1-1,k1-1);
               
       for (i=1; i<= nlstate ; i ++) {          for (i=1; i<= nlstate ; i ++) {
         if (i==cpt) fprintf(ficgp," %%lf (%%lf)");            if (i==cpt) fprintf(ficgp," %%lf (%%lf)");
         else        fprintf(ficgp," %%*lf (%%*lf)");            else        fprintf(ficgp," %%*lf (%%*lf)");
       }          }
       fprintf(ficgp,"\" t\"Period (stable) prevalence\" w l lt 0,\"%s\" every :::%d::%d u 1:($2+1.96*$3) \"%%lf",subdirf2(fileresu,"VPL_"),k1-1,k1-1);          fprintf(ficgp,"\" t\"Period (stable) prevalence\" w l lt 0,\"%s\" every :::%d::%d u 1:($2+1.96*$3) \"%%lf",subdirf2(fileresu,"VPL_"),k1-1,k1-1);
       for (i=1; i<= nlstate ; i ++) {          for (i=1; i<= nlstate ; i ++) {
         if (i==cpt) fprintf(ficgp," %%lf (%%lf)");            if (i==cpt) fprintf(ficgp," %%lf (%%lf)");
         else fprintf(ficgp," %%*lf (%%*lf)");            else fprintf(ficgp," %%*lf (%%*lf)");
       }           } 
       fprintf(ficgp,"\" t\"95%% CI\" w l lt 1,\"%s\" every :::%d::%d u 1:($2-1.96*$3) \"%%lf",subdirf2(fileresu,"VPL_"),k1-1,k1-1);           fprintf(ficgp,"\" t\"95%% CI\" w l lt 1,\"%s\" every :::%d::%d u 1:($2-1.96*$3) \"%%lf",subdirf2(fileresu,"VPL_"),k1-1,k1-1); 
       for (i=1; i<= nlstate ; i ++) {          for (i=1; i<= nlstate ; i ++) {
         if (i==cpt) fprintf(ficgp," %%lf (%%lf)");            if (i==cpt) fprintf(ficgp," %%lf (%%lf)");
         else fprintf(ficgp," %%*lf (%%*lf)");            else fprintf(ficgp," %%*lf (%%*lf)");
       }            }  
       fprintf(ficgp,"\" t\"\" w l lt 1,\"%s\" every :::%d::%d u 1:($%d) t\"Observed prevalence\" w l lt 2",subdirf2(fileresu,"P_"),k1-1,k1-1,2+4*(cpt-1));          fprintf(ficgp,"\" t\"\" w l lt 1,\"%s\" every :::%d::%d u 1:($%d) t\"Observed prevalence\" w l lt 2",subdirf2(fileresu,"P_"),k1-1,k1-1,2+4*(cpt-1));
       if(backcast==1){ /* We need to get the corresponding values of the covariates involved in this combination k1 */          if(backcast==1){ /* We need to get the corresponding values of the covariates involved in this combination k1 */
         /* fprintf(ficgp,",\"%s\" every :::%d::%d u 1:($%d) t\"Backward stable prevalence\" w l lt 3",subdirf2(fileresu,"PLB_"),k1-1,k1-1,1+cpt); */            /* fprintf(ficgp,",\"%s\" every :::%d::%d u 1:($%d) t\"Backward stable prevalence\" w l lt 3",subdirf2(fileresu,"PLB_"),k1-1,k1-1,1+cpt); */
         fprintf(ficgp,",\"%s\" u 1:((",subdirf2(fileresu,"PLB_")); /* Age is in 1 */            fprintf(ficgp,",\"%s\" u 1:((",subdirf2(fileresu,"PLB_")); /* Age is in 1 */
         if(cptcoveff ==0){            if(cptcoveff ==0){
           fprintf(ficgp,"$%d)) t 'Backward prevalence in state %d' with line ",  2+(cpt-1),  cpt );              fprintf(ficgp,"$%d)) t 'Backward prevalence in state %d' with line ",        2+(cpt-1),  cpt );
         }else{            }else{
           kl=0;              kl=0;
           for (k=1; k<=cptcoveff; k++){    /* For each combination of covariate  */              for (k=1; k<=cptcoveff; k++){    /* For each combination of covariate  */
             lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate value corresponding to k1 combination and kth covariate */                lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate value corresponding to k1 combination and kth covariate */
             /* decodtabm(1,1,4) = 1 because h=1  k= (1) 1  1  1 */                /* decodtabm(1,1,4) = 1 because h=1  k= (1) 1  1  1 */
             /* decodtabm(1,2,4) = 1 because h=1  k=  1 (1) 1  1 */                /* decodtabm(1,2,4) = 1 because h=1  k=  1 (1) 1  1 */
             /* decodtabm(13,3,4)= 2 because h=13 k=  1  1 (2) 2 */                /* decodtabm(13,3,4)= 2 because h=13 k=  1  1 (2) 2 */
             vlv= nbcode[Tvaraff[k]][lv];                vlv= nbcode[Tvaraff[k]][lv];
             kl++;  
             /* kl=6+(cpt-1)*(nlstate+1)+1+(i-1); /\* 6+(1-1)*(2+1)+1+(1-1)=7, 6+(2-1)(2+1)+1+(1-1)=10 *\/ */  
             /*6+(cpt-1)*(nlstate+1)+1+(i-1)+(nlstate+1)*nlstate; 6+(1-1)*(2+1)+1+(1-1) +(2+1)*2=13 */   
             /*6+1+(i-1)+(nlstate+1)*nlstate; 6+1+(1-1) +(2+1)*2=13 */   
             /* ''  u 6:(($1==1 && $2==0 && $3==2 && $4==0)? $9/(1.-$15) : 1/0):($5==2000? 3:2) t 'p.1' with line lc variable*/  
             if(k==cptcoveff){  
               fprintf(ficgp,"$%d==%d && $%d==%d)? $%d : 1/0) t 'Backward prevalence in state %d' ",kl+1, Tvaraff[k],kl+1+1,nbcode[Tvaraff[k]][lv], \  
                       4+(cpt-1),  cpt );  /* 4 or 6 ?*/  
             }else{  
               fprintf(ficgp,"$%d==%d && $%d==%d && ",kl+1, Tvaraff[k],kl+1+1,nbcode[Tvaraff[k]][lv]);  
               kl++;                kl++;
             }                /* kl=6+(cpt-1)*(nlstate+1)+1+(i-1); /\* 6+(1-1)*(2+1)+1+(1-1)=7, 6+(2-1)(2+1)+1+(1-1)=10 *\/ */
           } /* end covariate */                /*6+(cpt-1)*(nlstate+1)+1+(i-1)+(nlstate+1)*nlstate; 6+(1-1)*(2+1)+1+(1-1) +(2+1)*2=13 */ 
         } /* end if no covariate */                /*6+1+(i-1)+(nlstate+1)*nlstate; 6+1+(1-1) +(2+1)*2=13 */ 
       } /* end if backcast */                /* ''  u 6:(($1==1 && $2==0 && $3==2 && $4==0)? $9/(1.-$15) : 1/0):($5==2000? 3:2) t 'p.1' with line lc variable*/
       fprintf(ficgp,"\nset out \n");                if(k==cptcoveff){
                   fprintf(ficgp,"$%d==%d && $%d==%d)? $%d : 1/0) t 'Backward prevalence in state %d' ",kl+1, Tvaraff[k],kl+1+1,nbcode[Tvaraff[k]][lv], \
                           4+(cpt-1),  cpt );  /* 4 or 6 ?*/
                 }else{
                   fprintf(ficgp,"$%d==%d && $%d==%d && ",kl+1, Tvaraff[k],kl+1+1,nbcode[Tvaraff[k]][lv]);
                   kl++;
                 }
               } /* end covariate */
             } /* end if no covariate */
           } /* end if backcast */
           fprintf(ficgp,"\nset out \n");
         } /* nres */
     } /* k1 */      } /* k1 */
   } /* cpt */    } /* cpt */
   
       
   /*2 eme*/    /*2 eme*/
   for (k1=1; k1<= m ; k1 ++)      for (k1=1; k1<= m ; k1 ++){  
   for(nres=1; nres <= nresult; nres++){ /* For each resultline */      for(nres=1; nres <= nresult; nres++){ /* For each resultline */
     if(TKresult[nres]!= k1)        if(TKresult[nres]!= k1)
       continue;          continue;
     fprintf(ficgp,"\n# 2nd: Total life expectancy with CI: 't' files ");        fprintf(ficgp,"\n# 2nd: Total life expectancy with CI: 't' files ");
     for (k=1; k<=cptcoveff; k++){    /* For each covariate and each value */        for (k=1; k<=cptcoveff; k++){    /* For each covariate and each value */
       lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate number corresponding to k1 combination */  
       /* decodtabm(1,1,4) = 1 because h=1  k= (1) 1  1  1 */  
       /* decodtabm(1,2,4) = 1 because h=1  k=  1 (1) 1  1 */  
       /* decodtabm(13,3,4)= 2 because h=13 k=  1  1 (2) 2 */  
       vlv= nbcode[Tvaraff[k]][lv];  
       fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);  
     }  
     for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */  
       printf(" V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);  
       fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);  
     }  
     fprintf(ficgp,"\n#\n");  
     if(invalidvarcomb[k1]){  
       fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);   
       continue;  
     }  
                           
     fprintf(ficgp,"\nset out \"%s_%d.svg\" \n",subdirf2(optionfilefiname,"E_"),k1);  
     for(vpopbased=0; vpopbased <= popbased; vpopbased++){ /* Done for vpopbased=0 and vpopbased=1 if popbased==1*/  
       if(vpopbased==0)  
         fprintf(ficgp,"set ylabel \"Years\" \nset ter svg size 640, 480\nplot [%.f:%.f] ",ageminpar,fage);  
       else  
         fprintf(ficgp,"\nreplot ");  
       for (i=1; i<= nlstate+1 ; i ++) {  
         k=2*i;  
         fprintf(ficgp,"\"%s\" every :::%d::%d u 1:($2==%d && $4!=0 ?$4 : 1/0) \"%%lf %%lf %%lf",subdirf2(fileresu,"T_"),k1-1,k1-1, vpopbased);  
         for (j=1; j<= nlstate+1 ; j ++) {  
           if (j==i) fprintf(ficgp," %%lf (%%lf)");  
           else fprintf(ficgp," %%*lf (%%*lf)");  
         }     
         if (i== 1) fprintf(ficgp,"\" t\"TLE\" w l lt %d, \\\n",i);  
         else fprintf(ficgp,"\" t\"LE in state (%d)\" w l lt %d, \\\n",i-1,i+1);  
         fprintf(ficgp,"\"%s\" every :::%d::%d u 1:($2==%d && $4!=0 ? $4-$5*2 : 1/0) \"%%lf %%lf %%lf",subdirf2(fileresu,"T_"),k1-1,k1-1,vpopbased);  
         for (j=1; j<= nlstate+1 ; j ++) {  
           if (j==i) fprintf(ficgp," %%lf (%%lf)");  
           else fprintf(ficgp," %%*lf (%%*lf)");  
         }     
         fprintf(ficgp,"\" t\"\" w l lt 0,");  
         fprintf(ficgp,"\"%s\" every :::%d::%d u 1:($2==%d && $4!=0 ? $4+$5*2 : 1/0) \"%%lf %%lf %%lf",subdirf2(fileresu,"T_"),k1-1,k1-1,vpopbased);  
         for (j=1; j<= nlstate+1 ; j ++) {  
           if (j==i) fprintf(ficgp," %%lf (%%lf)");  
           else fprintf(ficgp," %%*lf (%%*lf)");  
         }     
         if (i== (nlstate+1)) fprintf(ficgp,"\" t\"\" w l lt 0");  
         else fprintf(ficgp,"\" t\"\" w l lt 0,\\\n");  
       } /* state */  
     } /* vpopbased */  
     fprintf(ficgp,"\nset out;set out \"%s_%d.svg\"; replot; set out; \n",subdirf2(optionfilefiname,"E_"),k1); /* Buggy gnuplot */  
   } /* k1 end 2 eme*/  
           
           
   /*3eme*/  
   for (k1=1; k1<= m ; k1 ++)   
   for(nres=1; nres <= nresult; nres++){ /* For each resultline */  
     if(TKresult[nres]!= k)  
       continue;  
   
     for (cpt=1; cpt<= nlstate ; cpt ++) {  
       fprintf(ficgp,"\n# 3d: Life expectancy with EXP_ files:  combination=%d state=%d",k1, cpt);  
       for (k=1; k<=cptcoveff; k++){    /* For each covariate dummy combination and each value */  
         lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate number corresponding to k1 combination */          lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate number corresponding to k1 combination */
         /* decodtabm(1,1,4) = 1 because h=1  k= (1) 1  1  1 */          /* decodtabm(1,1,4) = 1 because h=1  k= (1) 1  1  1 */
         /* decodtabm(1,2,4) = 1 because h=1  k=  1 (1) 1  1 */          /* decodtabm(1,2,4) = 1 because h=1  k=  1 (1) 1  1 */
Line 6431  plot [%.f:%.f] \"%s\" every :::%d::%d u Line 6384  plot [%.f:%.f] \"%s\" every :::%d::%d u
       }        }
       /* for(k=1; k <= ncovds; k++){ */        /* for(k=1; k <= ncovds; k++){ */
       for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */        for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
           printf(" V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
         fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);          fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
       }         }
       fprintf(ficgp,"\n#\n");        fprintf(ficgp,"\n#\n");
       if(invalidvarcomb[k1]){        if(invalidvarcomb[k1]){
         fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);           fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1); 
         continue;          continue;
       }        }
                                                   
       /*       k=2+nlstate*(2*cpt-2); */        fprintf(ficgp,"\nset out \"%s_%d.svg\" \n",subdirf2(optionfilefiname,"E_"),k1);
       k=2+(nlstate+1)*(cpt-1);        for(vpopbased=0; vpopbased <= popbased; vpopbased++){ /* Done for vpopbased=0 and vpopbased=1 if popbased==1*/
       fprintf(ficgp,"\nset out \"%s_%d%d.svg\" \n",subdirf2(optionfilefiname,"EXP_"),cpt,k1);          if(vpopbased==0)
       fprintf(ficgp,"set ter svg size 640, 480\n\            fprintf(ficgp,"set ylabel \"Years\" \nset ter svg size 640, 480\nplot [%.f:%.f] ",ageminpar,fage);
           else
             fprintf(ficgp,"\nreplot ");
           for (i=1; i<= nlstate+1 ; i ++) {
             k=2*i;
             fprintf(ficgp,"\"%s\" every :::%d::%d u 1:($2==%d && $4!=0 ?$4 : 1/0) \"%%lf %%lf %%lf",subdirf2(fileresu,"T_"),k1-1,k1-1, vpopbased);
             for (j=1; j<= nlstate+1 ; j ++) {
               if (j==i) fprintf(ficgp," %%lf (%%lf)");
               else fprintf(ficgp," %%*lf (%%*lf)");
             }   
             if (i== 1) fprintf(ficgp,"\" t\"TLE\" w l lt %d, \\\n",i);
             else fprintf(ficgp,"\" t\"LE in state (%d)\" w l lt %d, \\\n",i-1,i+1);
             fprintf(ficgp,"\"%s\" every :::%d::%d u 1:($2==%d && $4!=0 ? $4-$5*2 : 1/0) \"%%lf %%lf %%lf",subdirf2(fileresu,"T_"),k1-1,k1-1,vpopbased);
             for (j=1; j<= nlstate+1 ; j ++) {
               if (j==i) fprintf(ficgp," %%lf (%%lf)");
               else fprintf(ficgp," %%*lf (%%*lf)");
             }   
             fprintf(ficgp,"\" t\"\" w l lt 0,");
             fprintf(ficgp,"\"%s\" every :::%d::%d u 1:($2==%d && $4!=0 ? $4+$5*2 : 1/0) \"%%lf %%lf %%lf",subdirf2(fileresu,"T_"),k1-1,k1-1,vpopbased);
             for (j=1; j<= nlstate+1 ; j ++) {
               if (j==i) fprintf(ficgp," %%lf (%%lf)");
               else fprintf(ficgp," %%*lf (%%*lf)");
             }   
             if (i== (nlstate+1)) fprintf(ficgp,"\" t\"\" w l lt 0");
             else fprintf(ficgp,"\" t\"\" w l lt 0,\\\n");
           } /* state */
         } /* vpopbased */
         fprintf(ficgp,"\nset out;set out \"%s_%d.svg\"; replot; set out; \n",subdirf2(optionfilefiname,"E_"),k1); /* Buggy gnuplot */
       } /* end nres */
     } /* k1 end 2 eme*/
           
           
     /*3eme*/
     for (k1=1; k1<= m ; k1 ++){
       for(nres=1; nres <= nresult; nres++){ /* For each resultline */
         if(TKresult[nres]!= k)
           continue;
   
         for (cpt=1; cpt<= nlstate ; cpt ++) {
           fprintf(ficgp,"\n# 3d: Life expectancy with EXP_ files:  combination=%d state=%d",k1, cpt);
           for (k=1; k<=cptcoveff; k++){    /* For each covariate and each value */
             lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate number corresponding to k1 combination */
             /* decodtabm(1,1,4) = 1 because h=1  k= (1) 1  1  1 */
             /* decodtabm(1,2,4) = 1 because h=1  k=  1 (1) 1  1 */
             /* decodtabm(13,3,4)= 2 because h=13 k=  1  1 (2) 2 */
             vlv= nbcode[Tvaraff[k]][lv];
             fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);
           }
           for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
             fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
           }       
           fprintf(ficgp,"\n#\n");
           if(invalidvarcomb[k1]){
             fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1); 
             continue;
           }
                           
           /*       k=2+nlstate*(2*cpt-2); */
           k=2+(nlstate+1)*(cpt-1);
           fprintf(ficgp,"\nset out \"%s_%d%d.svg\" \n",subdirf2(optionfilefiname,"EXP_"),cpt,k1);
           fprintf(ficgp,"set ter svg size 640, 480\n\
 plot [%.f:%.f] \"%s\" every :::%d::%d u 1:%d t \"e%d1\" w l",ageminpar,fage,subdirf2(fileresu,"E_"),k1-1,k1-1,k,cpt);  plot [%.f:%.f] \"%s\" every :::%d::%d u 1:%d t \"e%d1\" w l",ageminpar,fage,subdirf2(fileresu,"E_"),k1-1,k1-1,k,cpt);
       /*fprintf(ficgp,",\"e%s\" every :::%d::%d u 1:($%d-2*$%d) \"\%%lf ",fileres,k1-1,k1-1,k,k+1);          /*fprintf(ficgp,",\"e%s\" every :::%d::%d u 1:($%d-2*$%d) \"\%%lf ",fileres,k1-1,k1-1,k,k+1);
         for (i=1; i<= nlstate*2 ; i ++) fprintf(ficgp,"\%%lf (\%%lf) ");            for (i=1; i<= nlstate*2 ; i ++) fprintf(ficgp,"\%%lf (\%%lf) ");
         fprintf(ficgp,"\" t \"e%d1\" w l",cpt);            fprintf(ficgp,"\" t \"e%d1\" w l",cpt);
         fprintf(ficgp,",\"e%s\" every :::%d::%d u 1:($%d+2*$%d) \"\%%lf ",fileres,k1-1,k1-1,k,k+1);            fprintf(ficgp,",\"e%s\" every :::%d::%d u 1:($%d+2*$%d) \"\%%lf ",fileres,k1-1,k1-1,k,k+1);
         for (i=1; i<= nlstate*2 ; i ++) fprintf(ficgp,"\%%lf (\%%lf) ");            for (i=1; i<= nlstate*2 ; i ++) fprintf(ficgp,"\%%lf (\%%lf) ");
         fprintf(ficgp,"\" t \"e%d1\" w l",cpt);            fprintf(ficgp,"\" t \"e%d1\" w l",cpt);
                                                                   
       */          */
       for (i=1; i< nlstate ; i ++) {          for (i=1; i< nlstate ; i ++) {
         fprintf(ficgp," ,\"%s\" every :::%d::%d u 1:%d t \"e%d%d\" w l",subdirf2(fileresu,"E_"),k1-1,k1-1,k+i,cpt,i+1);            fprintf(ficgp," ,\"%s\" every :::%d::%d u 1:%d t \"e%d%d\" w l",subdirf2(fileresu,"E_"),k1-1,k1-1,k+i,cpt,i+1);
         /*      fprintf(ficgp," ,\"%s\" every :::%d::%d u 1:%d t \"e%d%d\" w l",subdirf2(fileres,"e"),k1-1,k1-1,k+2*i,cpt,i+1);*/            /*    fprintf(ficgp," ,\"%s\" every :::%d::%d u 1:%d t \"e%d%d\" w l",subdirf2(fileres,"e"),k1-1,k1-1,k+2*i,cpt,i+1);*/
                                                                   
       }           } 
       fprintf(ficgp," ,\"%s\" every :::%d::%d u 1:%d t \"e%d.\" w l",subdirf2(fileresu,"E_"),k1-1,k1-1,k+nlstate,cpt);          fprintf(ficgp," ,\"%s\" every :::%d::%d u 1:%d t \"e%d.\" w l",subdirf2(fileresu,"E_"),k1-1,k1-1,k+nlstate,cpt);
     }        }
   }      } /* end nres */
     } /* end kl 3eme */
       
   /* 4eme */    /* 4eme */
   /* Survival functions (period) from state i in state j by initial state i */    /* Survival functions (period) from state i in state j by initial state i */
   for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */    for (k1=1; k1<=m; k1++){    /* For each covariate and each value */
   for(nres=1; nres <= nresult; nres++){ /* For each resultline */      for(nres=1; nres <= nresult; nres++){ /* For each resultline */
     if(TKresult[nres]!= k1)        if(TKresult[nres]!= k1)
       continue;  
   
     for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life state */  
       fprintf(ficgp,"\n#\n#\n# Survival functions in state j : 'LIJ_' files, cov=%d state=%d",k1, cpt);  
       for (k=1; k<=cptcoveff; k++){    /* For each covariate and each value */  
         lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate number corresponding to k1 combination */  
         /* decodtabm(1,1,4) = 1 because h=1  k= (1) 1  1  1 */  
         /* decodtabm(1,2,4) = 1 because h=1  k=  1 (1) 1  1 */  
         /* decodtabm(13,3,4)= 2 because h=13 k=  1  1 (2) 2 */  
         vlv= nbcode[Tvaraff[k]][lv];  
         fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);  
       }  
       for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */  
         fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);  
       }   
       fprintf(ficgp,"\n#\n");  
       if(invalidvarcomb[k1]){  
         fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);   
         continue;          continue;
       }        for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each life state cpt*/
                                   fprintf(ficgp,"\n#\n#\n# Survival functions in state j : 'LIJ_' files, cov=%d state=%d",k1, cpt);
       fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" \n",subdirf2(optionfilefiname,"LIJ_"),cpt,k1);          for (k=1; k<=cptcoveff; k++){    /* For each covariate and each value */
       fprintf(ficgp,"set xlabel \"Age\" \nset ylabel \"Probability to be alive\" \n\            lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate number corresponding to k1 combination */
 set ter svg size 640, 480\n                                                                                                                                                                                     \            /* decodtabm(1,1,4) = 1 because h=1  k= (1) 1  1  1 */
 unset log y\n                                                                                                                                                                                                                                           \            /* decodtabm(1,2,4) = 1 because h=1  k=  1 (1) 1  1 */
 plot [%.f:%.f]  ", ageminpar, agemaxpar);            /* decodtabm(13,3,4)= 2 because h=13 k=  1  1 (2) 2 */
       k=3;            vlv= nbcode[Tvaraff[k]][lv];
       for (i=1; i<= nlstate ; i ++){            fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);
         if(i==1){  
           fprintf(ficgp,"\"%s\"",subdirf2(fileresu,"PIJ_"));  
         }else{  
           fprintf(ficgp,", '' ");  
         }          }
         l=(nlstate+ndeath)*(i-1)+1;          for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
         fprintf(ficgp," u ($1==%d ? ($3):1/0):($%d/($%d",k1,k+l+(cpt-1),k+l);            fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
         for (j=2; j<= nlstate+ndeath ; j ++)          }       
           fprintf(ficgp,"+$%d",k+l+j-1);          fprintf(ficgp,"\n#\n");
         fprintf(ficgp,")) t \"l(%d,%d)\" w l",i,cpt);          if(invalidvarcomb[k1]){
       } /* nlstate */            fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1); 
       fprintf(ficgp,"\nset out\n");            continue;
     } /* end cpt state*/           }
   } /* end covariate */          
                   fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" \n",subdirf2(optionfilefiname,"LIJ_"),cpt,k1);
           fprintf(ficgp,"set xlabel \"Age\" \nset ylabel \"Probability to be alive\" \n\
   set ter svg size 640, 480\nunset log y\nplot [%.f:%.f]  ", ageminpar, agemaxpar);
           k=3;
           for (i=1; i<= nlstate ; i ++){
             if(i==1){
               fprintf(ficgp,"\"%s\"",subdirf2(fileresu,"PIJ_"));
             }else{
               fprintf(ficgp,", '' ");
             }
             l=(nlstate+ndeath)*(i-1)+1;
             fprintf(ficgp," u ($1==%d ? ($3):1/0):($%d/($%d",k1,k+l+(cpt-1),k+l);
             for (j=2; j<= nlstate+ndeath ; j ++)
               fprintf(ficgp,"+$%d",k+l+j-1);
             fprintf(ficgp,")) t \"l(%d,%d)\" w l",i,cpt);
           } /* nlstate */
           fprintf(ficgp,"\nset out\n");
         } /* end cpt state*/ 
       } /* end nres */
     } /* end covariate k1 */  
   
 /* 5eme */  /* 5eme */
   /* Survival functions (period) from state i in state j by final state j */    /* Survival functions (period) from state i in state j by final state j */
   for (k1=1; k1<= m ; k1 ++) /* For each covariate combination if any */    for (k1=1; k1<= m ; k1++){ /* For each covariate combination if any */
   for(nres=1; nres <= nresult; nres++){ /* For each resultline */      for(nres=1; nres <= nresult; nres++){ /* For each resultline */
     if(TKresult[nres]!= k1)        if(TKresult[nres]!= k1)
       continue;  
     for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each inital state  */  
                           
       fprintf(ficgp,"\n#\n#\n# Survival functions in state j and all livestates from state i by final state j: 'lij' files, cov=%d state=%d",k1, cpt);  
       for (k=1; k<=cptcoveff; k++){    /* For each covariate and each value */  
         lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate number corresponding to k1 combination */  
         /* decodtabm(1,1,4) = 1 because h=1  k= (1) 1  1  1 */  
         /* decodtabm(1,2,4) = 1 because h=1  k=  1 (1) 1  1 */  
         /* decodtabm(13,3,4)= 2 because h=13 k=  1  1 (2) 2 */  
         vlv= nbcode[Tvaraff[k]][lv];  
         fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);  
       }  
       for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */  
         fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);  
       }   
       fprintf(ficgp,"\n#\n");  
       if(invalidvarcomb[k1]){  
         fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1);   
         continue;          continue;
       }        for (cpt=1; cpt<=nlstate ; cpt ++) { /* For each inital state  */
           fprintf(ficgp,"\n#\n#\n# Survival functions in state j and all livestates from state i by final state j: 'lij' files, cov=%d state=%d",k1, cpt);
           for (k=1; k<=cptcoveff; k++){    /* For each covariate and each value */
             lv= decodtabm(k1,k,cptcoveff); /* Should be the covariate number corresponding to k1 combination */
             /* decodtabm(1,1,4) = 1 because h=1  k= (1) 1  1  1 */
             /* decodtabm(1,2,4) = 1 because h=1  k=  1 (1) 1  1 */
             /* decodtabm(13,3,4)= 2 because h=13 k=  1  1 (2) 2 */
             vlv= nbcode[Tvaraff[k]][lv];
             fprintf(ficgp," V%d=%d ",Tvaraff[k],vlv);
           }
           for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
             fprintf(ficgp," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
           }       
           fprintf(ficgp,"\n#\n");
           if(invalidvarcomb[k1]){
             fprintf(ficgp,"#Combination (%d) ignored because no cases \n",k1); 
             continue;
           }
               
       fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" \n",subdirf2(optionfilefiname,"LIJT_"),cpt,k1);          fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" \n",subdirf2(optionfilefiname,"LIJT_"),cpt,k1);
       fprintf(ficgp,"set xlabel \"Age\" \nset ylabel \"Probability to be alive\" \n\          fprintf(ficgp,"set xlabel \"Age\" \nset ylabel \"Probability to be alive\" \n\
 set ter svg size 640, 480\n                                             \  set ter svg size 640, 480\nunset log y\nplot [%.f:%.f]  ", ageminpar, agemaxpar);
 unset log y\n                                                           \          k=3;
 plot [%.f:%.f]  ", ageminpar, agemaxpar);          for (j=1; j<= nlstate ; j ++){ /* Lived in state j */
       k=3;            if(j==1)
       for (j=1; j<= nlstate ; j ++){ /* Lived in state j */              fprintf(ficgp,"\"%s\"",subdirf2(fileresu,"PIJ_"));
         if(j==1)            else
           fprintf(ficgp,"\"%s\"",subdirf2(fileresu,"PIJ_"));              fprintf(ficgp,", '' ");
         else            l=(nlstate+ndeath)*(cpt-1) +j;
           fprintf(ficgp,", '' ");            fprintf(ficgp," u (($1==%d && (floor($2)%%5 == 0)) ? ($3):1/0):($%d",k1,k+l);
         l=(nlstate+ndeath)*(cpt-1) +j;            /* for (i=2; i<= nlstate+ndeath ; i ++) */
         fprintf(ficgp," u (($1==%d && (floor($2)%%5 == 0)) ? ($3):1/0):($%d",k1,k+l);            /*   fprintf(ficgp,"+$%d",k+l+i-1); */
         /* for (i=2; i<= nlstate+ndeath ; i ++) */            fprintf(ficgp,") t \"l(%d,%d)\" w l",cpt,j);
         /*   fprintf(ficgp,"+$%d",k+l+i-1); */          } /* nlstate */
         fprintf(ficgp,") t \"l(%d,%d)\" w l",cpt,j);          fprintf(ficgp,", '' ");
       } /* nlstate */          fprintf(ficgp," u (($1==%d && (floor($2)%%5 == 0)) ? ($3):1/0):(",k1);
       fprintf(ficgp,", '' ");          for (j=1; j<= nlstate ; j ++){ /* Lived in state j */
       fprintf(ficgp," u (($1==%d && (floor($2)%%5 == 0)) ? ($3):1/0):(",k1);            l=(nlstate+ndeath)*(cpt-1) +j;
       for (j=1; j<= nlstate ; j ++){ /* Lived in state j */            if(j < nlstate)
         l=(nlstate+ndeath)*(cpt-1) +j;              fprintf(ficgp,"$%d +",k+l);
         if(j < nlstate)            else
           fprintf(ficgp,"$%d +",k+l);              fprintf(ficgp,"$%d) t\"l(%d,.)\" w l",k+l,cpt);
         else          }
           fprintf(ficgp,"$%d) t\"l(%d,.)\" w l",k+l,cpt);          fprintf(ficgp,"\nset out\n");
       }        } /* end cpt state*/ 
       fprintf(ficgp,"\nset out\n");      } /* end covariate */  
     } /* end cpt state*/     } /* end nres */
   } /* end covariate */    
       
 /* 6eme */  /* 6eme */
   /* CV preval stable (period) for each covariate */    /* CV preval stable (period) for each covariate */
Line 6593  plot [%.f:%.f]  ", ageminpar, agemaxpar) Line 6604  plot [%.f:%.f]  ", ageminpar, agemaxpar)
               
       fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" \n",subdirf2(optionfilefiname,"P_"),cpt,k1);        fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" \n",subdirf2(optionfilefiname,"P_"),cpt,k1);
       fprintf(ficgp,"set xlabel \"Age\" \nset ylabel \"Probability\" \n\        fprintf(ficgp,"set xlabel \"Age\" \nset ylabel \"Probability\" \n\
 set ter svg size 640, 480\n                                             \  set ter svg size 640, 480\nunset log y\nplot [%.f:%.f]  ", ageminpar, agemaxpar);
 unset log y\n                                                           \  
 plot [%.f:%.f]  ", ageminpar, agemaxpar);  
       k=3; /* Offset */        k=3; /* Offset */
       for (i=1; i<= nlstate ; i ++){        for (i=1; i<= nlstate ; i ++){
         if(i==1)          if(i==1)
Line 6641  plot [%.f:%.f]  ", ageminpar, agemaxpar) Line 6650  plot [%.f:%.f]  ", ageminpar, agemaxpar)
                   
         fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" \n",subdirf2(optionfilefiname,"PB_"),cpt,k1);          fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" \n",subdirf2(optionfilefiname,"PB_"),cpt,k1);
         fprintf(ficgp,"set xlabel \"Age\" \nset ylabel \"Probability\" \n\          fprintf(ficgp,"set xlabel \"Age\" \nset ylabel \"Probability\" \n\
 set ter svg size 640, 480\n                                             \  set ter svg size 640, 480\nunset log y\nplot [%.f:%.f]  ", ageminpar, agemaxpar);
 unset log y\n                                                           \  
 plot [%.f:%.f]  ", ageminpar, agemaxpar);  
         k=3; /* Offset */          k=3; /* Offset */
         for (i=1; i<= nlstate ; i ++){          for (i=1; i<= nlstate ; i ++){
           if(i==1)            if(i==1)
Line 6695  plot [%.f:%.f]  ", ageminpar, agemaxpar) Line 6702  plot [%.f:%.f]  ", ageminpar, agemaxpar)
         fprintf(ficgp,"# hpijx=probability over h years, hp.jx is weighted by observed prev\n ");          fprintf(ficgp,"# hpijx=probability over h years, hp.jx is weighted by observed prev\n ");
         fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" \n",subdirf2(optionfilefiname,"PROJ_"),cpt,k1);          fprintf(ficgp,"\nset out \"%s_%d-%d.svg\" \n",subdirf2(optionfilefiname,"PROJ_"),cpt,k1);
         fprintf(ficgp,"set xlabel \"Age\" \nset ylabel \"Prevalence\" \n\          fprintf(ficgp,"set xlabel \"Age\" \nset ylabel \"Prevalence\" \n\
 set ter svg size 640, 480\n                                             \  set ter svg size 640, 480\nunset log y\nplot [%.f:%.f]  ", ageminpar, agemaxpar);
 unset log y\n                                                           \  
 plot [%.f:%.f]  ", ageminpar, agemaxpar);  
         for (i=1; i<= nlstate+1 ; i ++){  /* nlstate +1 p11 p21 p.1 */          for (i=1; i<= nlstate+1 ; i ++){  /* nlstate +1 p11 p21 p.1 */
           /*#  V1  = 1  V2 =  0 yearproj age p11 p21 p.1 p12 p22 p.2 p13 p23 p.3*/            /*#  V1  = 1  V2 =  0 yearproj age p11 p21 p.1 p12 p22 p.2 p13 p23 p.3*/
           /*#   1    2   3    4    5      6  7   8   9   10   11 12  13   14  15 */               /*#   1    2   3    4    5      6  7   8   9   10   11 12  13   14  15 */   
Line 6764  plot [%.f:%.f]  ", ageminpar, agemaxpar) Line 6769  plot [%.f:%.f]  ", ageminpar, agemaxpar)
   } /* End if prevfcast */    } /* End if prevfcast */
       
       
   /* proba elementaires */    /* 9eme writing MLE parameters */
   fprintf(ficgp,"\n##############\n#MLE estimated parameters\n#############\n");    fprintf(ficgp,"\n##############\n#9eme MLE estimated parameters\n#############\n");
   for(i=1,jk=1; i <=nlstate; i++){    for(i=1,jk=1; i <=nlstate; i++){
     fprintf(ficgp,"# initial state %d\n",i);      fprintf(ficgp,"# initial state %d\n",i);
     for(k=1; k <=(nlstate+ndeath); k++){      for(k=1; k <=(nlstate+ndeath); k++){
Line 6782  plot [%.f:%.f]  ", ageminpar, agemaxpar) Line 6787  plot [%.f:%.f]  ", ageminpar, agemaxpar)
   fprintf(ficgp,"##############\n#\n");    fprintf(ficgp,"##############\n#\n");
       
   /*goto avoid;*/    /*goto avoid;*/
   fprintf(ficgp,"\n##############\n#Graphics of probabilities or incidences\n#############\n");    /* 10eme Graphics of probabilities or incidences using written MLE parameters */
     fprintf(ficgp,"\n##############\n#10eme Graphics of probabilities or incidences\n#############\n");
   fprintf(ficgp,"# logi(p12/p11)=a12+b12*age+c12age*age+d12*V1+e12*V1*age\n");    fprintf(ficgp,"# logi(p12/p11)=a12+b12*age+c12age*age+d12*V1+e12*V1*age\n");
   fprintf(ficgp,"# logi(p12/p11)=p1 +p2*age +p3*age*age+ p4*V1+ p5*V1*age\n");    fprintf(ficgp,"# logi(p12/p11)=p1 +p2*age +p3*age*age+ p4*V1+ p5*V1*age\n");
   fprintf(ficgp,"# logi(p13/p11)=a13+b13*age+c13age*age+d13*V1+e13*V1*age\n");    fprintf(ficgp,"# logi(p13/p11)=a13+b13*age+c13age*age+d13*V1+e13*V1*age\n");
Line 6797  plot [%.f:%.f]  ", ageminpar, agemaxpar) Line 6803  plot [%.f:%.f]  ", ageminpar, agemaxpar)
   fprintf(ficgp,"#       +exp(a14+b14*age+c14age*age+d14*V1+e14*V1*age)+...)\n");    fprintf(ficgp,"#       +exp(a14+b14*age+c14age*age+d14*V1+e14*V1*age)+...)\n");
   fprintf(ficgp,"#\n");    fprintf(ficgp,"#\n");
   for(ng=1; ng<=3;ng++){ /* Number of graphics: first is logit, 2nd is probabilities, third is incidences per year*/    for(ng=1; ng<=3;ng++){ /* Number of graphics: first is logit, 2nd is probabilities, third is incidences per year*/
     fprintf(ficgp,"#Number of graphics: first is logit, 2nd is probabilities, third is incidences per year \n");      fprintf(ficgp,"#Number of graphics: first is logit, 2nd is probabilities, third is incidences per year\n");
     fprintf(ficgp,"#model=%s \n",model);      fprintf(ficgp,"#model=%s \n",model);
     fprintf(ficgp,"# ng=%d\n",ng);      fprintf(ficgp,"# Type of graphic ng=%d\n",ng);
     fprintf(ficgp,"#   jk=1 to 2^%d=%d\n",cptcoveff,m);/* to be checked */      fprintf(ficgp,"#   jk=1 to 2^%d=%d\n",cptcoveff,m);/* to be checked */
     for(jk=1; jk <=m; jk++)  /* For each combination of covariate */      for(jk=1; jk <=m; jk++)  /* For each combination of covariate */
     for(nres=1; nres <= nresult; nres++){ /* For each resultline */      for(nres=1; nres <= nresult; nres++){ /* For each resultline */
Line 6856  plot [%.f:%.f]  ", ageminpar, agemaxpar) Line 6862  plot [%.f:%.f]  ", ageminpar, agemaxpar)
               /* printf("Tage[%d]=%d, j=%d\n", ij, Tage[ij], j); */                /* printf("Tage[%d]=%d, j=%d\n", ij, Tage[ij], j); */
               if(j==Tage[ij]) { /* Product by age */                if(j==Tage[ij]) { /* Product by age */
                 if(ij <=cptcovage) { /* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1, 2 V5 and V1 */                  if(ij <=cptcovage) { /* V5+V4+V3+V4*V3+V5*age+V2+V1*V2+V1*age+V1, 2 V5 and V1 */
                   if(Dummy[j]==0){                    if(DummyV[j]==0){
                     fprintf(ficgp,"+p%d*%d*x",i+j+2+nagesqr-1,Tinvresult[nres][Tvar[j]]);;                      fprintf(ficgp,"+p%d*%d*x",i+j+2+nagesqr-1,Tinvresult[nres][Tvar[j]]);;
                   }else{ /* quantitative */                    }else{ /* quantitative */
                     fprintf(ficgp,"+p%d*%f*x",i+j+2+nagesqr-1,Tqinvresult[nres][Tvar[j]]); /* Tqinvresult in decoderesult */                      fprintf(ficgp,"+p%d*%f*x",i+j+2+nagesqr-1,Tqinvresult[nres][Tvar[j]]); /* Tqinvresult in decoderesult */
Line 6867  plot [%.f:%.f]  ", ageminpar, agemaxpar) Line 6873  plot [%.f:%.f]  ", ageminpar, agemaxpar)
               }else if(j==Tprod[ijp]) { /* */                 }else if(j==Tprod[ijp]) { /* */ 
                 /* printf("Tprod[%d]=%d, j=%d\n", ij, Tprod[ijp], j); */                  /* printf("Tprod[%d]=%d, j=%d\n", ij, Tprod[ijp], j); */
                 if(ijp <=cptcovprod) { /* Product */                  if(ijp <=cptcovprod) { /* Product */
                   if(Dummy[Tvard[ijp][1]]==0){/* Vn is dummy */                    if(DummyV[Tvard[ijp][1]]==0){/* Vn is dummy */
                     if(Dummy[Tvard[ijp][2]]==0){/* Vn and Vm are dummy */                      if(DummyV[Tvard[ijp][2]]==0){/* Vn and Vm are dummy */
                       /* fprintf(ficgp,"+p%d*%d*%d",i+j+2+nagesqr-1,nbcode[Tvard[ijp][1]][codtabm(jk,j)],nbcode[Tvard[ijp][2]][codtabm(jk,j)]); */                        /* fprintf(ficgp,"+p%d*%d*%d",i+j+2+nagesqr-1,nbcode[Tvard[ijp][1]][codtabm(jk,j)],nbcode[Tvard[ijp][2]][codtabm(jk,j)]); */
                       fprintf(ficgp,"+p%d*%d*%d",i+j+2+nagesqr-1,Tinvresult[nres][Tvard[ijp][1]],Tinvresult[nres][Tvard[ijp][2]]);                        fprintf(ficgp,"+p%d*%d*%d",i+j+2+nagesqr-1,Tinvresult[nres][Tvard[ijp][1]],Tinvresult[nres][Tvard[ijp][2]]);
                     }else{ /* Vn is dummy and Vm is quanti */                      }else{ /* Vn is dummy and Vm is quanti */
Line 6876  plot [%.f:%.f]  ", ageminpar, agemaxpar) Line 6882  plot [%.f:%.f]  ", ageminpar, agemaxpar)
                       fprintf(ficgp,"+p%d*%d*%f",i+j+2+nagesqr-1,Tinvresult[nres][Tvard[ijp][1]],Tqinvresult[nres][Tvard[ijp][2]]);                        fprintf(ficgp,"+p%d*%d*%f",i+j+2+nagesqr-1,Tinvresult[nres][Tvard[ijp][1]],Tqinvresult[nres][Tvard[ijp][2]]);
                     }                      }
                   }else{ /* Vn*Vm Vn is quanti */                    }else{ /* Vn*Vm Vn is quanti */
                     if(Dummy[Tvard[ijp][2]]==0){                      if(DummyV[Tvard[ijp][2]]==0){
                       fprintf(ficgp,"+p%d*%d*%f",i+j+2+nagesqr-1,Tinvresult[nres][Tvard[ijp][2]],Tqinvresult[nres][Tvard[ijp][1]]);                        fprintf(ficgp,"+p%d*%d*%f",i+j+2+nagesqr-1,Tinvresult[nres][Tvard[ijp][2]],Tqinvresult[nres][Tvard[ijp][1]]);
                     }else{ /* Both quanti */                      }else{ /* Both quanti */
                       fprintf(ficgp,"+p%d*%f*%f",i+j+2+nagesqr-1,Tqinvresult[nres][Tvard[ijp][1]],Tqinvresult[nres][Tvard[ijp][2]]);                        fprintf(ficgp,"+p%d*%f*%f",i+j+2+nagesqr-1,Tqinvresult[nres][Tvard[ijp][1]],Tqinvresult[nres][Tvard[ijp][2]]);
                     }                      }
                   }                    }
                     ijp++;
                 }                  }
               } else{  /* simple covariate */                } else{  /* simple covariate */
                 /* fprintf(ficgp,"+p%d*%d",i+j+2+nagesqr-1,nbcode[Tvar[j]][codtabm(jk,j)]); /\* Valgrind bug nbcode *\/ */                  /* fprintf(ficgp,"+p%d*%d",i+j+2+nagesqr-1,nbcode[Tvar[j]][codtabm(jk,j)]); /\* Valgrind bug nbcode *\/ */
Line 7175  plot [%.f:%.f]  ", ageminpar, agemaxpar) Line 7182  plot [%.f:%.f]  ", ageminpar, agemaxpar)
       fprintf(ficresf," V%d (=) %d",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);        fprintf(ficresf," V%d (=) %d",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
     }      }
     for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */      for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
       printf(" V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);        fprintf(ficresf," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
       fprintf(ficlog," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);  
     }      }
     fprintf(ficresf," yearproj age");      fprintf(ficresf," yearproj age");
     for(j=1; j<=nlstate+ndeath;j++){       for(j=1; j<=nlstate+ndeath;j++){ 
Line 7822  int readdata(char datafile[], int firsto Line 7828  int readdata(char datafile[], int firsto
     /* Loops on waves */      /* Loops on waves */
     for (j=maxwav;j>=1;j--){      for (j=maxwav;j>=1;j--){
       for (iv=nqtv;iv>=1;iv--){  /* Loop  on time varying quantitative variables */        for (iv=nqtv;iv>=1;iv--){  /* Loop  on time varying quantitative variables */
                                 cutv(stra, strb, line, ' ');           cutv(stra, strb, line, ' '); 
                                 if(strb[0]=='.') { /* Missing value */          if(strb[0]=='.') { /* Missing value */
                                         lval=-1;            lval=-1;
                                         cotqvar[j][iv][i]=-1; /* 0.0/0.0 */            cotqvar[j][iv][i]=-1; /* 0.0/0.0 */
                                         cotvar[j][ntv+iv][i]=-1; /* For performance reasons */            cotvar[j][ntv+iv][i]=-1; /* For performance reasons */
                                         if(isalpha(strb[1])) { /* .m or .d Really Missing value */            if(isalpha(strb[1])) { /* .m or .d Really Missing value */
                                                 printf("Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be the %d th quantitative value out of %d measured at wave %d. If missing, you should remove this individual or impute a value.  Exiting.\n", strb, linei,i,line,iv, nqtv, j);              printf("Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be the %d th quantitative value out of %d measured at wave %d. If missing, you should remove this individual or impute a value.  Exiting.\n", strb, linei,i,line,iv, nqtv, j);
                                                 fprintf(ficlog,"Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be the %d th quantitative value out of %d measured at wave %d. If missing, you should remove this individual or impute a value.  Exiting.\n", strb, linei,i,line,iv, nqtv, j);fflush(ficlog);              fprintf(ficlog,"Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be the %d th quantitative value out of %d measured at wave %d. If missing, you should remove this individual or impute a value.  Exiting.\n", strb, linei,i,line,iv, nqtv, j);fflush(ficlog);
                                                 return 1;              return 1;
                                         }            }
                                 }else{          }else{
                                         errno=0;            errno=0;
                                         /* what_kind_of_number(strb); */            /* what_kind_of_number(strb); */
                                         dval=strtod(strb,&endptr);             dval=strtod(strb,&endptr); 
                                         /* if( strb[0]=='\0' || (*endptr != '\0')){ */            /* if( strb[0]=='\0' || (*endptr != '\0')){ */
                                         /* if(strb != endptr && *endptr == '\0') */            /* if(strb != endptr && *endptr == '\0') */
                                         /*    dval=dlval; */            /*    dval=dlval; */
                                         /* if (errno == ERANGE && (lval == LONG_MAX || lval == LONG_MIN)) */            /* if (errno == ERANGE && (lval == LONG_MAX || lval == LONG_MIN)) */
                                         if( strb[0]=='\0' || (*endptr != '\0')){            if( strb[0]=='\0' || (*endptr != '\0')){
                                                 printf("Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be the %d th quantitative value out of %d measured at wave %d. Setting maxwav=%d might be wrong.  Exiting.\n", strb, linei,i,line,iv, nqtv, j,maxwav);              printf("Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be the %d th quantitative value out of %d measured at wave %d. Setting maxwav=%d might be wrong.  Exiting.\n", strb, linei,i,line,iv, nqtv, j,maxwav);
                                                 fprintf(ficlog,"Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be the %d th quantitative value out of %d measured at wave %d. Setting maxwav=%d might be wrong.  Exiting.\n", strb, linei,i,line, iv, nqtv, j,maxwav);fflush(ficlog);              fprintf(ficlog,"Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be the %d th quantitative value out of %d measured at wave %d. Setting maxwav=%d might be wrong.  Exiting.\n", strb, linei,i,line, iv, nqtv, j,maxwav);fflush(ficlog);
                                                 return 1;              return 1;
                                         }            }
                                         cotqvar[j][iv][i]=dval;             cotqvar[j][iv][i]=dval; 
                                         cotvar[j][ntv+iv][i]=dval;             cotvar[j][ntv+iv][i]=dval; 
                                 }          }
                                 strcpy(line,stra);          strcpy(line,stra);
       }/* end loop ntqv */        }/* end loop ntqv */
               
       for (iv=ntv;iv>=1;iv--){  /* Loop  on time varying dummies */        for (iv=ntv;iv>=1;iv--){  /* Loop  on time varying dummies */
                                 cutv(stra, strb, line, ' ');           cutv(stra, strb, line, ' '); 
                                 if(strb[0]=='.') { /* Missing value */          if(strb[0]=='.') { /* Missing value */
                                         lval=-1;            lval=-1;
                                 }else{          }else{
                                         errno=0;            errno=0;
                                         lval=strtol(strb,&endptr,10);             lval=strtol(strb,&endptr,10); 
                                         /*      if (errno == ERANGE && (lval == LONG_MAX || lval == LONG_MIN))*/            /*    if (errno == ERANGE && (lval == LONG_MAX || lval == LONG_MIN))*/
                                         if( strb[0]=='\0' || (*endptr != '\0')){            if( strb[0]=='\0' || (*endptr != '\0')){
                                                 printf("Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be the %d th dummy covariate out of %d measured at wave %d. Setting maxwav=%d might be wrong.  Exiting.\n", strb, linei,i,line,iv, ntv, j,maxwav);              printf("Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be the %d th dummy covariate out of %d measured at wave %d. Setting maxwav=%d might be wrong.  Exiting.\n", strb, linei,i,line,iv, ntv, j,maxwav);
                                                 fprintf(ficlog,"Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be the %d dummy covariate out of %d measured wave %d. Setting maxwav=%d might be wrong.  Exiting.\n", strb, linei,i,line,iv, ntv,j,maxwav);fflush(ficlog);              fprintf(ficlog,"Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be the %d dummy covariate out of %d measured wave %d. Setting maxwav=%d might be wrong.  Exiting.\n", strb, linei,i,line,iv, ntv,j,maxwav);fflush(ficlog);
                                                 return 1;              return 1;
                                         }            }
                                 }          }
                                 if(lval <-1 || lval >1){          if(lval <-1 || lval >1){
                                         printf("Error reading data around '%ld' at line number %d for individual %d, '%s'\n \            printf("Error reading data around '%ld' at line number %d for individual %d, '%s'\n \
  Should be a value of %d(nth) covariate (0 should be the value for the reference and 1\n \   Should be a value of %d(nth) covariate (0 should be the value for the reference and 1\n \
  for the alternative. IMaCh does not build design variables automatically, do it yourself.\n \   for the alternative. IMaCh does not build design variables automatically, do it yourself.\n \
  For example, for multinomial values like 1, 2 and 3,\n                                                                 \   For example, for multinomial values like 1, 2 and 3,\n                 \
  build V1=0 V2=0 for the reference value (1),\n                                                                                                 \   build V1=0 V2=0 for the reference value (1),\n                         \
         V1=1 V2=0 for (2) \n                                                                                                                                                                            \          V1=1 V2=0 for (2) \n                                            \
  and V1=0 V2=1 for (3). V1=1 V2=1 should not exist and the corresponding\n \   and V1=0 V2=1 for (3). V1=1 V2=1 should not exist and the corresponding\n \
  output of IMaCh is often meaningless.\n                                                                                                                                \   output of IMaCh is often meaningless.\n                                \
  Exiting.\n",lval,linei, i,line,j);   Exiting.\n",lval,linei, i,line,j);
                                         fprintf(ficlog,"Error reading data around '%ld' at line number %d for individual %d, '%s'\n \            fprintf(ficlog,"Error reading data around '%ld' at line number %d for individual %d, '%s'\n \
  Should be a value of %d(nth) covariate (0 should be the value for the reference and 1\n \   Should be a value of %d(nth) covariate (0 should be the value for the reference and 1\n \
  for the alternative. IMaCh does not build design variables automatically, do it yourself.\n \   for the alternative. IMaCh does not build design variables automatically, do it yourself.\n \
  For example, for multinomial values like 1, 2 and 3,\n                                                                 \   For example, for multinomial values like 1, 2 and 3,\n                 \
  build V1=0 V2=0 for the reference value (1),\n                                                                                                 \   build V1=0 V2=0 for the reference value (1),\n                         \
         V1=1 V2=0 for (2) \n                                                                                                                                                                            \          V1=1 V2=0 for (2) \n                                            \
  and V1=0 V2=1 for (3). V1=1 V2=1 should not exist and the corresponding\n \   and V1=0 V2=1 for (3). V1=1 V2=1 should not exist and the corresponding\n \
  output of IMaCh is often meaningless.\n                                                                                                                                \   output of IMaCh is often meaningless.\n                                \
  Exiting.\n",lval,linei, i,line,j);fflush(ficlog);   Exiting.\n",lval,linei, i,line,j);fflush(ficlog);
                                         return 1;            return 1;
                                 }          }
                                 cotvar[j][iv][i]=(double)(lval);          cotvar[j][iv][i]=(double)(lval);
                                 strcpy(line,stra);          strcpy(line,stra);
       }/* end loop ntv */        }/* end loop ntv */
               
       /* Statuses  at wave */        /* Statuses  at wave */
       cutv(stra, strb, line, ' ');         cutv(stra, strb, line, ' '); 
       if(strb[0]=='.') { /* Missing value */        if(strb[0]=='.') { /* Missing value */
                                 lval=-1;          lval=-1;
       }else{        }else{
                                 errno=0;          errno=0;
                                 lval=strtol(strb,&endptr,10);           lval=strtol(strb,&endptr,10); 
                                 /*      if (errno == ERANGE && (lval == LONG_MAX || lval == LONG_MIN))*/          /*      if (errno == ERANGE && (lval == LONG_MAX || lval == LONG_MIN))*/
                                 if( strb[0]=='\0' || (*endptr != '\0')){          if( strb[0]=='\0' || (*endptr != '\0')){
                                         printf("Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be a status of wave %d. Setting maxwav=%d might be wrong.  Exiting.\n", strb, linei,i,line,j,maxwav);            printf("Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be a status of wave %d. Setting maxwav=%d might be wrong.  Exiting.\n", strb, linei,i,line,j,maxwav);
                                         fprintf(ficlog,"Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be a status of wave %d. Setting maxwav=%d might be wrong.  Exiting.\n", strb, linei,i,line,j,maxwav);fflush(ficlog);            fprintf(ficlog,"Error reading data around '%s' at line number %d for individual %d, '%s'\nShould be a status of wave %d. Setting maxwav=%d might be wrong.  Exiting.\n", strb, linei,i,line,j,maxwav);fflush(ficlog);
                                         return 1;            return 1;
                                 }          }
       }        }
               
       s[j][i]=lval;        s[j][i]=lval;
Line 8207  int decodemodel( char model[], int lasto Line 8213  int decodemodel( char model[], int lasto
         * - Tvard[k]  p Tvard[1][1]@4 {7, 8, 5, 6} for V7*V8 and V5*V6 .          * - Tvard[k]  p Tvard[1][1]@4 {7, 8, 5, 6} for V7*V8 and V5*V6 .
         */          */
 {  {
   int i, j, k, ks;    int i, j, k, ks, v;
   int  j1, k1, k2, k3, k4;    int  j1, k1, k2, k3, k4;
   char modelsav[80];    char modelsav[80];
   char stra[80], strb[80], strc[80], strd[80],stre[80];    char stra[80], strb[80], strc[80], strd[80],stre[80];
Line 8414  Typevar: 0 for simple covariate (dummy, Line 8420  Typevar: 0 for simple covariate (dummy,
 Fixed[k] 0=fixed (product or simple), 1 varying, 2 fixed with age product, 3 varying with age product \n\  Fixed[k] 0=fixed (product or simple), 1 varying, 2 fixed with age product, 3 varying with age product \n\
 Dummy[k] 0=dummy (0 1), 1 quantitative (single or product without age), 2 dummy with age product, 3 quant with age product\n",model);  Dummy[k] 0=dummy (0 1), 1 quantitative (single or product without age), 2 dummy with age product, 3 quant with age product\n",model);
   
     for(v=1; v <=ncovcol;v++){
       DummyV[v]=0;
       FixedV[v]=0;
     }
     for(v=ncovcol+1; v <=ncovcol+nqv;v++){
       DummyV[v]=1;
       FixedV[v]=0;
     }
     for(v=ncovcol+nqv+1; v <=ncovcol+nqv+ntv;v++){
       DummyV[v]=0;
       FixedV[v]=1;
     }
     for(v=ncovcol+nqv+ntv+1; v <=ncovcol+nqv+ntv+nqtv;v++){
       DummyV[v]=1;
       FixedV[v]=1;
     }
     for(v=1; v <=ncovcol+nqv+ntv+nqtv;v++){
       printf("Decodemodel: V%d, Dummy(V%d)=%d, FixedV(V%d)=%d\n",v,v,DummyV[v],v,FixedV[v]);
       fprintf(ficlog,"Decodemodel: V%d, Dummy(V%d)=%d, FixedV(V%d)=%d\n",v,v,DummyV[v],v,FixedV[v]);
     }
   for(k=1, ncovf=0, nsd=0, nsq=0, ncovv=0, ncova=0, ncoveff=0, nqfveff=0, ntveff=0, nqtveff=0;k<=cptcovt; k++){ /* or cptocvt */    for(k=1, ncovf=0, nsd=0, nsq=0, ncovv=0, ncova=0, ncoveff=0, nqfveff=0, ntveff=0, nqtveff=0;k<=cptcovt; k++){ /* or cptocvt */
     if (Tvar[k] <=ncovcol && Typevar[k]==0 ){ /* Simple fixed dummy (<=ncovcol) covariates */      if (Tvar[k] <=ncovcol && Typevar[k]==0 ){ /* Simple fixed dummy (<=ncovcol) covariates */
       Fixed[k]= 0;        Fixed[k]= 0;
Line 9028  int prevalence_limit(double *p, double * Line 9054  int prevalence_limit(double *p, double *
   i1=pow(2,cptcoveff); /* Number of combination of dummy covariates */    i1=pow(2,cptcoveff); /* Number of combination of dummy covariates */
   if (cptcovn < 1){i1=1;}    if (cptcovn < 1){i1=1;}
   
   for(nres=1; nres <= nresult; nres++) /* For each resultline */    for(k=1; k<=i1;k++){ /* For each combination k of dummy covariates in the model */
   for(k=1; k<=i1;k++){      for(nres=1; nres <= nresult; nres++){ /* For each resultline */
     if(TKresult[nres]!= k)        if(TKresult[nres]!= k)
       continue;          continue;
   
   /* for(cptcov=1,k=0;cptcov<=i1;cptcov++){ */        /* for(cptcov=1,k=0;cptcov<=i1;cptcov++){ */
     /* for(cptcov=1,k=0;cptcov<=1;cptcov++){ */        /* for(cptcov=1,k=0;cptcov<=1;cptcov++){ */
     //for(cptcod=1;cptcod<=ncodemax[cptcov];cptcod++){        //for(cptcod=1;cptcod<=ncodemax[cptcov];cptcod++){
     /* k=k+1; */        /* k=k+1; */
     /* to clean */        /* to clean */
     //printf("cptcov=%d cptcod=%d codtab=%d\n",cptcov, cptcod,codtabm(cptcod,cptcov));        //printf("cptcov=%d cptcod=%d codtab=%d\n",cptcov, cptcod,codtabm(cptcod,cptcov));
     fprintf(ficrespl,"#******");        fprintf(ficrespl,"#******");
     printf("#******");        printf("#******");
     fprintf(ficlog,"#******");        fprintf(ficlog,"#******");
     for(j=1;j<=cptcoveff ;j++) {/* all covariates */        for(j=1;j<=cptcoveff ;j++) {/* all covariates */
       fprintf(ficrespl," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]); /* Here problem for varying dummy*/          fprintf(ficrespl," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]); /* Here problem for varying dummy*/
       printf(" V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);          printf(" V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
       fprintf(ficlog," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);          fprintf(ficlog," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
     }        }
     for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */        for (k4=1; k4<= nsq; k4++){ /* For each selected (single) quantitative value */
       printf(" V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);          printf(" V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
       fprintf(ficlog," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);          fprintf(ficrespl," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
     }          fprintf(ficlog," V%d=%f ",Tvqresult[nres][k4],Tqresult[nres][k4]);
     fprintf(ficrespl,"******\n");        }
     printf("******\n");        fprintf(ficrespl,"******\n");
     fprintf(ficlog,"******\n");        printf("******\n");
     if(invalidvarcomb[k]){        fprintf(ficlog,"******\n");
       printf("\nCombination (%d) ignored because no case \n",k);         if(invalidvarcomb[k]){
       fprintf(ficrespl,"#Combination (%d) ignored because no case \n",k);           printf("\nCombination (%d) ignored because no case \n",k); 
       fprintf(ficlog,"\nCombination (%d) ignored because no case \n",k);           fprintf(ficrespl,"#Combination (%d) ignored because no case \n",k); 
                                                 continue;          fprintf(ficlog,"\nCombination (%d) ignored because no case \n",k); 
     }          continue;
         }
   
     fprintf(ficrespl,"#Age ");        fprintf(ficrespl,"#Age ");
     for(j=1;j<=cptcoveff;j++) {        for(j=1;j<=cptcoveff;j++) {
       fprintf(ficrespl,"V%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);          fprintf(ficrespl,"V%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
     }        }
     for(i=1; i<=nlstate;i++) fprintf(ficrespl,"  %d-%d   ",i,i);        for(i=1; i<=nlstate;i++) fprintf(ficrespl,"  %d-%d   ",i,i);
     fprintf(ficrespl,"Total Years_to_converge\n");        fprintf(ficrespl,"Total Years_to_converge\n");
           
     for (age=agebase; age<=agelim; age++){        for (age=agebase; age<=agelim; age++){
       /* for (age=agebase; age<=agebase; age++){ */          /* for (age=agebase; age<=agebase; age++){ */
       prevalim(prlim, nlstate, p, age, oldm, savm, ftolpl, ncvyearp, k, nres);          prevalim(prlim, nlstate, p, age, oldm, savm, ftolpl, ncvyearp, k, nres);
       fprintf(ficrespl,"%.0f ",age );          fprintf(ficrespl,"%.0f ",age );
       for(j=1;j<=cptcoveff;j++)          for(j=1;j<=cptcoveff;j++)
         fprintf(ficrespl,"%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);            fprintf(ficrespl,"%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
       tot=0.;          tot=0.;
       for(i=1; i<=nlstate;i++){          for(i=1; i<=nlstate;i++){
         tot +=  prlim[i][i];            tot +=  prlim[i][i];
         fprintf(ficrespl," %.5f", prlim[i][i]);            fprintf(ficrespl," %.5f", prlim[i][i]);
       }          }
       fprintf(ficrespl," %.3f %d\n", tot, *ncvyearp);          fprintf(ficrespl," %.3f %d\n", tot, *ncvyearp);
     } /* Age */        } /* Age */
     /* was end of cptcod */        /* was end of cptcod */
   } /* cptcov */      } /* cptcov */
     } /* nres */
   return 0;    return 0;
 }  }
   
Line 9123  int back_prevalence_limit(double *p, dou Line 9151  int back_prevalence_limit(double *p, dou
   i1=pow(2,cptcoveff);    i1=pow(2,cptcoveff);
   if (cptcovn < 1){i1=1;}    if (cptcovn < 1){i1=1;}
       
   for(nres=1; nres <= nresult; nres++) /* For each resultline */    for(nres=1; nres <= nresult; nres++){ /* For each resultline */
   for(k=1; k<=i1;k++){ /* For any combination of dummy covariates, fixed and varying */      for(k=1; k<=i1;k++){ /* For any combination of dummy covariates, fixed and varying */
     if(TKresult[nres]!= k)        if(TKresult[nres]!= k)
       continue;          continue;
     //printf("cptcov=%d cptcod=%d codtab=%d\n",cptcov, cptcod,codtabm(cptcod,cptcov));        //printf("cptcov=%d cptcod=%d codtab=%d\n",cptcov, cptcod,codtabm(cptcod,cptcov));
     fprintf(ficresplb,"#******");        fprintf(ficresplb,"#******");
     printf("#******");        printf("#******");
     fprintf(ficlog,"#******");        fprintf(ficlog,"#******");
     for(j=1;j<=cptcoveff ;j++) {/* all covariates */        for(j=1;j<=cptcoveff ;j++) {/* all covariates */
       fprintf(ficresplb," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);          fprintf(ficresplb," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
       printf(" V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);          printf(" V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
       fprintf(ficlog," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);          fprintf(ficlog," V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
     }        }
     for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */        for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */
       printf(" V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);          printf(" V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
       fprintf(ficresplb," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);          fprintf(ficresplb," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
       fprintf(ficlog," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);          fprintf(ficlog," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
     }        }
     fprintf(ficresplb,"******\n");        fprintf(ficresplb,"******\n");
     printf("******\n");        printf("******\n");
     fprintf(ficlog,"******\n");        fprintf(ficlog,"******\n");
     if(invalidvarcomb[k]){        if(invalidvarcomb[k]){
       printf("\nCombination (%d) ignored because no cases \n",k);           printf("\nCombination (%d) ignored because no cases \n",k); 
       fprintf(ficresplb,"#Combination (%d) ignored because no cases \n",k);           fprintf(ficresplb,"#Combination (%d) ignored because no cases \n",k); 
       fprintf(ficlog,"\nCombination (%d) ignored because no cases \n",k);           fprintf(ficlog,"\nCombination (%d) ignored because no cases \n",k); 
       continue;          continue;
     }        }
           
     fprintf(ficresplb,"#Age ");        fprintf(ficresplb,"#Age ");
     for(j=1;j<=cptcoveff;j++) {        for(j=1;j<=cptcoveff;j++) {
       fprintf(ficresplb,"V%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);          fprintf(ficresplb,"V%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
     }        }
     for(i=1; i<=nlstate;i++) fprintf(ficresplb,"  %d-%d   ",i,i);        for(i=1; i<=nlstate;i++) fprintf(ficresplb,"  %d-%d   ",i,i);
     fprintf(ficresplb,"Total Years_to_converge\n");        fprintf(ficresplb,"Total Years_to_converge\n");
           
           
     for (age=agebase; age<=agelim; age++){        for (age=agebase; age<=agelim; age++){
       /* for (age=agebase; age<=agebase; age++){ */          /* for (age=agebase; age<=agebase; age++){ */
       if(mobilavproj > 0){          if(mobilavproj > 0){
         /* bprevalim(bprlim, mobaverage, nlstate, p, age, ageminpar, agemaxpar, oldm, savm, doldm, dsavm, ftolpl, ncvyearp, k); */            /* bprevalim(bprlim, mobaverage, nlstate, p, age, ageminpar, agemaxpar, oldm, savm, doldm, dsavm, ftolpl, ncvyearp, k); */
         /* bprevalim(bprlim, mobaverage, nlstate, p, age, oldm, savm, dnewm, doldm, dsavm, ftolpl, ncvyearp, k); */            /* bprevalim(bprlim, mobaverage, nlstate, p, age, oldm, savm, dnewm, doldm, dsavm, ftolpl, ncvyearp, k); */
         bprevalim(bprlim, mobaverage, nlstate, p, age, ftolpl, ncvyearp, k);            bprevalim(bprlim, mobaverage, nlstate, p, age, ftolpl, ncvyearp, k);
       }else if (mobilavproj == 0){          }else if (mobilavproj == 0){
         printf("There is no chance to get back prevalence limit if data aren't non zero and summing to 1, please try a non null mobil_average(=%d) parameter or mobil_average=-1 if you want to try at your own risk.\n",mobilavproj);            printf("There is no chance to get back prevalence limit if data aren't non zero and summing to 1, please try a non null mobil_average(=%d) parameter or mobil_average=-1 if you want to try at your own risk.\n",mobilavproj);
         fprintf(ficlog,"There is no chance to get back prevalence limit if data aren't non zero and summing to 1, please try a non null mobil_average(=%d) parameter or mobil_average=-1 if you want to try at your own risk.\n",mobilavproj);            fprintf(ficlog,"There is no chance to get back prevalence limit if data aren't non zero and summing to 1, please try a non null mobil_average(=%d) parameter or mobil_average=-1 if you want to try at your own risk.\n",mobilavproj);
         exit(1);            exit(1);
       }else{          }else{
         /* bprevalim(bprlim, probs, nlstate, p, age, oldm, savm, dnewm, doldm, dsavm, ftolpl, ncvyearp, k); */            /* bprevalim(bprlim, probs, nlstate, p, age, oldm, savm, dnewm, doldm, dsavm, ftolpl, ncvyearp, k); */
         bprevalim(bprlim, probs, nlstate, p, age, ftolpl, ncvyearp, k);            bprevalim(bprlim, probs, nlstate, p, age, ftolpl, ncvyearp, k);
       }          }
       fprintf(ficresplb,"%.0f ",age );          fprintf(ficresplb,"%.0f ",age );
       for(j=1;j<=cptcoveff;j++)          for(j=1;j<=cptcoveff;j++)
         fprintf(ficresplb,"%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);            fprintf(ficresplb,"%d %d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
       tot=0.;          tot=0.;
       for(i=1; i<=nlstate;i++){          for(i=1; i<=nlstate;i++){
         tot +=  bprlim[i][i];            tot +=  bprlim[i][i];
         fprintf(ficresplb," %.5f", bprlim[i][i]);            fprintf(ficresplb," %.5f", bprlim[i][i]);
       }          }
       fprintf(ficresplb," %.3f %d\n", tot, *ncvyearp);          fprintf(ficresplb," %.3f %d\n", tot, *ncvyearp);
     } /* Age */        } /* Age */
     /* was end of cptcod */        /* was end of cptcod */
   } /* cptcov */      } /* end of any combination */
       } /* end of nres */  
   /* hBijx(p, bage, fage); */    /* hBijx(p, bage, fage); */
   /* fclose(ficrespijb); */    /* fclose(ficrespijb); */
       
Line 9277  int hPijx(double *p, int bage, int fage) Line 9305  int hPijx(double *p, int bage, int fage)
         int ageminl;          int ageminl;
   int hstepm;    int hstepm;
   int nhstepm;    int nhstepm;
   int h, i, i1, j, k;    int h, i, i1, j, k, nres;
                   
   double agedeb;    double agedeb;
   double ***p3mat;    double ***p3mat;
Line 9305  int hPijx(double *p, int bage, int fage) Line 9333  int hPijx(double *p, int bage, int fage)
   /* for(cptcov=1,k=0;cptcov<=i1;cptcov++){ */    /* for(cptcov=1,k=0;cptcov<=i1;cptcov++){ */
   /*    /\*for(cptcod=1;cptcod<=ncodemax[cptcov];cptcod++){*\/ */    /*    /\*for(cptcod=1;cptcod<=ncodemax[cptcov];cptcod++){*\/ */
   /*    k=k+1;  */    /*    k=k+1;  */
   for (k=1; k <= (int) pow(2,cptcoveff); k++){    for(nres=1; nres <= nresult; nres++){ /* For each resultline */
     fprintf(ficrespijb,"\n#****** ");      for(k=1; k<=i1;k++){ /* For any combination of dummy covariates, fixed and varying */
     for(j=1;j<=cptcoveff;j++)        if(TKresult[nres]!= k)
       fprintf(ficrespijb,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);          continue;
     fprintf(ficrespijb,"******\n");        fprintf(ficrespijb,"\n#****** ");
     if(invalidvarcomb[k]){        for(j=1;j<=cptcoveff;j++)
       fprintf(ficrespijb,"\n#Combination (%d) ignored because no cases \n",k);           fprintf(ficrespijb,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
       continue;        for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */
     }          fprintf(ficrespijb," V%d=%f ",Tvqresult[nres][j],Tqresult[nres][j]);
             }
     /* for (agedeb=fage; agedeb>=bage; agedeb--){ /\* If stepm=6 months *\/ */        fprintf(ficrespijb,"******\n");
     for (agedeb=bage; agedeb<=fage; agedeb++){ /* If stepm=6 months and estepm=24 (2 years) */        if(invalidvarcomb[k]){
       /* nhstepm=(int) rint((agelim-agedeb)*YEARM/stepm); /\* Typically 20 years = 20*12/6=40 *\/ */          fprintf(ficrespijb,"\n#Combination (%d) ignored because no cases \n",k); 
       nhstepm=(int) rint((agedeb-ageminl)*YEARM/stepm); /* Typically 20 years = 20*12/6=40 */          continue;
       nhstepm = nhstepm/hstepm; /* Typically 40/4=10, because estepm=24 stepm=6 => hstepm=24/6=4 */        }
               
       /*          nhstepm=nhstepm*YEARM; aff par mois*/        /* for (agedeb=fage; agedeb>=bage; agedeb--){ /\* If stepm=6 months *\/ */
               for (agedeb=bage; agedeb<=fage; agedeb++){ /* If stepm=6 months and estepm=24 (2 years) */
       p3mat=ma3x(1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm);          /* nhstepm=(int) rint((agelim-agedeb)*YEARM/stepm); /\* Typically 20 years = 20*12/6=40 *\/ */
       /* oldm=oldms;savm=savms; */          nhstepm=(int) rint((agedeb-ageminl)*YEARM/stepm); /* Typically 20 years = 20*12/6=40 */
       /* hbxij(p3mat,nhstepm,agedeb,hstepm,p,nlstate,stepm,oldm,savm, k);   */          nhstepm = nhstepm/hstepm; /* Typically 40/4=10, because estepm=24 stepm=6 => hstepm=24/6=4 */
       hbxij(p3mat,nhstepm,agedeb,hstepm,p,prevacurrent,nlstate,stepm, k);          
       /* hbxij(p3mat,nhstepm,agedeb,hstepm,p,prevacurrent,nlstate,stepm,oldm,savm, dnewm, doldm, dsavm, k); */          /*        nhstepm=nhstepm*YEARM; aff par mois*/
       fprintf(ficrespijb,"# Cov Agex agex-h hpijx with i,j=");          
       for(i=1; i<=nlstate;i++)          p3mat=ma3x(1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm);
         for(j=1; j<=nlstate+ndeath;j++)          /* oldm=oldms;savm=savms; */
           fprintf(ficrespijb," %1d-%1d",i,j);          /* hbxij(p3mat,nhstepm,agedeb,hstepm,p,nlstate,stepm,oldm,savm, k);   */
       fprintf(ficrespijb,"\n");          hbxij(p3mat,nhstepm,agedeb,hstepm,p,prevacurrent,nlstate,stepm, k);
       for (h=0; h<=nhstepm; h++){          /* hbxij(p3mat,nhstepm,agedeb,hstepm,p,prevacurrent,nlstate,stepm,oldm,savm, dnewm, doldm, dsavm, k); */
         /*agedebphstep = agedeb + h*hstepm/YEARM*stepm;*/          fprintf(ficrespijb,"# Cov Agex agex-h hpijx with i,j=");
         fprintf(ficrespijb,"%d %3.f %3.f",k, agedeb, agedeb - h*hstepm/YEARM*stepm );  
         /* fprintf(ficrespijb,"%d %3.f %3.f",k, agedeb, agedeb + h*hstepm/YEARM*stepm ); */  
         for(i=1; i<=nlstate;i++)          for(i=1; i<=nlstate;i++)
           for(j=1; j<=nlstate+ndeath;j++)            for(j=1; j<=nlstate+ndeath;j++)
             fprintf(ficrespijb," %.5f", p3mat[i][j][h]);              fprintf(ficrespijb," %1d-%1d",i,j);
         fprintf(ficrespijb,"\n");          fprintf(ficrespijb,"\n");
       }          for (h=0; h<=nhstepm; h++){
       free_ma3x(p3mat,1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm);            /*agedebphstep = agedeb + h*hstepm/YEARM*stepm;*/
       fprintf(ficrespijb,"\n");            fprintf(ficrespijb,"%d %3.f %3.f",k, agedeb, agedeb - h*hstepm/YEARM*stepm );
     }            /* fprintf(ficrespijb,"%d %3.f %3.f",k, agedeb, agedeb + h*hstepm/YEARM*stepm ); */
     /*}*/            for(i=1; i<=nlstate;i++)
   }              for(j=1; j<=nlstate+ndeath;j++)
                 fprintf(ficrespijb," %.5f", p3mat[i][j][h]);
             fprintf(ficrespijb,"\n");
           }
           free_ma3x(p3mat,1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm);
           fprintf(ficrespijb,"\n");
         } /* end age deb */
       } /* end combination */
     } /* end nres */
   return 0;    return 0;
  } /*  hBijx */   } /*  hBijx */
   
Line 9975  Please run with mle=-1 to get a correct Line 10009  Please run with mle=-1 to get a correct
   Typevar=ivector(-1,NCOVMAX); /* -1 to 2 */    Typevar=ivector(-1,NCOVMAX); /* -1 to 2 */
   Fixed=ivector(-1,NCOVMAX); /* -1 to 3 */    Fixed=ivector(-1,NCOVMAX); /* -1 to 3 */
   Dummy=ivector(-1,NCOVMAX); /* -1 to 3 */    Dummy=ivector(-1,NCOVMAX); /* -1 to 3 */
     DummyV=ivector(1,NCOVMAX); /* 1 to 3 */
     FixedV=ivector(1,NCOVMAX); /* 1 to 3 */
   /*  V2+V1+V4+age*V3 is a model with 4 covariates (3 plus signs).     /*  V2+V1+V4+age*V3 is a model with 4 covariates (3 plus signs). 
       For each model-covariate stores the data-covariate id. Tvar[1]=2, Tvar[2]=1, Tvar[3]=4,         For each model-covariate stores the data-covariate id. Tvar[1]=2, Tvar[2]=1, Tvar[3]=4, 
       Tvar[4=age*V3] is 3 and 'age' is recorded in Tage.        Tvar[4=age*V3] is 3 and 'age' is recorded in Tage.
Line 10751  Please run with mle=-1 to get a correct Line 10787  Please run with mle=-1 to get a correct
       ungetc(c,ficpar);        ungetc(c,ficpar);
       fgets(line, MAXLINE, ficpar);        fgets(line, MAXLINE, ficpar);
       fputs(line,stdout);        fputs(line,stdout);
         fputs(line,ficres);
       fputs(line,ficparo);        fputs(line,ficparo);
     }      }
     ungetc(c,ficpar);      ungetc(c,ficpar);
Line 10767  Please run with mle=-1 to get a correct Line 10804  Please run with mle=-1 to get a correct
       fgets(line, MAXLINE, ficpar);        fgets(line, MAXLINE, ficpar);
       fputs(line,stdout);        fputs(line,stdout);
       fputs(line,ficparo);        fputs(line,ficparo);
         fputs(line,ficres);
     }      }
     ungetc(c,ficpar);      ungetc(c,ficpar);
           
Line 10785  Please run with mle=-1 to get a correct Line 10823  Please run with mle=-1 to get a correct
         fputs(line,stdout);          fputs(line,stdout);
         fputs(line,ficparo);          fputs(line,ficparo);
         fputs(line,ficlog);          fputs(line,ficlog);
           fputs(line,ficres);
         continue;          continue;
       }else        }else
         break;          break;
Line 10803  Please run with mle=-1 to get a correct Line 10842  Please run with mle=-1 to get a correct
         goto end;          goto end;
       }        }
       decoderesult(resultline, nresult); /* Fills TKresult[nresult] combination and Tresult[nresult][k4+1] combination values */        decoderesult(resultline, nresult); /* Fills TKresult[nresult] combination and Tresult[nresult][k4+1] combination values */
         fprintf(ficparo,"result: %s\n",resultline);
         fprintf(ficres,"result: %s\n",resultline);
         fprintf(ficlog,"result: %s\n",resultline);
       while(fgets(line, MAXLINE, ficpar)) {        while(fgets(line, MAXLINE, ficpar)) {
         /* If line starts with a # it is a comment */          /* If line starts with a # it is a comment */
         if (line[0] == '#') {          if (line[0] == '#') {
           numlinepar++;            numlinepar++;
           fputs(line,stdout);            fputs(line,stdout);
           fputs(line,ficparo);            fputs(line,ficparo);
             fputs(line,ficres);
           fputs(line,ficlog);            fputs(line,ficlog);
           continue;            continue;
         }else          }else
Line 10953  Please run with mle=-1 to get a correct Line 10996  Please run with mle=-1 to get a correct
     }      }
     printf("Computing Health Expectancies: result on file '%s' ...", filerese);fflush(stdout);      printf("Computing Health Expectancies: result on file '%s' ...", filerese);fflush(stdout);
     fprintf(ficlog,"Computing Health Expectancies: result on file '%s' ...", filerese);fflush(ficlog);      fprintf(ficlog,"Computing Health Expectancies: result on file '%s' ...", filerese);fflush(ficlog);
   
       pstamp(ficreseij);
                                   
     i1=pow(2,cptcoveff); /* Number of combination of dummy covariates */      i1=pow(2,cptcoveff); /* Number of combination of dummy covariates */
     if (cptcovn < 1){i1=1;}      if (cptcovn < 1){i1=1;}
Line 11065  Please run with mle=-1 to get a correct Line 11110  Please run with mle=-1 to get a correct
       fprintf(ficrescveij,"******\n");        fprintf(ficrescveij,"******\n");
               
       fprintf(ficresvij,"\n#****** ");        fprintf(ficresvij,"\n#****** ");
         /* pstamp(ficresvij); */
       for(j=1;j<=cptcoveff;j++)         for(j=1;j<=cptcoveff;j++) 
         fprintf(ficresvij,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);          fprintf(ficresvij,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]);
       for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */        for (j=1; j<= nsq; j++){ /* For each selected (single) quantitative value */
Line 11242  Please run with mle=-1 to get a correct Line 11288  Please run with mle=-1 to get a correct
   free_ivector(ncodemaxwundef,1,NCOVMAX);    free_ivector(ncodemaxwundef,1,NCOVMAX);
   free_ivector(Dummy,-1,NCOVMAX);    free_ivector(Dummy,-1,NCOVMAX);
   free_ivector(Fixed,-1,NCOVMAX);    free_ivector(Fixed,-1,NCOVMAX);
     free_ivector(DummyV,1,NCOVMAX);
     free_ivector(FixedV,1,NCOVMAX);
   free_ivector(Typevar,-1,NCOVMAX);    free_ivector(Typevar,-1,NCOVMAX);
   free_ivector(Tvar,1,NCOVMAX);    free_ivector(Tvar,1,NCOVMAX);
   free_ivector(TvarsQ,1,NCOVMAX);    free_ivector(TvarsQ,1,NCOVMAX);

Removed from v.1.237  
changed lines
  Added in v.1.238


FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>