cudaconfidencemapsdemo.cpp 20.8 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
// ================================================================================================
// 
// This file is part of the CAMPVis Software Framework.
// 
// If not explicitly stated otherwise: Copyright (C) 2012-2014, all rights reserved,
//      Christian Schulte zu Berge <christian.szb@in.tum.de>
//      Chair for Computer Aided Medical Procedures
//      Technische Universitaet Muenchen
//      Boltzmannstr. 3, 85748 Garching b. Muenchen, Germany
// 
// For a full list of authors and contributors, please refer to the file "AUTHORS.txt".
// 
// Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file 
// except in compliance with the License. You may obtain a copy of the License at
// 
// http://www.apache.org/licenses/LICENSE-2.0
// 
// Unless required by applicable law or agreed to in writing, software distributed under the 
// License is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, 
// either express or implied. See the License for the specific language governing permissions 
// and limitations under the License.
// 
// ================================================================================================

#include "cudaconfidencemapsdemo.h"

27
28
29
30
#include <iomanip>
#include <sstream>
#include <QDateTime>
#include <QDir>
31
32
#include <QClipboard>
#include <QApplication>
33
34
35
36
37
38

#ifdef CAMPVIS_HAS_MODULE_DEVIL
#include <IL/il.h>
#include <IL/ilu.h>
#endif

39
#include "core/datastructures/imagedata.h"
40
#include "core/datastructures/imagerepresentationgl.h"
41
42
#include "core/classification/geometry1dtransferfunction.h"
#include "core/classification/tfgeometry1d.h"
43
#include "cgt/event/keyevent.h"
44
45
46

