src/Tools/ctfInfo/ctfDeterminationByBayes/src/ctfDeterminationByBayes.c

   1 /*
   2 # ctfDeterminationByBayes : $Revision$
   3 # $Date$
   4 # Created by $Author$
   5 # Usage : ctfDeterminationByBayes
   6 # Attention
   7 #   $Loccker$
   8 #       $State$
   9 #
  10 */
  11 #include <stdio.h>
  12 #include <stdlib.h>
  13 #include <string.h>
  14 #include <math.h>
  15 #define GLOBAL_DECLARATION
  16 #include "../inc/config.h"
  17 #include "eosBayes.h"
  18 #include "ctfInfo.h"
  19 #include "ctffuncforbayes.h"
  20
  21
  22 #define DEBUG
  23 #include "genUtil.h"
  24 #include "ctfDeterminationByBayes.h"
  25 #include "String.h"
  26
  27 /*
  28 Example:
  29 typedef struct lctfDeterminationByBayesInfo {
  30         float a;
  31         int   b;
  32 } lctfDeterminationByBayesInfo;
  33
  34 typedef enum lctfDeterminationByBayesMode {
  35         a=0,
  36         b=1
  37 } lctfDeterminationByBayesMode;
  38 */
  39
  40 void DefocusEstimationLoop(ctfDeterminationByBayesInfo info, eosBayes* out, ctfInfo* ctfinfo, ctfforbayes* ctf, intensityforbayes* intensity, defocuslikelihoodforbayes* likelihood, int count, int mode);
  41
  42 int
  43 main(int argc, char* argv[])
  44 {
  45
  46   int i;
  47   int j;
  48
  49   ctfDeterminationByBayesInfo info;
  50   ctfInfo ctfinfo;
  51   eosBayesDistributionInfo llinfo;
  52
  53   eosBayes out;
  54
  55   ctfforbayes ctf;
  56   intensityforbayes intensity;
  57   defocuslikelihoodforbayes likelihood;
  58
  59   float Rmax;
  60   float Rmin;
  61   float Rintermediate;
  62
  63   init0(&info);
  64   argCheck(&info, argc, argv);
  65   init1(&info);
  66
  67   DEBUGPRINT("Program Start\n");
  68
  69   ctfinfo.mode = 33;
  70
  71   /*For Prior*/
  72   leosBayesInit(&out, info.numDistribution, info.numLikelihood, info.mode);
  73   eosBayesCondition(info.fptIn, &out, &llinfo, 0);
  74   eosBayesPriorWrite(&out,0);
  75   DistributionChangeCheck(&out, 0, 0);
  76
  77
  78   /*For CTF and intensity*/
  79   ctfforbayesInit(info.fptIn2,&ctf,0);
  80   intensityforbayesInit(&intensity,ctf.n,info.MTFmode,info.Ainmode,0);
  81   spatialfrequencyRead(info.fptIn2,&ctf,info.rmax,info.rmin,0);
  82   ctfInfoRead(info.fptIn3,&ctfinfo,"",0);
  83   DEBUGPRINT6("kV:%f  Cs:%f  Ain:%f Cc:%f MTF:%f ctfMode:%ld\n",ctfinfo.kV,ctfinfo.Cs,ctfinfo.Ain,ctfinfo.Cc,ctfinfo.BofMTF,ctfinfo.mode);
  84   DEBUGPRINT1("flagIn4:%ld\n",info.flagIn4);
  85
  86   Rmax = ctf.rmax;
  87   Rmin = ctf.rmin;
  88   Rintermediate = info.rintermediate;
  89
  90   /*global carse search*/
  91   DEBUGPRINT("carse search start\n");
  92   ctf.rmax = Rmax;
  93   ctf.rmin = Rmin;
  94   DEBUGPRINT2("set R range:%f~%f\n",ctf.rmin,ctf.rmax);
  95   DefocusEstimationLoop(info,&out,&ctfinfo,&ctf,&intensity,&likelihood,0,0);
  96   DEBUGPRINT("carse search finish\n");
  97   ctfforbayescheckIntensityOfMaxPosterior(ctf,out.posterior,out.numDistribution,&ctfinfo,intensity.MTFmode,intensity.Ainmode,1,1);
  98
  99   /*MTF,B serch start\n*/
 100   DEBUGPRINT("MTF,B search start\n");
 101   ctf.rmax = Rmax;
 102   ctf.rmin = Rintermediate;
 103   DEBUGPRINT2("set R range:%f~%f\n",ctf.rmin,ctf.rmax);
 104
 105   DEBUGPRINT("set maxposteriorflag of df, A, k\n");
 106   eosBayesMaxposteriorflagSet(&out,0);
 107   eosBayesMaxposteriorflagSet(&out,1);
 108   eosBayesMaxposteriorflagSet(&out,3);
 109
 110   DEBUGPRINT("----maxflagcheck----\n");
 111   #ifdef DEBUG
 112   for(i=0;i<out.numDistribution;i++){
 113         printf("dist%d:%d\n",i,(out.posterior[i].maxposteriorflag));
 114   }
 115   #endif
 116
 117   DEBUGPRINT("set canstantflag of Ain");
 118   out.posterior[5].constantflag = 1;
 119
 120   DEBUGPRINT("set resetflag of \n");
 121   out.prior[4].resetflag = 1;
 122   eosBayesProbabilityReset(info.fptIn,&out,10,0);
 123
 124   DEBUGPRINT("set rangechangeflag of B\n");
 125   out.posterior[2].rangechangeflag = 1;
 126   eosBayesProbabilityRangeChange(&out,100,0.95,0);
 127
 128   DistributionChangeCheck(&out, 1, 0);
 129
 130   DefocusEstimationLoop(info,&out,&ctfinfo,&ctf,&intensity,&likelihood,1,1);
 131   DEBUGPRINT("MTF,B search finish\n");
 132
 133   DEBUGPRINT("----maxflagcheck----\n");
 134   #ifdef DEBUG
 135   for(i=0;i<out.numDistribution;i++){
 136         printf("dist%d:%d\n",i,(out.posterior[i].maxposteriorflag));
 137   }
 138   #endif
 139
 140   ctfforbayescheckIntensityOfMaxPosterior(ctf,out.posterior,out.numDistribution,&ctfinfo,intensity.MTFmode,intensity.Ainmode,1,1);
 141
 142
 143   /*df,A B serch*/
 144   DEBUGPRINT("df, A, k search start\n");
 145   ctf.rmax = Rintermediate;
 146   ctf.rmin = Rmin;
 147   DEBUGPRINT2("set R range:%f~%f\n",ctf.rmin,ctf.rmax);
 148
 149   DEBUGPRINT("set maxposteriorflag of MTF,B\n");
 150   eosBayesMaxposteriorflagSet(&out,4);
 151
 152   DEBUGPRINT("set canstantflag of Ain");
 153   out.posterior[3].constantflag = 1;
 154   out.posterior[5].constantflag = 1;
 155
 156   DEBUGPRINT("set resetflag of df,k\n");
 157   out.prior[0].resetflag = 1;
 158   out.prior[1].resetflag = 1;
 159   out.prior[2].resetflag = 1;
 160   eosBayesProbabilityReset(info.fptIn,&out,1,0);
 161
 162   DistributionChangeCheck(&out, 4, 0);
 163   DEBUGPRINT("----maxflagcheck----\n");
 164   #ifdef DEBUG
 165   for(i=0;i<out.numDistribution;i++){
 166         printf("dist%d:%d\n",i,(out.posterior[i].maxposteriorflag));
 167   }
 168   #endif
 169
 170   DefocusEstimationLoop(info,&out,&ctfinfo,&ctf,&intensity,&likelihood,4,1);
 171   DEBUGPRINT("A,df,k search finish\n");
 172   ctfforbayescheckIntensityOfMaxPosterior(ctf,out.posterior,out.numDistribution,&ctfinfo,intensity.MTFmode,intensity.Ainmode,1,1);
 173
 174
 175   /*df,A serch*/
 176
 177   DEBUGPRINT("df, A, k search start\n");
 178   ctf.rmax = Rintermediate;
 179   ctf.rmin = Rmin;
 180   DEBUGPRINT2("set R range:%f~%f\n",ctf.rmin,ctf.rmax);
 181
 182   DEBUGPRINT("set maxposteriorflag of MTF,B\n");
 183   eosBayesMaxposteriorflagSet(&out,4);
 184
 185   DEBUGPRINT("set canstantflag of Ain");
 186   out.posterior[5].constantflag = 1;
 187   out.posterior[3].constantflag = 1;
 188
 189   DEBUGPRINT("set resetflag of df,k\n");
 190   out.prior[0].resetflag = 1;
 191   eosBayesProbabilityReset(info.fptIn,&out,5,0);
 192
 193
 194   DEBUGPRINT("set rangechangeflag of A\n");
 195   out.posterior[1].rangechangeflag = 1;
 196   out.posterior[2].rangechangeflag = 1;
 197   eosBayesProbabilityRangeChange(&out,20,0.95,0);
 198
 199   DistributionChangeCheck(&out, 5, 0);
 200   DEBUGPRINT("----maxflagcheck----\n");
 201   #ifdef DEBUG
 202   for(i=0;i<out.numDistribution;i++){
 203         printf("dist%d:%d\n",i,(out.posterior[i].maxposteriorflag));
 204   }
 205   #endif
 206
 207   DefocusEstimationLoop(info,&out,&ctfinfo,&ctf,&intensity,&likelihood,5,1);
 208   DEBUGPRINT("A,df,k search finish\n");
 209   ctfforbayescheckIntensityOfMaxPosterior(ctf,out.posterior,out.numDistribution,&ctfinfo,intensity.MTFmode,intensity.Ainmode,1,1);
 210
 211   DEBUGPRINT("defocus refinement start\n");
 212   ctf.rmax = Rintermediate;
 213   ctf.rmin = Rmin;
 214   DEBUGPRINT2("set R range:%f~%f\n",ctf.rmin,ctf.rmax);
 215
 216   DEBUGPRINT("set maxposteriorflag of defocus,A,B,MTF\n");
 217   out.posterior[0].maxposteriorflag = 0;
 218   out.posterior[1].maxposteriorflag = 0;
 219   out.posterior[3].maxposteriorflag = 0;
 220   out.posterior[5].maxposteriorflag = 0;
 221   eosBayesMaxposteriorflagSet(&out,4);
 222
 223   DEBUGPRINT("set canstantflag of Ain");
 224   out.posterior[3].constantflag = 0;
 225   out.posterior[5].constantflag = 0;
 226
 227   DEBUGPRINT("set resetflag of df,k\n");
 228   out.prior[3].resetflag = 1;
 229   out.prior[5].resetflag = 1;
 230   eosBayesProbabilityReset(info.fptIn,&out,1,0);
 231
 232   DEBUGPRINT("set rangechangeflag of A\n");
 233   out.posterior[0].rangechangeflag = 1;
 234   out.posterior[1].rangechangeflag = 1;
 235   out.posterior[2].rangechangeflag = 1;
 236   eosBayesProbabilityRangeChange(&out,20,0.95,0);
 237
 238
 239   DistributionChangeCheck(&out, 6, 0);
 240   defocusLikelihoodInit(&likelihood,&out,1);
 241
 242   DefocusEstimationLoop(info,&out,&ctfinfo,&ctf,&intensity,&likelihood,6,1);
 243   DEBUGPRINT("defocus refinement finish\n");
 244   ctfforbayescheckIntensityOfMaxPosterior(ctf,out.posterior,out.numDistribution,&ctfinfo,intensity.MTFmode,intensity.Ainmode,1,1);
 245
 246   /*defocus refinement*/
 247 /*
 248   DEBUGPRINT("defocus refinement start\n");
 249   ctf.rmax = Rintermediate;
 250   ctf.rmin = Rmin;
 251   DEBUGPRINT2("set R range:%f~%f\n",ctf.rmin,ctf.rmax);
 252
 253   DEBUGPRINT("set maxposteriorflag of defocus,A,B,MTF\n");
 254   out.posterior[0].maxposteriorflag = 0;
 255   out.posterior[1].maxposteriorflag = 0;
 256   out.posterior[3].maxposteriorflag = 0;
 257   out.posterior[5].maxposteriorflag = 0;
 258   eosBayesMaxposteriorflagSet(&out,4);
 259
 260   DEBUGPRINT("set canstantflag of Ain");
 261   out.posterior[3].constantflag = 0;
 262   out.posterior[5].constantflag = 0;
 263
 264   DEBUGPRINT("set refmode\n");
 265   out.posterior[0].refmode = 1;
 266   out.posterior[1].refmode = 1;
 267   out.posterior[2].refmode = 1;
 268   out.posterior[3].refmode = 0;
 269   out.posterior[4].refmode = 0;
 270   out.posterior[5].refmode = 0;
 271
 272   eosBayesRefinementCondition(info.fptIn,&out,&llinfo,5,0);
 273   DistributionChangeCheck(&out, 6, 0);
 274   defocusLikelihoodInit(&likelihood,&out,1);
 275
 276   DefocusEstimationLoop(info,&out,&ctfinfo,&ctf,&intensity,&likelihood,6,1);
 277   DEBUGPRINT("defocus refinement finish\n");
 278   ctfforbayescheckIntensityOfMaxPosterior(ctf,out.posterior,out.numDistribution,&ctfinfo,intensity.MTFmode,intensity.Ainmode,1,1);
 279 */
 280   eosBayesPosteriorWrite(&out,0);
 281
 282   ctfforbayesFree(&ctf,&intensity,&likelihood,out.numLikelihood);
 283
 284   eosBayesFree(&out);
 285
 286   exit(EXIT_SUCCESS);
 287 }
 288
 289 void DefocusEstimationLoop(ctfDeterminationByBayesInfo info, eosBayes* out, ctfInfo* ctfinfo, ctfforbayes* ctf, intensityforbayes* intensity, defocuslikelihoodforbayes* likelihood, int count,int mode){
 290   int i;
 291
 292   defocusLikelihoodInit(likelihood,out,mode);
 293   DEBUGPRINT("In DefocusEstimationLoop\n");
 294   for(i=0;i<info.flagIn4;i++){
 295
 296         intensityRead(info.In4[i],intensity,i,0);
 297         DEBUGPRINT("ctf start\n");
 298         ctfFunctionforbayes(ctf,intensity,out,ctfinfo,likelihood,count,1);
 299
 300         DEBUGPRINT("estimate start\n");
 301         eosBayesEstimation(out,0);
 302         DistributionChangeCheck(out,count,1);
 303
 304         DEBUGPRINT("change start\n");
 305         eosBayesPosteriortoPrior(out,0);
 306
 307   }
 308
 309 }
 310
 311 void
 312 additionalUsage()
 313 {
 314   fprintf(stderr, "----- Additional Usage -----\n");
 315 }
 316