src/Tools/filter/moi2mrc/src/moi2mrc.c

   1 /*
   2 # %M% %Y% %I%
   3 # The latest update : %G% at %U%
   4 #
   5 #%Z% moi2mrc ver %I%
   6 #%Z% Created by
   7 #%Z%
   8 #%Z% Usage : moi2mrc
   9 #%Z% Attention
  10 #%Z%
  11 */
  12 static char __sccs_id[] = "%Z%moi2mrc ver%I%; Date:%D% %Z%";
  13 #include <stdio.h>
  14 #include <stdlib.h>
  15 #include <string.h>
  16 #include <math.h>
  17 #define GLOBAL_DECLARATION
  18 #include "../inc/config.h"
  19 #define DEBUG
  20 #include "genUtil.h"
  21 #include "Memory.h"
  22 #include "File.h"
  23 #include "eosString.h"
  24 #include "mrcImage.h"
  25
  26 typedef int            MedOpticsInteger;
  27 typedef float          MedOpticsReal;
  28 typedef short          MedOpticsImageInteger;
  29 #define MedOpticsCCDPixelSize (24)
  30
  31 typedef struct MedOpticsScan {
  32         MedOpticsInteger prescan;
  33         MedOpticsInteger underscan;
  34         MedOpticsInteger origin;
  35         MedOpticsInteger read;
  36         MedOpticsInteger postscan;
  37         MedOpticsInteger overscan;
  38         MedOpticsInteger bin;
  39 } MedOpticsScan;
  40
  41 typedef struct MedOpticsCameraParam {
  42         MedOpticsScan ser;
  43         MedOpticsScan par;
  44         MedOpticsInteger gain;
  45         MedOpticsInteger mpp;
  46         MedOpticsInteger openDelay;
  47         MedOpticsInteger closeDelay;
  48         MedOpticsInteger expTime;
  49         MedOpticsInteger clears;
  50         MedOpticsInteger cclear;
  51         MedOpticsInteger frame;
  52         MedOpticsInteger parShift;
  53         MedOpticsInteger numImages;
  54         MedOpticsInteger imDelay;
  55         MedOpticsInteger expInc;
  56         MedOpticsInteger crab;
  57         MedOpticsInteger offset;
  58         MedOpticsInteger setpointTemp;
  59         MedOpticsInteger spec3;
  60         MedOpticsInteger spec2;
  61         MedOpticsInteger spec1;
  62 } MedOpticsCameraParam;
  63
  64 typedef struct MedOpticsROI {
  65         MedOpticsInteger x;
  66         MedOpticsInteger y;
  67         MedOpticsInteger xsz;
  68         MedOpticsInteger ysz;
  69 } MedOpticsROI;
  70
  71 typedef struct MedOpticsLabels {
  72         char* name[4];
  73         char* comment[4];
  74         char* filter[4];
  75         char* lens[4];
  76         char* aperture[4];
  77         char* history[4];
  78 } MedOpticsLabels;
  79
  80 typedef struct MedOpticsHeader {
  81         char* title;
  82         char* vers;
  83         MedOpticsCameraParam params;
  84         MedOpticsInteger serDim;
  85         MedOpticsInteger parDim;
  86         MedOpticsInteger Shutter;
  87         char* startTime;
  88         char* endTime;
  89         MedOpticsInteger Temperature;
  90         char* type;
  91         MedOpticsInteger numMean;
  92         MedOpticsInteger exclude;
  93         MedOpticsLabels  labels;
  94         MedOpticsInteger min;
  95         MedOpticsInteger max;
  96         MedOpticsReal    mean;
  97         MedOpticsReal    stdev;
  98         MedOpticsROI     roi;
  99         MedOpticsInteger wmin;
 100         MedOpticsInteger wmax;
 101 } MedOpticsHeader;
 102
 103 typedef struct MedOpticsImage {
 104         MedOpticsHeader        Header;
 105         MedOpticsImageInteger* Image;
 106 } MedOpticsImage;
 107
 108 typedef struct lmoi2mrcInfo {
 109         mrcImageParaTypeReal ADUperElectron;
 110         mrcImageParaTypeReal Magnification;
 111 } lmoi2mrcInfo;
 112
 113 extern void MedOpticsImageFileRead(MedOpticsImage* moi, char* filename, int mode);
 114 extern void MedOpticsImageHeaderPrint(FILE* fpt, MedOpticsHeader* moi, int mode);
 115 extern void lmoi2mrc(MedOpticsImage* moi, mrcImage*  mrc, lmoi2mrcInfo linfo, int mode);
 116
 117 int
 118 main(int argc, char* argv[])
 119 {
 120         MedOpticsImage moi;
 121         mrcImage       mrc;
 122         moi2mrcInfo    info;
 123         lmoi2mrcInfo   linfo;
 124
 125         init0(&info);
 126     argCheck(&info, argc, argv);
 127     init1(&info);
 128
 129         linfo.ADUperElectron = info.ADUperE;
 130         linfo.Magnification  = info.Mag*info.Mag2*1000;
 131
 132         DEBUGPRINT("Program Start\n");
 133         MedOpticsImageFileRead(&moi, info.In, 0);
 134
 135         switch(info.Info) {
 136                 case 1: {
 137                         MedOpticsImageHeaderPrint(stdout, &(moi.Header), 0);
 138                         break;
 139                 }
 140                 default: {
 141                         break;
 142                 }
 143         }
 144
 145         lmoi2mrc(&moi, &mrc, linfo, info.mode);
 146
 147         if(info.flagOut) {
 148                 mrcFileWrite(&mrc, info.Out, "in main", 0);
 149         }
 150         exit(EXIT_SUCCESS);
 151 }
 152
 153 void
 154 additionalUsage()
 155 {
 156         fprintf(stderr, "----- Additional Usage -----\n");
 157         fprintf(stderr, ">>>> -info Integer\n");
 158         fprintf(stderr, " 1: All Header Info List \n");
 159 }
 160
 161
 162 void
 163 lmoi2mrc(MedOpticsImage* moi, mrcImage*  mrc, lmoi2mrcInfo linfo, int mode)
 164 {
 165         mrcDefaultHeaderValueSet(mrc);
 166         mrc->HeaderN.x = moi->Header.serDim;
 167         mrc->HeaderN.y = moi->Header.parDim;
 168         mrc->HeaderN.z = 1;
 169         mrc->HeaderMode = mode;
 170         mrc->HeaderLength.x = MedOpticsCCDPixelSize*10000/linfo.Magnification;
 171         mrc->HeaderLength.y = mrc->HeaderLength.x;
 172         mrc->HeaderLength.z = mrc->HeaderLength.x;
 173         mrcInit(mrc, NULL);
 174
 175         switch(mode) {
 176                 case mrcCharImage: {
 177                         int n = mrc->HeaderN.x*mrc->HeaderN.y;
 178                         int i;
 179                         int count;
 180                         double sum  = 0;
 181                         double sum2 = 0;
 182                         double mean;
 183                         double stdev;
 184                         double stdev3;
 185
 186                         fprintf(stderr, "mean+-3 sd\n");
 187                         count = 0;
 188                         for(i=0; i<n; i++) {
 189                                 if(0<=moi->Image[i]) {
 190                                         sum  +=     ((double)moi->Image[i]);
 191                                         sum2 +=  SQR((double)moi->Image[i]);
 192                                         count++;
 193                                 }
 194                         }
 195                         mean  = sum/count;
 196                         stdev = sqrt(sum2/count - mean*mean);
 197                         stdev3 = stdev*3;
 198
 199                         for(i=0; i<n; i++) {
 200                                 if(moi->Image[i] <= mean - stdev3) {
 201                                         mrc->ImageCharImage[i] = 0;
 202                                 } else if(moi->Image[i] >= mean + stdev3) {
 203                                         mrc->ImageCharImage[i] = 255;
 204                                 } else {
 205                                         mrc->ImageCharImage[i] =  (((double)moi->Image[i]) - mean)/stdev3*127 + 127;
 206                                 }
 207                         }
 208                         break;
 209                 }
 210                 case mrcShortImage: {
 211                         int i;
 212                         int n = mrc->HeaderN.x*mrc->HeaderN.y;
 213
 214                         fprintf(stderr, "Direct(Positive only)\n");
 215                         for(i=0; i<n; i++) {
 216                                 if(moi->Image[i] < 0) {
 217                                         mrc->ImageShortImage[i] = 0;
 218                                 } else {
 219                                         mrc->ImageShortImage[i] = moi->Image[i];
 220                                 }
 221 #ifdef DEBUG
 222                                 if(i%2048==1024) {
 223                                                 fprintf(stderr, "%d %d\n", mrc->ImageShortImage[i], moi->Image[i]);
 224                                 }
 225 #endif
 226                         }
 227                         break;
 228                 }
 229                 case mrcFloatImage: {
 230                         int n = mrc->HeaderN.x*mrc->HeaderN.y;
 231                         int i;
 232
 233                         fprintf(stderr, "Electon Number\n");
 234                         for(i=0; i<n; i++) {
 235                                 if(moi->Image[i] < 0) {
 236                                         mrc->ImageFloatImage[i] = 0;
 237                                 } else {
 238                                         mrc->ImageFloatImage[i] = ((double)moi->Image[i])/linfo.ADUperElectron;
 239                                 }
 240                         }
 241                         break;
 242                 }
 243                 default: {
 244                         fprintf(stderr, "Not supported  mode in lmoi2mrc: %d\n", mode);
 245                         exit(EXIT_FAILURE);
 246                 }
 247         }
 248 }
 249
 250 void
 251 MedOpticsImageFileReadScan(MedOpticsScan* scan, FILE* fpt, int mode)
 252 {
 253         fread(scan, sizeof(MedOpticsInteger), 7, fpt);
 254 #ifdef M68
 255         memoryByteSwap(scan, sizeof(MedOpticsInteger), 7);
 256 #endif
 257 }
 258
 259 void
 260 MedOpticsImageScanPrint(FILE* fpt, MedOpticsScan* scan, int mode)
 261 {
 262         fprintf(fpt, "Scan ----\n");
 263         fprintf(fpt, "prescan:   %d\n", scan->prescan);
 264         fprintf(fpt, "underscan: %d\n", scan->underscan);
 265         fprintf(fpt, "origin:    %d\n", scan->origin);
 266         fprintf(fpt, "read:      %d\n", scan->read);
 267         fprintf(fpt, "postscan:  %d\n", scan->postscan);
 268         fprintf(fpt, "overscan:  %d\n", scan->overscan);
 269         fprintf(fpt, "bin:       %d\n", scan->bin);
 270 }
 271
 272 void
 273 MedOpticsImageFileReadCameraParam(MedOpticsCameraParam* cam, FILE* fpt, int mode)
 274 {
 275         MedOpticsImageFileReadScan(&(cam->ser), fpt, mode);
 276         MedOpticsImageFileReadScan(&(cam->par), fpt, mode);
 277         fread(&(cam->gain), sizeof(MedOpticsInteger), 18, fpt);
 278 #ifdef M68
 279         memoryByteSwap(&(cam->gain), sizeof(MedOpticsInteger), 18);
 280 #endif
 281 }
 282
 283 void
 284 MedOpticsImageCameraParamPrint(FILE* fpt, MedOpticsCameraParam* cam, int mode)
 285 {
 286         fprintf(fpt, "Camera Parameters ----\n");
 287         fprintf(fpt, "---- Serial ");
 288         MedOpticsImageScanPrint(fpt, &(cam->ser), mode);
 289         fprintf(fpt, "---- Parallel ");
 290         MedOpticsImageScanPrint(fpt, &(cam->par), mode);
 291         fprintf(fpt, "gain:       %d\n", cam->gain);
 292         fprintf(fpt, "mpp:        %d\n", cam->mpp);
 293         fprintf(fpt, "openDelay:  %d\n", cam->openDelay);
 294         fprintf(fpt, "closeDelay: %d\n", cam->closeDelay);
 295         fprintf(fpt, "expTime:    %d\n", cam->expTime);
 296         fprintf(fpt, "clears:     %d\n", cam->clears);
 297         fprintf(fpt, "cclear:     %d\n", cam->cclear);
 298         fprintf(fpt, "frame:      %d\n", cam->frame);
 299         fprintf(fpt, "parShift:   %d\n", cam->parShift);
 300         fprintf(fpt, "numImages:  %d\n", cam->numImages);
 301         fprintf(fpt, "imDelay:    %d\n", cam->imDelay);
 302         fprintf(fpt, "expInc:     %d\n", cam->expInc);
 303         fprintf(fpt, "crab:       %d\n", cam->crab);
 304         fprintf(fpt, "offset:     %d\n", cam->offset);
 305         fprintf(fpt, "setTemp:    %d\n", cam->setpointTemp);
 306         fprintf(fpt, "spec3:      %d\n", cam->spec3);
 307         fprintf(fpt, "spec2:      %d\n", cam->spec2);
 308         fprintf(fpt, "spec1:      %d\n", cam->spec1);
 309 }
 310
 311 void
 312 MedOpticsImageFileReadROI(MedOpticsROI* roi, FILE* fpt, int mode)
 313 {
 314         fread(&(roi->x), sizeof(MedOpticsInteger), 4, fpt);
 315 #ifdef M68
 316         memoryByteSwap(&(roi->x), sizeof(MedOpticsInteger), 4);
 317 #endif
 318 }
 319
 320 void
 321 MedOpticsImageROIPrint(FILE* fpt, MedOpticsROI* roi,  int mode)
 322 {
 323         fprintf(fpt, "ROI ----\n");
 324         fprintf(fpt, "x:      %d\n", roi->x);
 325         fprintf(fpt, "y:      %d\n", roi->y);
 326         fprintf(fpt, "xsz:    %d\n", roi->xsz);
 327         fprintf(fpt, "ysz:    %d\n", roi->ysz);
 328 }
 329
 330 void
 331 MedOpticsImageFileReadLabels(MedOpticsLabels* labels, FILE* fpt, int mode)
 332 {
 333         char s[1024];
 334         int i;
 335
 336         for(i=0; i<4; i++) {
 337                 stringGetFromFile(s, "MedOpticsImageFileRead", fpt, stdout, 3); labels->name[i] = stringDuplicate(s);
 338         }
 339         for(i=0; i<4; i++) {
 340                 stringGetFromFile(s, "MedOpticsImageFileRead", fpt, stdout, 3); labels->comment[i] = stringDuplicate(s);
 341         }
 342         for(i=0; i<4; i++) {
 343                 stringGetFromFile(s, "MedOpticsImageFileRead", fpt, stdout, 3); labels->filter[i] = stringDuplicate(s);
 344         }
 345         for(i=0; i<4; i++) {
 346                 stringGetFromFile(s, "MedOpticsImageFileRead", fpt, stdout, 3); labels->lens[i] = stringDuplicate(s);
 347         }
 348         for(i=0; i<4; i++) {
 349                 stringGetFromFile(s, "MedOpticsImageFileRead", fpt, stdout, 3); labels->aperture[i] = stringDuplicate(s);
 350         }
 351         for(i=0; i<4; i++) {
 352                 stringGetFromFile(s, "MedOpticsImageFileRead", fpt, stdout, 3); labels->history[i] = stringDuplicate(s);
 353         }
 354 }
 355
 356 void
 357 MedOpticsImageLabelsPrint(FILE* fpt, MedOpticsLabels* labels, int mode)
 358 {
 359         int i;
 360
 361         fprintf(fpt, "Labels ----\n");
 362         fprintf(fpt, "---- Name ----\n");
 363         for(i=0; i<4; i++) {
 364                 fprintf(fpt, "%s\n", labels->name[i]);
 365         }
 366         fprintf(fpt, "---- Comment ----\n");
 367         for(i=0; i<4; i++) {
 368                 fprintf(fpt, "%s\n", labels->comment[i]);
 369         }
 370         fprintf(fpt, "---- Filter ----\n");
 371         for(i=0; i<4; i++) {
 372                 fprintf(fpt, "%s\n", labels->filter[i]);
 373         }
 374         fprintf(fpt, "---- Lens ----\n");
 375         for(i=0; i<4; i++) {
 376                 fprintf(fpt, "%s\n", labels->lens[i]);
 377         }
 378         fprintf(fpt, "---- Aperture ----\n");
 379         for(i=0; i<4; i++) {
 380                 fprintf(fpt, "%s\n", labels->aperture[i]);
 381         }
 382         fprintf(fpt, "---- History ----\n");
 383         for(i=0; i<4; i++) {
 384                 fprintf(fpt, "%s\n", labels->history[i]);
 385         }
 386 }
 387
 388 void
 389 MedOpticsImageFileRead(MedOpticsImage* MOI, char* filename, int mode)
 390 {
 391         char s[1024];
 392         FILE* fpt;
 393         MedOpticsHeader* moi;
 394         int n;
 395
 396         fpt = fileOpen(filename, "r");
 397
 398         moi = &(MOI->Header);
 399         stringGetFromFile(s, "MedOpticsImageFileRead", fpt, stdout, 3); moi->title = stringDuplicate(s);
 400         stringGetFromFile(s, "MedOpticsImageFileRead", fpt, stdout, 3); moi->vers  = stringDuplicate(s);
 401         MedOpticsImageFileReadCameraParam(&(moi->params), fpt, mode);
 402         fread(&(moi->serDim), sizeof(MedOpticsInteger), 3, fpt);
 403 #ifdef M68
 404         memoryByteSwap(&(moi->serDim), sizeof(MedOpticsInteger), 3);
 405 #endif
 406         stringGetFromFile(s, "MedOpticsImageFileRead", fpt, stdout, 3); moi->startTime = stringDuplicate(s);
 407         stringGetFromFile(s, "MedOpticsImageFileRead", fpt, stdout, 3); moi->endTime   = stringDuplicate(s);
 408         fread(&(moi->Temperature), sizeof(MedOpticsInteger), 1, fpt);
 409 #ifdef M68
 410         memoryByteSwap(&(moi->Temperature), sizeof(MedOpticsInteger), 1);
 411 #endif
 412         stringGetFromFile(s, "MedOpticsImageFileRead", fpt, stdout, 3); moi->type = stringDuplicate(s);
 413         fread(&(moi->numMean), sizeof(MedOpticsInteger), 2, fpt);
 414 #ifdef M68
 415         memoryByteSwap(&(moi->numMean), sizeof(MedOpticsInteger), 2);
 416 #endif
 417         MedOpticsImageFileReadLabels(&(moi->labels), fpt, mode);
 418         fread(&(moi->min), sizeof(MedOpticsInteger), 4, fpt);
 419 #ifdef M68
 420         memoryByteSwap(&(moi->min), sizeof(MedOpticsInteger), 4);
 421 #endif
 422         MedOpticsImageFileReadROI(&(moi->roi), fpt, mode);
 423         fread(&(moi->wmin), sizeof(MedOpticsInteger), 2, fpt);
 424 #ifdef M68
 425         memoryByteSwap(&(moi->wmin), sizeof(MedOpticsInteger), 2);
 426 #endif
 427         MOI->Image = (MedOpticsImageInteger*)memoryAllocate(sizeof(MedOpticsImageInteger)*moi->serDim*moi->parDim, "MedOptics");
 428         n = fread(MOI->Image, sizeof(MedOpticsImageInteger), moi->serDim*moi->parDim, fpt);
 429         if(n != moi->serDim*moi->parDim) {
 430                 fprintf(stderr, "All Data Cannot Read: %d(present) != %d(predict) at %ld\n", n, moi->serDim*moi->parDim, ftell(fpt));
 431         }
 432 #ifdef M68
 433         memoryByteSwap(MOI->Image, sizeof(MedOpticsImageInteger), moi->serDim*moi->parDim);
 434 #endif
 435
 436 }
 437
 438 void
 439 MedOpticsImageHeaderPrint(FILE* fpt, MedOpticsHeader* moi, int mode)
 440 {
 441         fprintf(fpt, "title: %s\n", moi->title);
 442         fprintf(fpt, "vers:  %s\n", moi->vers);
 443         fprintf(fpt, "---- Camera Params ");
 444         MedOpticsImageCameraParamPrint(fpt, &(moi->params), mode);
 445         fprintf(fpt, "serDim:    %d\n", moi->serDim);
 446         fprintf(fpt, "parDim:    %d\n", moi->parDim);
 447         fprintf(fpt, "Shutter:   %d\n", moi->Shutter);
 448         fprintf(fpt, "startTime:  %s\n", moi->startTime);
 449         fprintf(fpt, "endTime:  %s\n", moi->endTime);
 450         fprintf(fpt, "Temperature: %d\n", moi->Temperature);
 451         fprintf(fpt, "type:  %s\n", moi->type);
 452         fprintf(fpt, "numMean:  %d\n", moi->numMean);
 453         fprintf(fpt, "exclude:  %d\n", moi->exclude);
 454         MedOpticsImageLabelsPrint(fpt, &(moi->labels), mode);
 455         fprintf(fpt, "min:  %d\n", moi->min);
 456         fprintf(fpt, "max:  %d\n", moi->max);
 457         fprintf(fpt, "mean:  %15.6f\n", moi->mean);
 458         fprintf(fpt, "stdev:  %15.6f\n", moi->stdev);
 459         MedOpticsImageROIPrint(fpt, &(moi->roi), mode);
 460         fprintf(fpt, "wmin:  %d\n", moi->wmin);
 461         fprintf(fpt, "wmax:  %d\n", moi->wmax);
 462 }
 463