GCC Code Coverage Report
Directory: ./ Exec Total Coverage
File: src/unrrdu/dering.c Lines: 29 81 35.8 %
Date: 2017-05-26 Branches: 1 50 2.0 %

Line Branch Exec Source
1
/*
2
  Teem: Tools to process and visualize scientific data and images             .
3
  Copyright (C) 2013, 2012, 2011, 2010, 2009  University of Chicago
4
  Copyright (C) 2008, 2007, 2006, 2005  Gordon Kindlmann
5
  Copyright (C) 2004, 2003, 2002, 2001, 2000, 1999, 1998  University of Utah
6
7
  This library is free software; you can redistribute it and/or
8
  modify it under the terms of the GNU Lesser General Public License
9
  (LGPL) as published by the Free Software Foundation; either
10
  version 2.1 of the License, or (at your option) any later version.
11
  The terms of redistributing and/or modifying this software also
12
  include exceptions to the LGPL that facilitate static linking.
13
14
  This library is distributed in the hope that it will be useful,
15
  but WITHOUT ANY WARRANTY; without even the implied warranty of
16
  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
17
  Lesser General Public License for more details.
18
19
  You should have received a copy of the GNU Lesser General Public License
20
  along with this library; if not, write to Free Software Foundation, Inc.,
21
  51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
22
*/
23
24
#include "unrrdu.h"
25
#include "privateUnrrdu.h"
26
27
#define INFO "Ring removal for CT"
28
static const char *_unrrdu_deringInfoL =
29
(INFO
30
 ". Should be considered a work-in-progress. ");