namespace campvis {

47
48
49
50
    CudaConfidenceMapsDemo::CudaConfidenceMapsDemo(DataContainer* dc)
        : AutoEvaluationPipeline(dc)
        , _usIgtlReader()
        , _usCropFilter(&_canvasSize)
51
        , _usBlurFilter(&_canvasSize)
52
        , _usResampler(&_canvasSize)
53
        , _usMapsSolver()
54
        , _usFusion(&_canvasSize)
55
        , _usFanRenderer(&_canvasSize)
56
57
58
        , p_useFixedIterationCount("UseFixedIterationCount", "Use Fixed Iteration Count", false)
        , p_millisecondBudget("MillisecondBudget", "(P)CG Milliseconds per frame", 24.0f, 10.0f, 60.0f)
        , p_iterationBudget("IterationBudget", "(P)CG Iteration Count", 100, 0, 1000)
59
        , p_connectDisconnectButton("ConnectToIGTLink", "Connect/Disconnect")
Declara Denis's avatar
Declara Denis committed
60
61
        , p_resamplingScale("ResampleScale", "Resample Scale", 0.5f, 0.01f, 1.0f)
        , p_beta("Beta", "Beta", 80.0f, 1.0f, 200.0f)
62
63
        , p_collectStatistics("CollectStatistics", "Collect Statistics", false)
        , p_copyStatisticsToClipboard("CopyStatisticsToClipboard", "Copy Statistics To Clipboard as CSV")
64
65
66
        , p_showAdvancedOptions("ShowAdvancedOptions", "Advanced options...", false)
        , p_useAlphaBetaFilter("UseAlphaBetaFilter", "Alpha-Beta-Filter", true)
        , p_gaussianFilterSize("GaussianSigma", "Blur amount", 2.5f, 1.0f, 10.0f)
Declara Denis's avatar
Declara Denis committed
67
        , p_gradientScaling("GradientScaling", "Scaling factor for gradients", 2.0f, 0.001, 10)
68
69
        , p_alpha("Alpha", "Alpha", 2.0f, 0.0f, 10.0f)
        , p_gamma("Gamma", "Gamma", 0.03f, 0.0f, 0.4f, 0.001, 4)
70
        , p_fanHalfAngle("FanHalfAngle", "Fan Half Angle", 28.0f, 1.0f, 90.0f)
71
        , p_fanInnerRadius("FanInnerRadius", "Fan Inner Radius", 0.222f, 0.001f, 0.999f)
72
        , p_useSpacingEncodedFanGeometry("UseSpacingEncodedFanGeomtry", "Use spacing encoded fan geometry", true)
73
74
75
        , p_recordingDirectory("RecordingDirectory", "Recording output direcotry", "D:\\us_acquisitions\\")
        , p_enableRecording("EnableRecording", "Enable recording", false)
        , _recordedFrames(0)
Declara Denis's avatar
Declara Denis committed
76
        , _statisticsLastUpdateTime()
77
    {
78
79
80
        // Calculate file prefix, using date, hour, minute and second of when the pipeline was created
        _filePrefix = QDateTime::currentDateTime().toString("yyyyMMdd_hhmmss_").toStdString();

81
82
        addProcessor(&_usIgtlReader);
        addProcessor(&_usCropFilter);
83
        addProcessor(&_usBlurFilter);
84
        addProcessor(&_usResampler);
85
        addProcessor(&_usMapsSolver);
86
        addProcessor(&_usFusion);
87
        addProcessor(&_usFanRenderer);
88

89
        addProperty(p_useFixedIterationCount);
90
        addProperty(p_millisecondBudget);
91
        addProperty(p_iterationBudget);
92
        addProperty(p_connectDisconnectButton);
93
        addProperty(p_resamplingScale);
94
        addProperty(p_beta);
95
96
        addProperty(p_collectStatistics);
        addProperty(p_copyStatisticsToClipboard);
97
98
99
100
101

        addProperty(p_showAdvancedOptions);

        addProperty(p_useAlphaBetaFilter);
        addProperty(p_gaussianFilterSize);
Declara Denis's avatar
Declara Denis committed
102
        addProperty(p_gradientScaling);
103
104
105
106
        addProperty(p_alpha);
        addProperty(p_gamma);
        addProperty(p_fanHalfAngle);
        addProperty(p_fanInnerRadius);
107
108
109
        addProperty(p_recordingDirectory);
        addProperty(p_enableRecording);

110
111
        addProperty(p_useSpacingEncodedFanGeometry);

112
        setAdvancedPropertiesVisibility(false);
113

114
115
116
        _canvasSize.setVisible(false);
        _renderTargetID.setVisible(false);

117
118
119
        // Reserve memory for statistics, so that (hopefully) no reallocation happens at runtime
        _statistics.reserve(1000);
        _objectCreationTime = tbb::tick_count::now();
120
121
122
123
124
125
126
127
    }

    CudaConfidenceMapsDemo::~CudaConfidenceMapsDemo() {
    }

