[BACK]Return to imach.c CVS log [TXT] [DIR] Up to [Local Repository] / imach / src

Diff for /imach/src/imach.c between versions 1.118 and 1.119

version 1.118, 2006/03/14 18:20:07 version 1.119, 2006/03/15 17:42:26
Line 1 Line 1
 /* $Id$  /* $Id$
   $State$    $State$
   $Log$    $Log$
     Revision 1.119  2006/03/15 17:42:26  brouard
     (Module): Bug if status = -2, the loglikelihood was
     computed as likelihood omitting the logarithm. Version O.98e
   
   Revision 1.118  2006/03/14 18:20:07  brouard    Revision 1.118  2006/03/14 18:20:07  brouard
   (Module): varevsij Comments added explaining the second    (Module): varevsij Comments added explaining the second
   table of variances if popbased=1 .    table of variances if popbased=1 .
Line 325  extern int errno; Line 329  extern int errno;
 /* $Id$ */  /* $Id$ */
 /* $State$ */  /* $State$ */
   
 char version[]="Imach version 0.98d, March 2006, INED-EUROREVES-Institut de longevite ";  char version[]="Imach version 0.98e, March 2006, INED-EUROREVES-Institut de longevite ";
 char fullversion[]="$Revision$ $Date$";   char fullversion[]="$Revision$ $Date$"; 
 char strstart[80];  char strstart[80];
 char optionfilext[10], optionfilefiname[FILENAMELENGTH];  char optionfilext[10], optionfilefiname[FILENAMELENGTH];