31
32
int
33
unrrdu_deringMain(int argc, const char **argv, const char *me,
34
                  hestParm *hparm) {
35
2
  hestOpt *opt = NULL;
36
1
  char *out, *err;
37
1
  Nrrd *nin, *nout;
38
  airArray *mop;
39
  int pret;
40
41
1
  Nrrd *nmask;
42
  NrrdDeringContext *drc;
43
1
  double center[2], radScale, clampPerc[2], backval;
44
1
  int verbose, linterp, vertSeam;
45
1
  unsigned int thetaNum;
46
1
  NrrdKernelSpec *rkspec, *tkspec;
47
48
  /* HEY is this needed? (to display -rk and -tk kernels) */
49
1
  hparm->elideSingleOtherDefault = AIR_FALSE;
50
51
1
  hestOptAdd(&opt, "c,center", "x y", airTypeDouble, 2, 2, center, NULL,
52
             "center of rings, in index space of fastest two axes");
53
1
  hestOptAdd(&opt, "v,verbose", "v", airTypeInt, 1, 1, &verbose, "0",
54
             "verbosity level");
55
1
  hestOptAdd(&opt, "li,linterp", "bool", airTypeBool, 1, 1, &linterp, "false",
56
             "whether to use linear interpolation during polar transform");
57
1
  hestOptAdd(&opt, "vs,vertseam", "bool", airTypeBool, 1, 1, &vertSeam, "false",
58
             "whether to dering left and right sides separately "
59
             "(requires an even value for -tn thetanum)");
60
1
  hestOptAdd(&opt, "tn,thetanum", "# smpls", airTypeUInt, 1, 1, &thetaNum,
61
             "20", "# of theta samples");
62
1
  hestOptAdd(&opt, "rs,radscale", "scale", airTypeDouble, 1, 1, &radScale,
63
             "1.0", "scaling on radius in polar transform");
64
1
  hestOptAdd(&opt, "rk,radiuskernel", "kern", airTypeOther, 1, 1, &rkspec,
65
             "gauss:3,4",
66
             "kernel for high-pass filtering along radial direction",
67
1
             NULL, NULL, nrrdHestKernelSpec);
68
1
  hestOptAdd(&opt, "tk,thetakernel", "kern", airTypeOther, 1, 1, &tkspec,
69
             "box",
70
             "kernel for blurring along theta direction.",
71
1
             NULL, NULL, nrrdHestKernelSpec);
72
1
  hestOptAdd(&opt, "cp,clampperc", "lo hi", airTypeDouble, 2, 2, clampPerc,
73
             "0.0 0.0",
74
             "when clamping values as part of ring estimation, the "
75
             "clamping range is set to exclude this percent of values "
76
             "from the low and high end of the data range");
77
1
  hestOptAdd(&opt, "m,mask", "mask", airTypeOther, 1, 1, &nmask, "",
78
             "optional: after deringing, output undergoes a lerp, "
79
             "parameterized by this array, from the background value "
80
             "(via \"-b\") where mask=0 to the original deringing "
81
             "output where mask=1.  This lerp is effectively the same "
82
             "as a \"unu 3op lerp\", so this should either be match the "
83
             "input in size, or match its slices along the slowest axis.",
84
1
             NULL, NULL, nrrdHestNrrd);
85
1
  hestOptAdd(&opt, "b,back", "val", airTypeDouble, 1, 1, &backval, "0.0",
86
             "when using a mask (\"-m\"), the background value to "
87
             "lerp with.");
88
1
  OPT_ADD_NIN(nin, "input nrrd");
89
1
  OPT_ADD_NOUT(out, "output nrrd");
90
91
1
  mop = airMopNew();
92
1
  airMopAdd(mop, opt, (airMopper)hestOptFree, airMopAlways);
93
94
2
  USAGE(_unrrdu_deringInfoL);
95
  PARSE();
96
  airMopAdd(mop, opt, (airMopper)hestParseFree, airMopAlways);
97
98
  nout = nrrdNew();
99
  airMopAdd(mop, nout, (airMopper)nrrdNuke, airMopAlways);
100
101
  if (nmask) {
102
    if (!(2 == nmask->dim
103
          && nrrdTypeBlock != nmask->type
104
          && nmask->axis[0].size == nin->axis[0].size
105
          && nmask->axis[1].size == nin->axis[1].size)) {
106
      fprintf(stderr, "%s: given mask not 2-D %u-by-%u array of scalar type",
107
              me, AIR_CAST(unsigned int, nin->axis[0].size),
108
              AIR_CAST(unsigned int, nin->axis[1].size));
109
      airMopError(mop);
110
      return 1;
111
    }
112
  }
113
114
  drc = nrrdDeringContextNew();
115
  airMopAdd(mop, drc, (airMopper)nrrdDeringContextNix, airMopAlways);
116
  if (nrrdDeringVerboseSet(drc, verbose)
117
      || nrrdDeringLinearInterpSet(drc, linterp)
118
      || nrrdDeringVerticalSeamSet(drc, vertSeam)
119
      || nrrdDeringInputSet(drc, nin)
120
      || nrrdDeringCenterSet(drc, center[0], center[1])
121
      || nrrdDeringRadiusScaleSet(drc, radScale)
122
      || nrrdDeringThetaNumSet(drc, thetaNum)
123
      || nrrdDeringRadialKernelSet(drc, rkspec->kernel, rkspec->parm)
124
      || nrrdDeringThetaKernelSet(drc, tkspec->kernel, tkspec->parm)
125
      || nrrdDeringClampPercSet(drc, clampPerc[0], clampPerc[1])
126
      || nrrdDeringExecute(drc, nout)) {
127
    airMopAdd(mop, err = biffGetDone(NRRD), airFree, airMopAlways);
128
    fprintf(stderr, "%s: error deringing:\n%s", me, err);
129
    airMopError(mop);
130
    return 1;
131
  }
132
133
  if (nmask) {
134
    NrrdIter *nitout, *nitmask, *nitback;
135
    Nrrd *ntmp;
136
    nitout = nrrdIterNew();
137
    airMopAdd(mop, nitout, (airMopper)nrrdIterNix, airMopAlways);
138
    nitmask = nrrdIterNew();
139
    airMopAdd(mop, nitmask, (airMopper)nrrdIterNix, airMopAlways);
140
    nitback = nrrdIterNew();
141
    airMopAdd(mop, nitback, (airMopper)nrrdIterNix, airMopAlways);
142
    nrrdIterSetValue(nitback, backval);
143
144
    ntmp = nrrdNew();
145
    airMopAdd(mop, ntmp, (airMopper)nrrdNuke, airMopAlways);
146
147
    nrrdIterSetNrrd(nitout, nout);
148
    nrrdIterSetNrrd(nitmask, nmask);
149
    if (nrrdArithIterTernaryOpSelect(ntmp, nrrdTernaryOpLerp,
150
                                     nitmask, nitback, nitout, 2)
151
        || nrrdCopy(nout, ntmp)) {
152
      airMopAdd(mop, err = biffGetDone(NRRD), airFree, airMopAlways);
153
      fprintf(stderr, "%s: error masking:\n%s", me, err);
154
      airMopError(mop);
155
      return 1;
156
    }
157
  }
158
159
  SAVE(out, nout, NULL);
160
161
  airMopOkay(mop);
162
  return 0;
163
1
}
164
165
UNRRDU_CMD_HIDE(dering, INFO);