    void CudaConfidenceMapsDemo::init() {
        AutoEvaluationPipeline::init();

128
        // Set intial options
129
130
131
        _usIgtlReader.p_receiveImages.setValue(true);
        _usIgtlReader.p_receiveTransforms.setValue(false);
        _usIgtlReader.p_receivePositions.setValue(false);
132
        _usResampler.p_resampleScale.setValue(0.25f);
133
134
135
136
        // Set transfer function
        Geometry1DTransferFunction* tf = new Geometry1DTransferFunction(256);
        tf->addGeometry(TFGeometry1D::createQuad(cgt::vec2(0.0f, 0.5f), cgt::col4(0, 0, 0, 255), cgt::col4(0, 0, 0, 0)));
        _usFusion.p_confidenceTF.replaceTF(tf);
137
138

        // Create connectors
139
140
        _usIgtlReader.p_targetImagePrefix.setValue("us.igtl.");

141
142
        _usCropFilter.p_inputImage.setValue("us.igtl.CAMPUS");
        _usCropFilter.p_outputImage.setValue("us");
143
144

        _usBlurFilter.p_inputImage.setValue("us");
145
        _usBlurFilter.p_outputImage.setValue("us.blurred");
146
147
        _usBlurFilter.p_outputImage.addSharedProperty(&_usResampler.p_inputImage);
        _usBlurFilter.p_outputImage.addSharedProperty(&_usFusion.p_blurredImageId);
148

149
150
        _usResampler.p_outputImage.setValue("us.resampled");
        _usResampler.p_outputImage.addSharedProperty(&_usMapsSolver.p_inputImage);
151

152
        _usMapsSolver.p_outputConfidenceMap.setValue("us.confidence");
153
154
        _usMapsSolver.p_outputConfidenceMap.addSharedProperty(&_usFusion.p_confidenceImageID);

155
        _usFusion.p_usImageId.setValue("us");
156
157
        _usFusion.p_targetImageID.setValue("us.fusion");
        _usFusion.p_view.setValue(12);
158
159
        _usFusion.p_renderToTexture.setValue(true);
        _usFusion.p_targetImageID.addSharedProperty(&_usFanRenderer.p_inputImage);
160
161
        _usFusion.p_transferFunction.setAutoFitWindowToData(false);
        _usFusion.p_confidenceTF.setAutoFitWindowToData(false);
162

163
164
165
166
167
        _usFanRenderer.p_renderTargetID.setValue("us.fused_fan");
        _usFanRenderer.p_innerRadius.setValue(120.0f/540.0f);
        _usFanRenderer.p_halfAngle.setValue(37);

        _renderTargetID.setValue("us.fused_fan");
168

169
        // Bind buttons to event handlers
170
        p_connectDisconnectButton.s_clicked.connect(this, &CudaConfidenceMapsDemo::toggleIGTLConnection);
171
172
173
174
175
176
        p_copyStatisticsToClipboard.s_clicked.connect(this, &CudaConfidenceMapsDemo::copyStatisticsToClipboard);

        // Bind pipeline proeprties to processor properties
        p_useFixedIterationCount.addSharedProperty(&_usMapsSolver.p_useFixedIterationCount);
        p_millisecondBudget.addSharedProperty(&_usMapsSolver.p_millisecondBudget);
        p_iterationBudget.addSharedProperty(&_usMapsSolver.p_iterationBudget);
177

178
179
        p_gaussianFilterSize.addSharedProperty(&_usBlurFilter.p_sigma);
        p_resamplingScale.addSharedProperty(&_usResampler.p_resampleScale);
Declara Denis's avatar
Declara Denis committed
180
        p_gradientScaling.addSharedProperty(&_usMapsSolver.p_gradientScaling);
181
182
183
184
185
186
        p_alpha.addSharedProperty(&_usMapsSolver.p_paramAlpha);
        p_beta.addSharedProperty(&_usMapsSolver.p_paramBeta);
        p_gamma.addSharedProperty(&_usMapsSolver.p_paramGamma);
        p_useAlphaBetaFilter.addSharedProperty(&_usMapsSolver.p_useAlphaBetaFilter);
        p_fanHalfAngle.addSharedProperty(&_usFanRenderer.p_halfAngle);
        p_fanInnerRadius.addSharedProperty(&_usFanRenderer.p_innerRadius);
187
188
189
190
191
192
    }