Line 1383  double func( double *x) Line 1387  double func( double *x)
         } else if  (s2==-2) {          } else if  (s2==-2) {
           for (j=1,survp=0. ; j<=nlstate; j++)             for (j=1,survp=0. ; j<=nlstate; j++) 
             survp += out[s1][j];              survp += out[s1][j];
           lli= survp;            lli= log(survp);
         }          }
                   
         else if  (s2==-4) {  /*      else if  (s2==-4) { */
           for (j=3,survp=0. ; j<=nlstate; j++)   /*        for (j=3,survp=0. ; j<=nlstate; j++)  */
             survp += out[s1][j];  /*          survp += out[s1][j]; */
           lli= survp;  /*        lli= survp; */
         }  /*      } */
                   
         else if  (s2==-5) {  /*      else if  (s2==-5) { */
           for (j=1,survp=0. ; j<=2; j++)   /*        for (j=1,survp=0. ; j<=2; j++)  */
             survp += out[s1][j];  /*          survp += out[s1][j]; */
           lli= survp;  /*        lli= survp; */
         }  /*      } */
   
   
         else{          else{
Line 1593  double funcone( double *x) Line 1597  double funcone( double *x)
        */         */
       if( s2 > nlstate && (mle <5) ){  /* Jackson */        if( s2 > nlstate && (mle <5) ){  /* Jackson */
         lli=log(out[s1][s2] - savm[s1][s2]);          lli=log(out[s1][s2] - savm[s1][s2]);
       } else if (mle==1){        } else if  (s2==-2) {
           for (j=1,survp=0. ; j<=nlstate; j++) 
             survp += out[s1][j];
           lli= log(survp);
         }else if (mle==1){
         lli= log((1.+bbh)*out[s1][s2]- bbh*savm[s1][s2]); /* linear interpolation */          lli= log((1.+bbh)*out[s1][s2]- bbh*savm[s1][s2]); /* linear interpolation */
       } else if(mle==2){        } else if(mle==2){
         lli= (savm[s1][s2]>(double)1.e-8 ?log((1.+bbh)*out[s1][s2]- bbh*savm[s1][s2]):log((1.+bbh)*out[s1][s2])); /* linear interpolation */          lli= (savm[s1][s2]>(double)1.e-8 ?log((1.+bbh)*out[s1][s2]- bbh*savm[s1][s2]):log((1.+bbh)*out[s1][s2])); /* linear interpolation */
Line 1609  double funcone( double *x) Line 1617  double funcone( double *x)
       ll[s[mw[mi][i]][i]] += 2*weight[i]*lli;        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]); */  /*       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]); */
       if(globpr){        if(globpr){
         fprintf(ficresilk,"%9d %6d %1d %1d %1d %1d %3d %10.6f %6.4f\          fprintf(ficresilk,"%9d %6d %2d %2d %1d %1d %3d %11.6f %8.4f\
  %10.6f %10.6f %10.6f ", \   %11.6f %11.6f %11.6f ", \
                 num[i],i,s1,s2,mi,mw[mi][i],dh[mi][i],exp(lli),weight[i],                  num[i],i,s1,s2,mi,mw[mi][i],dh[mi][i],exp(lli),weight[i],
                 2*weight[i]*lli,out[s1][s2],savm[s1][s2]);                  2*weight[i]*lli,out[s1][s2],savm[s1][s2]);
         for(k=1,llt=0.,l=0.; k<=nlstate; k++){          for(k=1,llt=0.,l=0.; k<=nlstate; k++){
Line 2475  void evsij(char fileres[], double ***eij Line 2483  void evsij(char fileres[], double ***eij
   for (age=bage; age<=fage; age ++){ /* If stepm=6 months */    for (age=bage; age<=fage; age ++){ /* If stepm=6 months */
     /* Computed by stepm unit matrices, product of hstepm matrices, stored      /* Computed by stepm unit matrices, product of hstepm matrices, stored
        in an array of nhstepm length: nhstepm=10, hstepm=4, stepm=6 months */         in an array of nhstepm length: nhstepm=10, hstepm=4, stepm=6 months */
       
     hpxij(p3mat,nhstepm,age,hstepm,x,nlstate,stepm,oldm, savm, cij);        hpxij(p3mat,nhstepm,age,hstepm,x,nlstate,stepm,oldm, savm, cij);  
        
     hf=hstepm*stepm/YEARM;  /* Duration of hstepm expressed in year unit. */      hf=hstepm*stepm/YEARM;  /* Duration of hstepm expressed in year unit. */
       
     /* Computing  Variances of health expectancies */      printf("%d|",(int)age);fflush(stdout);
     /* Gradient is computed with plus gp and minus gm. Code is duplicated in order to      fprintf(ficlog,"%d|",(int)age);fflush(ficlog);
        decrease memory allocation */      
      printf("%d|",(int)age);fflush(stdout);  
      fprintf(ficlog,"%d|",(int)age);fflush(ficlog);  
     /* Computing expectancies */      /* Computing expectancies */
     for(i=1; i<=nlstate;i++)      for(i=1; i<=nlstate;i++)
       for(j=1; j<=nlstate;j++)        for(j=1; j<=nlstate;j++)
         for (h=0, eij[i][j][(int)age]=0; h<=nhstepm-1; h++){          for (h=0, eij[i][j][(int)age]=0; h<=nhstepm-1; h++){
           eij[i][j][(int)age] += (p3mat[i][j][h]+p3mat[i][j][h+1])/2.0*hf;            eij[i][j][(int)age] += (p3mat[i][j][h]+p3mat[i][j][h+1])/2.0*hf;
                       
 /* if((int)age==70)printf("i=%2d,j=%2d,h=%2d,age=%3d,%9.4f,%9.4f,%9.4f\n",i,j,h,(int)age,p3mat[i][j][h],hf,eij[i][j][(int)age]);*/            /* if((int)age==70)printf("i=%2d,j=%2d,h=%2d,age=%3d,%9.4f,%9.4f,%9.4f\n",i,j,h,(int)age,p3mat[i][j][h],hf,eij[i][j][(int)age]);*/
   
         }          }
   
Line 2505  void evsij(char fileres[], double ***eij Line 2511  void evsij(char fileres[], double ***eij
       fprintf(ficreseij,"%9.4f", eip );        fprintf(ficreseij,"%9.4f", eip );
     }      }
     fprintf(ficreseij,"\n");      fprintf(ficreseij,"\n");
       
   }    }
   free_ma3x(p3mat,1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm);    free_ma3x(p3mat,1,nlstate+ndeath,1, nlstate+ndeath, 0,nhstepm);
   printf("\n");    printf("\n");
   fprintf(ficlog,"\n");    fprintf(ficlog,"\n");
     
 }  }
   
 void cvevsij(char fileres[], double ***eij, double x[], int nlstate, int stepm, int bage, int fage, double **oldm, double **savm, int cij, int estepm,double delti[],double **matcov,char strstart[] )  void cvevsij(char fileres[], double ***eij, double x[], int nlstate, int stepm, int bage, int fage, double **oldm, double **savm, int cij, int estepm,double delti[],double **matcov,char strstart[] )
Line 4411  int main(int argc, char *argv[]) Line 4417  int main(int argc, char *argv[])
       while ((val = strsep(&tok, "\"" )) != NULL && *val == '\0');        while ((val = strsep(&tok, "\"" )) != NULL && *val == '\0');
       printf("val= |%s| pathr=%s\n",val,pathr);        printf("val= |%s| pathr=%s\n",val,pathr);
       strcpy (pathtot, val);        strcpy (pathtot, val);
       if(pathr[0] == '\0') break; /* Un peu sale */        if(pathr[0] == '\0') break; /* Dirty */
     }      }
   }    }
   else{    else{
Line 4430  int main(int argc, char *argv[]) Line 4436  int main(int argc, char *argv[])
   /* Split argv[1]=pathtot, parameter file name to get path, optionfile, extension and name */    /* Split argv[1]=pathtot, parameter file name to get path, optionfile, extension and name */
   split(pathtot,path,optionfile,optionfilext,optionfilefiname);    split(pathtot,path,optionfile,optionfilext,optionfilefiname);
   printf("\npathtot=%s,\npath=%s,\noptionfile=%s \noptionfilext=%s \noptionfilefiname=%s\n",pathtot,path,optionfile,optionfilext,optionfilefiname);    printf("\npathtot=%s,\npath=%s,\noptionfile=%s \noptionfilext=%s \noptionfilefiname=%s\n",pathtot,path,optionfile,optionfilext,optionfilefiname);
   chdir(path);    chdir(path); /* Can be a relative path */
     if(getcwd(pathcd,MAXLINE) > 0) /* So pathcd is the full path */
       printf("Current directory %s!\n",pathcd);
   strcpy(command,"mkdir ");    strcpy(command,"mkdir ");
   strcat(command,optionfilefiname);    strcat(command,optionfilefiname);
   if((outcmd=system(command)) != 0){    if((outcmd=system(command)) != 0){
Line 5263  Interval (in months) between two waves: Line 5271  Interval (in months) between two waves:
       printf("%d %.0lf %lf %.0lf %.0lf %.0lf %lf\n",k,lsurv[k],p[1]*exp(p[2]*(k-agegomp)),(p[1]*exp(p[2]*(k-agegomp)))*lsurv[k],lpop[k],tpop[k],tpop[k]/lsurv[k]);        printf("%d %.0lf %lf %.0lf %.0lf %.0lf %lf\n",k,lsurv[k],p[1]*exp(p[2]*(k-agegomp)),(p[1]*exp(p[2]*(k-agegomp)))*lsurv[k],lpop[k],tpop[k],tpop[k]/lsurv[k]);
           
           
     replace_back_to_slash(pathc,path); /* Even gnuplot wants a / */      replace_back_to_slash(pathc,pathcd); /* Even gnuplot wants a / */
     printinggnuplotmort(fileres, optionfilefiname,ageminpar,agemaxpar,fage, pathc,p);      printinggnuplotmort(fileres, optionfilefiname,ageminpar,agemaxpar,fage, pathc,p);
           
     printinghtmlmort(fileres,title,datafile, firstpass, lastpass, \      printinghtmlmort(fileres,title,datafile, firstpass, lastpass, \
Line 5500  Interval (in months) between two waves: Line 5508  Interval (in months) between two waves:
     /*  freqsummary(fileres, agemin, agemax, s, agev, nlstate, imx,Tvaraff,nbcode, ncodemax,mint,anint);*/      /*  freqsummary(fileres, agemin, agemax, s, agev, nlstate, imx,Tvaraff,nbcode, ncodemax,mint,anint);*/
     /*,dateprev1,dateprev2,jprev1, mprev1,anprev1,jprev2, mprev2,anprev2);*/      /*,dateprev1,dateprev2,jprev1, mprev1,anprev1,jprev2, mprev2,anprev2);*/
           
     replace_back_to_slash(pathc,path); /* Even gnuplot wants a / */      replace_back_to_slash(pathc,pathcd); /* Even gnuplot wants a / */
     printinggnuplot(fileres, optionfilefiname,ageminpar,agemaxpar,fage, pathc,p);      printinggnuplot(fileres, optionfilefiname,ageminpar,agemaxpar,fage, pathc,p);
           
     printinghtml(fileres,title,datafile, firstpass, lastpass, stepm, weightopt,\      printinghtml(fileres,title,datafile, firstpass, lastpass, stepm, weightopt,\
Line 5895  Interval (in months) between two waves: Line 5903  Interval (in months) between two waves:
   fclose(ficlog);    fclose(ficlog);
   /*------ End -----------*/    /*------ End -----------*/
   
   chdir(path);  
      printf("Before Current directory %s!\n",pathcd);
      if(chdir(pathcd) != 0)
       printf("Can't move to directory %s!\n",path);
     if(getcwd(pathcd,MAXLINE) > 0)
       printf("Current directory %s!\n",pathcd);
   /*strcat(plotcmd,CHARSEPARATOR);*/    /*strcat(plotcmd,CHARSEPARATOR);*/
   sprintf(plotcmd,"gnuplot");    sprintf(plotcmd,"gnuplot");
 #ifndef UNIX  #ifndef UNIX
Removed from v.1.118  
changed lines
  Added in v.1.119

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