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 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,2);
 151   eosBayesMaxposteriorflagSet(&out,4);
 152
 153   DEBUGPRINT("set canstantflag of Ain");
 154 //  out.posterior[3].constantflag = 1;
 155   out.posterior[5].constantflag = 1;
 156
 157   DEBUGPRINT("set resetflag of df,k\n");
 158   out.prior[0].resetflag = 1;
 159   out.prior[1].resetflag = 1;
 160   out.prior[3].resetflag = 1;
 161   out.prior[5].resetflag = 1;
 162   eosBayesProbabilityReset(info.fptIn,&out,1,0);
 163
 164
 165   DEBUGPRINT("set rangechangeflag of A\n");
 166   out.posterior[2].rangechangeflag = 1;
 167   eosBayesProbabilityRangeChange(&out,30,0.98,0);
 168
 169   DistributionChangeCheck(&out, 4, 0);
 170   DEBUGPRINT("----maxflagcheck----\n");
 171   #ifdef DEBUG
 172   for(i=0;i<out.numDistribution;i++){
 173         printf("dist%d:%d\n",i,(out.posterior[i].maxposteriorflag));
 174   }
 175   #endif
 176
 177   DefocusEstimationLoop(info,&out,&ctfinfo,&ctf,&intensity,&likelihood,4,1);
 178   DEBUGPRINT("A,df,k search finish\n");
 179   ctfforbayescheckIntensityOfMaxPosterior(ctf,out.posterior,out.numDistribution,&ctfinfo,intensity.MTFmode,intensity.Ainmode,1,1);
 180
 181
 182   /*df,A serch*/
 183 /*
 184   DEBUGPRINT("df, A, k search start\n");
 185   ctf.rmax = Rintermediate;
 186   ctf.rmin = Rmin;
 187   DEBUGPRINT2("set R range:%f~%f\n",ctf.rmin,ctf.rmax);
 188
 189   DEBUGPRINT("set maxposteriorflag of MTF,B\n");
 190 //  eosBayesMaxposteriorflagSet(&out,2);
 191   eosBayesMaxposteriorflagSet(&out,4);
 192
 193   DEBUGPRINT("set canstantflag of Ain");
 194   out.posterior[5].constantflag = 1;
 195
 196   DEBUGPRINT("set resetflag of df,k\n");
 197   out.prior[0].resetflag = 1;
 198 //  out.prior[1].resetflag = 1;
 199 //  out.prior[2].resetflag = 1;
 200   out.prior[3].resetflag = 1;
 201   eosBayesProbabilityReset(info.fptIn,&out,1,0);
 202
 203
 204   DEBUGPRINT("set rangechangeflag of A\n");
 205   out.posterior[1].rangechangeflag = 1;
 206   out.posterior[2].rangechangeflag = 1;
 207   eosBayesProbabilityRangeChange(&out,40,0.98,0);
 208
 209   DistributionChangeCheck(&out, 2, 0);
 210   DEBUGPRINT("----maxflagcheck----\n");
 211   #ifdef DEBUG
 212   for(i=0;i<out.numDistribution;i++){
 213         printf("dist%d:%d\n",i,(out.posterior[i].maxposteriorflag));
 214   }
 215   #endif
 216
 217   DefocusEstimationLoop(info,&out,&ctfinfo,&ctf,&intensity,&likelihood,2,1);
 218   DEBUGPRINT("A,df,k search finish\n");
 219   ctfforbayescheckIntensityOfMaxPosterior(ctf,out.posterior,out.numDistribution,&ctfinfo,intensity.MTFmode,intensity.Ainmode,1,1);
 220 */
 221
 222
 223   /*df serch*/
 224 /*
 225   DEBUGPRINT("df serch start\n");
 226   ctf.rmax = Rintermediate;
 227   ctf.rmin = Rmin;
 228   DEBUGPRINT2("set R range:%f~%f\n",ctf.rmin,ctf.rmax);
 229
 230   DEBUGPRINT("set maxposteriorflag of defocus,A,B\n");
 231   eosBayesMaxposteriorflagSet(&out,1);
 232   eosBayesMaxposteriorflagSet(&out,2);
 233   eosBayesMaxposteriorflagSet(&out,3);
 234   eosBayesMaxposteriorflagSet(&out,4);
 235
 236   DEBUGPRINT("set canstantflag of Ain");
 237   out.posterior[5].constantflag = 1;
 238
 239   DEBUGPRINT("set resetflag of df\n");
 240   out.prior[0].resetflag = 1;
 241   eosBayesProbabilityReset(info.fptIn,&out,5,0);
 242   DistributionChangeCheck(&out, 3, 0);
 243
 244   DEBUGPRINT("----maxflagcheck----\n");
 245   #ifdef DEBUG
 246   for(i=0;i<out.numDistribution;i++){
 247         printf("dist%d:%d\n",i,(out.posterior[i].maxposteriorflag));
 248   }
 249   #endif
 250
 251   DefocusEstimationLoop(info,&out,&ctfinfo,&ctf,&intensity,&likelihood,3,1);
 252   DEBUGPRINT("df search finish\n");
 253   ctfforbayescheckIntensityOfMaxPosterior(ctf,out.posterior,out.numDistribution,&ctfinfo,intensity.MTFmode,intensity.Ainmode,1,1);
 254   */
 255
 256   /*MTF and B refinement*/
 257 /*
 258   DEBUGPRINT("MTF and B refinement\n");
 259   ctf.rmax = Rmax;
 260   ctf.rmin = Rintermediate;
 261   DEBUGPRINT2("set R range:%f~%f\n",ctf.rmin,ctf.rmax);
 262
 263   DEBUGPRINT("set maxposteriorflag of defocus,A,B,k\n");
 264   eosBayesMaxposteriorflagSet(&out,0);
 265   eosBayesMaxposteriorflagSet(&out,1);
 266   eosBayesMaxposteriorflagSet(&out,3);
 267
 268   DEBUGPRINT("set canstantflag of Ain");
 269   out.posterior[5].constantflag = 1;
 270
 271   DEBUGPRINT("set resetflag of MTF\n");
 272   out.prior[4].resetflag = 1;
 273   eosBayesProbabilityReset(info.fptIn,&out,5,0);
 274
 275   DEBUGPRINT("set rangechangeflag of B\n");
 276   out.posterior[2].rangechangeflag = 1;
 277   eosBayesProbabilityRangeChange(&out,100,0.9,0);
 278
 279   DistributionChangeCheck(&out, 4, 0);
 280
 281   DEBUGPRINT("----maxflagcheck----\n");
 282   #ifdef DEBUG
 283   for(i=0;i<out.numDistribution;i++){
 284         printf("dist%d:%d\n",i,(out.posterior[i].maxposteriorflag));
 285   }
 286   #endif
 287
 288   DefocusEstimationLoop(info,&out,&ctfinfo,&ctf,&intensity,&likelihood,4,1);
 289   DEBUGPRINT("MTF refinement finish\n");
 290   ctfforbayescheckIntensityOfMaxPosterior(ctf,out.posterior,out.numDistribution,&ctfinfo,intensity.MTFmode,intensity.Ainmode,1,1);
 291 <<<<<<< HEAD
 292
 293
 294 =======
 295 */
 296 >>>>>>> 2a62a3046b5d765cfdcc331395154d3f4a9e081f
 297   /*df serch 2*/
 298 /*
 299   DEBUGPRINT("df serch 2 start\n");
 300   ctf.rmax = Rintermediate;
 301   ctf.rmin = Rmin;
 302   DEBUGPRINT2("set R range:%f~%f\n",ctf.rmin,ctf.rmax);
 303
 304   DEBUGPRINT("set maxposteriorflag of defocus,A,B\n");
 305   eosBayesMaxposteriorflagSet(&out,1);
 306   eosBayesMaxposteriorflagSet(&out,2);
 307   eosBayesMaxposteriorflagSet(&out,3);
 308   eosBayesMaxposteriorflagSet(&out,4);
 309
 310   DEBUGPRINT("set canstantflag of Ain");
 311   out.posterior[5].constantflag = 1;
 312
 313   DEBUGPRINT("set resetflag of df\n");
 314   out.prior[0].resetflag = 1;
 315   eosBayesProbabilityReset(info.fptIn,&out,5,0);
 316   DistributionChangeCheck(&out, 5, 0);
 317
 318   DEBUGPRINT("----maxflagcheck----\n");
 319   #ifdef DEBUG
 320   for(i=0;i<out.numDistribution;i++){
 321         printf("dist%d:%d\n",i,(out.posterior[i].maxposteriorflag));
 322   }
 323   #endif
 324
 325   DefocusEstimationLoop(info,&out,&ctfinfo,&ctf,&intensity,&likelihood,5,1);
 326   DEBUGPRINT("df search2 finish\n");
 327   ctfforbayescheckIntensityOfMaxPosterior(ctf,out.posterior,out.numDistribution,&ctfinfo,intensity.MTFmode,intensity.Ainmode,1,1);
 328 */
 329   /*defocus refinement*/
 330 /*
 331   DEBUGPRINT("defocus refinement start\n");
 332   ctf.rmax = Rintermediate;
 333   ctf.rmin = Rmin;
 334   DEBUGPRINT2("set R range:%f~%f\n",ctf.rmin,ctf.rmax);
 335
 336   DEBUGPRINT("set maxposteriorflag of defocus,A,B,MTF\n");
 337   out.posterior[0].maxposteriorflag = 0;
 338   out.posterior[1].maxposteriorflag = 0;
 339   out.posterior[3].maxposteriorflag = 0;
 340 // out.posterior[5].maxposteriorflag = 0;
 341 // eosBayesMaxposteriorflagSet(&out,2);
 342   eosBayesMaxposteriorflagSet(&out,4);
 343 //  eosBayesMaxposteriorflagSet(&out,5);
 344
 345   DEBUGPRINT("set canstantflag of Ain");
 346   out.posterior[5].constantflag = 0;
 347
 348   DEBUGPRINT("set refmode\n");
 349   out.posterior[0].refmode = 1;
 350   out.posterior[1].refmode = 1;
 351   out.posterior[2].refmode = 1;
 352   out.posterior[3].refmode = 1;
 353   out.posterior[4].refmode = 0;
 354   out.posterior[5].refmode = 0;
 355
 356   eosBayesRefinementCondition(info.fptIn,&out,&llinfo,5,0);
 357   DistributionChangeCheck(&out, 6, 0);
 358   defocusLikelihoodInit(&likelihood,&out,1);
 359
 360   DefocusEstimationLoop(info,&out,&ctfinfo,&ctf,&intensity,&likelihood,6,1);
 361   DEBUGPRINT("defocus refinement finish\n");
 362   ctfforbayescheckIntensityOfMaxPosterior(ctf,out.posterior,out.numDistribution,&ctfinfo,intensity.MTFmode,intensity.Ainmode,1,1);
 363 */
 364   eosBayesPosteriorWrite(&out,0);
 365
 366   ctfforbayesFree(&ctf,&intensity,&likelihood,out.numLikelihood);
 367
 368   eosBayesFree(&out);
 369
 370   exit(EXIT_SUCCESS);
 371 }
 372
 373 void DefocusEstimationLoop(ctfDeterminationByBayesInfo info, eosBayes* out, ctfInfo* ctfinfo, ctfforbayes* ctf, intensityforbayes* intensity, defocuslikelihoodforbayes* likelihood, int count,int mode){
 374   int i;
 375
 376   defocusLikelihoodInit(likelihood,out,mode);
 377   DEBUGPRINT("In DefocusEstimationLoop\n");
 378   for(i=0;i<info.flagIn4;i++){
 379
 380         intensityRead(info.In4[i],intensity,i,0);
 381         DEBUGPRINT("ctf start\n");
 382         ctfFunctionforbayes(ctf,intensity,out,ctfinfo,likelihood,count,1);
 383
 384         DEBUGPRINT("estimate start\n");
 385         eosBayesEstimation(out,0);
 386         DistributionChangeCheck(out,count,1);
 387
 388         DEBUGPRINT("change start\n");
 389         eosBayesPosteriortoPrior(out,0);
 390
 391   }
 392
 393 }
 394
 395 void
 396 additionalUsage()
 397 {
 398   fprintf(stderr, "----- Additional Usage -----\n");
 399 }
 400