    void CudaConfidenceMapsDemo::deinit() {
        AutoEvaluationPipeline::deinit();
    }

193
    void CudaConfidenceMapsDemo::executePipeline() {
194
195
196
197
198
199
        // Only launch the pipeline if the IgtlReader has recieved new data
        // FIXME: It would be better to check if a new image actually arrived, instead
        // of just checking the invaildation state of the IGTLReader
        if (!_usIgtlReader.isValid()) {
            auto startTime = tbb::tick_count::now();

200
            // Make sure that the whole pipeline gets invalidated and executed
201
202
203
204
205
206
207
208
209
210
211
            _usBlurFilter.invalidate(AbstractProcessor::INVALID_RESULT);
            _usCropFilter.invalidate(AbstractProcessor::INVALID_RESULT);
            _usResampler.invalidate(AbstractProcessor::INVALID_RESULT);
            _usMapsSolver.invalidate(AbstractProcessor::INVALID_RESULT);
            _usFusion.invalidate(AbstractProcessor::INVALID_RESULT);
            executeProcessorAndCheckOpenGLState(&_usIgtlReader);
            executeProcessorAndCheckOpenGLState(&_usCropFilter);
            executeProcessorAndCheckOpenGLState(&_usBlurFilter) ;
            executeProcessorAndCheckOpenGLState(&_usResampler);
            executeProcessorAndCheckOpenGLState(&_usMapsSolver);

212
213
214
215
216
217
218
219
220
221
222
            // Read fan geomtry from encoded image...
            if (p_useSpacingEncodedFanGeometry.getValue()) {
                ImageRepresentationGL::ScopedRepresentation img(*_data, _usCropFilter.p_inputImage.getValue());
                auto image = reinterpret_cast<const ImageData*>(_data->getData(_usCropFilter.p_inputImage.getValue()).getData());
                if (image != nullptr) {
                    cgt::vec3 encodedData = image->getMappingInformation().getVoxelSize();
                    p_fanHalfAngle.setValue(encodedData.x / 2.0f);
                    p_fanInnerRadius.setValue(encodedData.y);
                }
            }

223
224
            executeProcessorAndCheckOpenGLState(&_usFusion);
            executeProcessorAndCheckOpenGLState(&_usFanRenderer);
225
                
226
227
228
229
230
231
            auto endTime = tbb::tick_count::now();

            if ((startTime - _statisticsLastUpdateTime).seconds() > 0.5f) {
                _statisticsLastUpdateTime = startTime;

                auto ms = (endTime - startTime).seconds() * 1000.0f;
232
                auto solverMs = _usMapsSolver.getActualSolverExecutionTime();
233
                std::stringstream string;
234
                string << "Mode: " << _usFusion.p_view.getOptionValue() << std::endl;
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
                string << "Execution time: " << static_cast<int>(ms) << "ms" << std::endl;
                string << "Solver time: " << static_cast<int>(solverMs) << "ms" << std::endl;
                string << "CG Iterations: " << _usMapsSolver.getActualConjugentGradientIterations() << std::endl;
                string << "Error: " << _usMapsSolver.getResidualNorm() << std::endl;
                _usFanRenderer.p_text.setValue(string.str());
            }

            if (p_enableRecording.getValue()) {
                _recordedFrames++;
                QDir dir(QString::fromStdString(p_recordingDirectory.getValue()));
                std::stringstream filename;
                filename << _filePrefix << std::setw(4) << std::setfill('0') << _recordedFrames << ".png";
                auto url = dir.absoluteFilePath(QString::fromStdString(filename.str())).toStdString();

                ScopedTypedData<ImageData> rd(*_data, "us");
                const ImageRepresentationGL *rep = rd->getRepresentation<ImageRepresentationGL>();
                if (rep != 0) {
#ifdef CAMPVIS_HAS_MODULE_DEVIL
                    if (!cgt::FileSystem::dirExists(dir.absolutePath().toStdString()))
                        cgt::FileSystem::createDirectory(dir.absolutePath().toStdString());

                    // get color buffer content
                    GLubyte* colorBuffer = rep->getTexture()->downloadTextureToBuffer(GL_RED, GL_UNSIGNED_BYTE);
                    cgt::ivec2 size = rep->getSize().xy();

                    // create Devil image from image data and write it to file
                    ILuint img;
                    ilGenImages(1, &img);
                    ilBindImage(img);

                    // put pixels into IL-Image
                    ilTexImage(size.x, size.y, 1, 1, GL_LUMINANCE, IL_UNSIGNED_BYTE, colorBuffer);
                    ilEnable(IL_FILE_OVERWRITE);
                    ilResetWrite();
                    ILboolean success = ilSaveImage(url.c_str());
                    ilDeleteImages(1, &img);

                    delete[] colorBuffer;

                    if (!success) {
                        LERROR("Could not save image to file: " << ilGetError());
                    }
#endif
278
279
                }
            }
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312

            // Collect statistics
            if (p_collectStatistics.getValue()) {
                StatisticsEntry entry;
                entry.time = (startTime - _objectCreationTime).seconds() * 1000.0f;

                ImageRepresentationGL::ScopedRepresentation originalImage(*_data, _usCropFilter.p_outputImage.getValue());
                ImageRepresentationGL::ScopedRepresentation downsampledImage(*_data, _usResampler.p_outputImage.getValue());

                if (originalImage && downsampledImage) {
                    entry.originalWidth = originalImage->getSize().x;
                    entry.originalHeight = originalImage->getSize().y;
                    entry.downsampledWidth = downsampledImage->getSize().x;
                    entry.downsampledHeight = downsampledImage->getSize().y;
                } else {
                    LWARNING("Could not read image size");
                    entry.originalWidth = entry.originalHeight = -1;
                    entry.downsampledWidth = entry.downsampledHeight = -1;
                }

                entry.gaussianKernelSize = _usBlurFilter.p_sigma.getValue();
                entry.scalingFactor = _usResampler.p_resampleScale.getValue();
                entry.alpha = _usMapsSolver.p_paramAlpha.getValue();
                entry.beta = _usMapsSolver.p_paramBeta.getValue();
                entry.gamma = _usMapsSolver.p_paramGamma.getValue();
                entry.gradientScaling = _usMapsSolver.p_gradientScaling.getValue();
                entry.iterations = _usMapsSolver.getActualConjugentGradientIterations();
                entry.solverExecutionTime = _usMapsSolver.getActualSolverExecutionTime();
                entry.totalExecutionTime = (endTime - startTime).seconds() * 1000.0f;
                entry.solverError = _usMapsSolver.getResidualNorm();

                _statistics.push_back(entry);
            }
313
314
315
        }
    }

316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
    void CudaConfidenceMapsDemo::onEvent(cgt::Event* e) {
        // Allow for rapid switching between different visualizations
        //  F1: Ultrasound only
        //  F2: Sharpness
        //  F3: LAB
        //  F4: Color overlay
        //  F5: CM only
        if (typeid(*e) == typeid(cgt::KeyEvent)) {
            cgt::KeyEvent *keyEvent = reinterpret_cast<cgt::KeyEvent*>(e);
            if (keyEvent == nullptr) return;

            if (keyEvent->pressed()) {
                bool eventHandled = true;
                switch (keyEvent->keyCode()) {
                case cgt::KeyEvent::K_F1:
                    _usFusion.p_view.setValue(0); // US only
                    break;
                case cgt::KeyEvent::K_F2:
                    _usFusion.p_view.setValue(10); // Sharpness
                    {
                        Geometry1DTransferFunction* tf = new Geometry1DTransferFunction(256);
                        tf->addGeometry(TFGeometry1D::createQuad(cgt::vec2(0.0f, 1.0f), cgt::col4(0, 0, 0, 255), cgt::col4(0, 0, 0, 0)));
                        _usFusion.p_confidenceTF.replaceTF(tf);
                    }
                    break;
                case cgt::KeyEvent::K_F3:
                    _usFusion.p_view.setValue(8); // LAB
                    {
                        Geometry1DTransferFunction* tf = new Geometry1DTransferFunction(256);
                        tf->addGeometry(TFGeometry1D::createQuad(cgt::vec2(0.0f, 0.5f), cgt::col4(0, 0, 0, 255), cgt::col4(0, 0, 0, 0)));
                        _usFusion.p_confidenceTF.replaceTF(tf);
                        _usFusion.p_hue.setValue(0.23f);
                    }
                    break;
                case cgt::KeyEvent::K_F4:
                    _usFusion.p_view.setValue(12); // Color overlay
                    {
                        Geometry1DTransferFunction* tf = new Geometry1DTransferFunction(256);
                        tf->addGeometry(TFGeometry1D::createQuad(cgt::vec2(0.0f, 0.5f), cgt::col4(0, 0, 0, 255), cgt::col4(0, 0, 0, 0)));
                        _usFusion.p_confidenceTF.replaceTF(tf);
                        _usFusion.p_hue.setValue(0.15f);
                    }
                    break;
                case cgt::KeyEvent::K_F5:
                    _usFusion.p_view.setValue(2); // CM Only
                    break;
                default:
                    eventHandled = false;
                };
                if (eventHandled) {
                    e->accept();
                    // Force HUD statistics to be updated
                    _statisticsLastUpdateTime = tbb::tick_count();
                }
            }
        }
    }

374
375
376
377
378
379
380
381
382
    void CudaConfidenceMapsDemo::onPropertyChanged(const AbstractProperty* prop) {
        if (prop == &p_showAdvancedOptions) {
            setAdvancedPropertiesVisibility(p_showAdvancedOptions.getValue());
            return;
        }
        else if (prop == &_renderTargetID && _renderTargetID.getValue() != "us.fused_fan") {
            // Prevent the program from changing the render target
            _renderTargetID.setValue("us.fused_fan");
        }
383
384
385
386
387

        bool useFixedIterationCount = p_useFixedIterationCount.getValue();
        p_millisecondBudget.setVisible(!useFixedIterationCount);
        p_iterationBudget.setVisible(useFixedIterationCount);

388
389
390
391
392
393
394
395
396
397
        AutoEvaluationPipeline::onPropertyChanged(prop);
    }

