1 struct Foo
2 {
3 double length() const;
4 };
5
6 std::vector<Foo> foos;
7
8 // now we want to remove all foos being shorter than 5.0
9 std::remove_if(foos.begin(), foos.end(), boost::bind( std::less_equal<double>(), boost::bind(&Foo::length, _1), 5.0));
10
11 // well, all this in a single line. it works but what does it do?
12 // 1) first a functor is declared that takes a single double and returns a boolean value.
13 // this is done by binding the value 5.0 to the second parameter of the function
14 // std::less_equal<double>(double lhs, double rhs)
15 boost::function< bool(double) > le_five = boost::bind( std::less_equal<double>(), _1, 5.0);
16
17 // 2) the second two binds, bind the return value of Foo::length to the first parameter of
18 // the functor created by the first bind operation.
19 std::remove_if(foo.begin(), foo.end(), boost::bind(le_five, boost::bind(&Foo::length, _1))); |
1 typedef float ValueType;
2 typedef itk::Image<ValueType, 2> Image2D;
3 typedef itk::Statistics::ScalarImageToHistogramGenerator<Image2D> HistogramGeneratorType;
4
5 const unsigned int bins = 64;
6
7 HistogramGeneratorType::RealPixelType low = -0.5/(bins-1);
8 HistogramGeneratorType::RealPixelType high = (bins-0.5)/(bins-1);
9
10 HistogramGeneratorType::Pointer hist_gen = HistogramGeneratorType::New();
11
12 hist_gen->SetNumberOfBins(bins);
13 hist_gen->SetHistogramMin(low);
14 hist_gen->SetHistogramMax(high);
15 hist_gen->SetInput(frangi_maxima);
16 hist_gen->Compute();
17
18 const HistogramGeneratorType::HistogramType* histo = hist_gen->GetOutput(); |
1 HistogramGeneratorType::HistogramType::FrequencyType sum = 0;
2
3 HistogramGeneratorType::HistogramType::FrequencyType element_count = histo->GetTotalFrequency();
4 HistogramGeneratorType::HistogramType::FrequencyType upper = 0.7*element_count;
5
6 HistogramGeneratorType::HistogramType::ConstIterator hit = histo->Begin();
7 HistogramGeneratorType::HistogramType::ConstIterator hit_end = histo->End();
8 for (; hit != hit_end; ++hit)
9 {
10 sum += hit.GetFrequency(); // access the number of elements in the current bin
11 if (sum > upper) break;
12 }
13
14 const unsigned int bin = hit.GetInstanceIdentifier();
15 const HistogramGeneratorType::RealPixelType bin_value = static_cast<HistogramGeneratorType::RealPixelType>(instance+1)/64; |
1 function deconvolution()
2 f1 = [0 -1 1];
3
4 I = rand(512, 512);
5
6 F = convmtx2(f1, size(I));
7 R = reshape(F*I(:),size(f1)+size(I)-1);
8
9 tic
10 I_rec = reshape(F\R(:), size(I));
11 toc
12
13 sum(sum(I_rec-I))
14 end |
1 function deconvolution( image, method, epsilon )
2 %% argument checks...
3 if nargin~=3
4 epsilon = eps('double');
5 end
6
7 %% data setup
8 if strcmp(method,'disk')
9 % symmteric 9-by-9 filter
10 hs = fspecial('disk',4);
11 elseif strcmp(method, 'laplacian')
12 % 3x3 laplacian
13 hs = fspecial('laplacian', 0.0);
14 elseif strcmp(method, 'assym_gradient1')
15 % assymmetric filter
16 f = [0 -1 1];
17 hs = f'*f;
18 elseif strcmp(method, 'assym_gradient2')
19 f = [0 -1 1];
20 hs = padarray(f, [(size(f',1)-size(f,1))/2 0]) + ...
21 padarray(f', [0 (size(f,2)-size(f',2))/2]);
22 end
23
24 % convert the image to double for fftn (im2double would also work, but
25 % i don't want the normalization to [0;1])...
26 cam = double(image);
27 [h w] = size(cam);
28
29 % reconstruction will always be done based on twice the input image
30 % size - otherwise it does simply not work (proven empirically, no idea
31 % why? is it due to shannon's theorem???)
32 filter_size = 2.*size(cam);
33
34 hf = fftn(hs, filter_size);
35 cam_filtered = real(ifft2(hf.*fftn(cam, filter_size)));
36
37 % we consider two cased
38 % a) trim the output of hs*I to the resolution of I
39 % b) keep the output of hs*I at the resolution of 2*size(I)
40 figure('Name','Input Data','NumberTitle','off')
41
42 subplot(1,3,1);
43 imagesc(cam);
44 axis equal tight;
45 xlabel('input image')
46
47 subplot(1,3,2);
48 imagesc(cam_filtered(1:h,1:w));
49 axis equal tight;
50 xlabel('filtered')
51
52 subplot(1,3,3);
53 imagesc(cam_filtered);
54 axis equal tight;
55 xlabel('filtered (double resolution)')
56
57 colormap gray;
58
59 %% computation of the (pseudo-)inverse filter
60 % computation of the inverse filter taken from the brilliant blog of
61 % Steve - check it out!
62 % http://blogs.mathworks.com/steve/2007/08/13/image-deblurring-introduction/
63
64 hf_inv = conj(hf)./(abs(hf).^2 + epsilon);
65
66 %% full resolution - full information
67 figure('Name','Double Resolution Reconstruction','NumberTitle','off');
68
69 cam_pinv = real(ifft2(fftn(cam_filtered).*hf_inv));
70
71 subplot(1,2,1);
72
73 % compute the mean and variance of the difference image. they are
74 % expressed as a fraction of the maximal intensity value of the input
75 % image, so we can cope with different types of image normalizations
76 diff = abs(cam_pinv(1:h,1:w)-cam)/max(cam(:));
77 diff_mean = mean(diff(:))*100;
78 diff_var = var(diff(:))*100;
79 label = sprintf('absolute difference image\n(mean: %2.2f %%, var: %2.2f %%)', diff_mean, diff_var);
80
81 imagesc(diff);
82 axis equal tight;
83 xlabel(label)
84
85 subplot(1,2,2);
86 imagesc(cam_pinv(1:h,1:w));
87 axis equal tight;
88 label = sprintf('reconstruction\n(trimmed to input resolution)');
89 xlabel(label)
90
91 colormap gray;
92
93 %% trimmed to original image resolution
94 figure('Name','Single Resolution Reconstruction','NumberTitle','off')
95
96 cam_filtered = cam_filtered(1:h, 1:w);
97 cam_pinv2 = real(ifft2(fftn(cam_filtered, filter_size).*hf_inv));
98
99 subplot(1,2,1);
100
101 % see the comment from above...
102 diff = abs(cam_pinv2(1:h,1:w)-cam)/max(cam(:));
103 diff_mean = mean(diff(:))*100;
104 diff_var = var(diff(:))*100;
105 label = sprintf('absolute difference image\n(mean: %2.2f %%, var: %2.2f %%)', diff_mean, diff_var);
106
107 imagesc(diff);
108 axis equal tight;
109 xlabel(label)
110
111 subplot(1,2,2);
112 imagesc(cam_pinv2(1:h,1:w));
113 axis equal tight;
114 label = sprintf('reconstruction\n(trimmed to input resolution)');
115 xlabel(label)
116
117 colormap gray;
118
119 end |