    void CudaConfidenceMapsDemo::toggleIGTLConnection() {
        // Simulate a click on the currently visible button of the IGTL reader,
        // causing it to either connect or disconnect from a server.
        if (_usIgtlReader.p_connect.isVisible())
            _usIgtlReader.p_connect.click();
        else
            _usIgtlReader.p_disconnect.click();
398
    }
399

400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
    void CudaConfidenceMapsDemo::copyStatisticsToClipboard() {
        // Copy statistics to the clipboard in CSV format
        std::cout << "Have " << _statistics.size() << " stats items" << std::endl;

        std::stringstream stream;
        stream << "time, originalWidth, originalHeight, downsampledWidth, downsampledHeight, gaussianKernelSize, scalingFactor, alpha, beta, gamma, gradientScaling, iterations, solverExecutionTime, totalExecutionTime, solverError" << std::endl;

        for (auto& item : _statistics) {
            stream << item.time << ", ";
            stream << item.originalWidth << ", ";
            stream << item.originalHeight << ", ";
            stream << item.downsampledWidth << ", ";
            stream << item.downsampledHeight << ", ";
            stream << item.gaussianKernelSize << ", ";
            stream << item.scalingFactor << ", ";
            stream << item.alpha << ", ";
            stream << item.beta << ", ";
            stream << item.gamma << ", ";
            stream << item.gradientScaling << ", ";
            stream << item.iterations << ", ";
            stream << item.solverExecutionTime << ", ";
            stream << item.totalExecutionTime << ", ";
            stream << item.solverError << std::endl;
        }

        // Copy CSV stream to clipboard
        QClipboard *clipboard = QApplication::clipboard();
        QString text = QString::fromStdString(stream.str());
        clipboard->setText(text);

        _statistics.clear();
        _statistics.reserve(1000);
    }

434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
    void CudaConfidenceMapsDemo::setAdvancedPropertiesVisibility(bool visible) {
        if (p_showAdvancedOptions.getValue() != visible) {
            p_showAdvancedOptions.setValue(visible);
        }

        p_useAlphaBetaFilter.setVisible(visible);
        p_gaussianFilterSize.setVisible(visible);
        p_gradientScaling.setVisible(visible);
        p_alpha.setVisible(visible);
        p_gamma.setVisible(visible);
        p_fanHalfAngle.setVisible(visible);
        p_fanInnerRadius.setVisible(visible);
        p_recordingDirectory.setVisible(visible);
        p_enableRecording.setVisible(visible);
    }

450
}