C++ pointer::getImageAccess方法代码示例

void executeAlgorithmOnHost(Image::pointer input, Image::pointer output) {
    ImageAccess inputAccess = input->getImageAccess(ACCESS_READ);
    ImageAccess outputAccess = output->getImageAccess(ACCESS_READ_WRITE);

    T * inputData = (T*)inputAccess.get();
    T * outputData = (T*)outputAccess.get();

    unsigned int nrOfElements = input->getWidth()*input->getHeight()*input->getDepth()*input->getNrOfComponents();
    for(unsigned int i = 0; i < nrOfElements; i++) {
        outputData[i] = 2.0*inputData[i];
    }
}

void SeededRegionGrowing::executeOnHost(T* input, Image::pointer output) {
    ImageAccess::pointer outputAccess = output->getImageAccess(ACCESS_READ_WRITE);
    uchar* outputData = (uchar*)outputAccess->get();
    // initialize output to all zero
    memset(outputData, 0, output->getWidth()*output->getHeight()*output->getDepth());
    std::stack<Vector3ui> queue;

    // Add seeds to queue
    for(int i = 0; i < mSeedPoints.size(); i++) {
        Vector3ui pos = mSeedPoints[i];

        // Check if seed point is in bounds
        if(pos.x() < 0 || pos.y() < 0 || pos.z() < 0 ||
            pos.x() >= output->getWidth() || pos.y() >= output->getHeight() || pos.z() >= output->getDepth())
            throw Exception("One of the seed points given to SeededRegionGrowing was out of bounds.");

        queue.push(pos);
    }

    // Process queue
    while(!queue.empty()) {
        Vector3ui pos = queue.top();
        queue.pop();

        // Add neighbors to queue
        for(int a = -1; a < 2; a++) {
        for(int b = -1; b < 2; b++) {
        for(int c = -1; c < 2; c++) {
            if(abs(a)+abs(b)+abs(c) != 1) // connectivity
                continue;
            Vector3ui neighbor(pos.x()+a,pos.y()+b,pos.z()+c);
            // Check for out of bounds
            if(neighbor.x() < 0 || neighbor.y() < 0 || neighbor.z() < 0 ||
                neighbor.x() >= output->getWidth() || neighbor.y() >= output->getHeight() || neighbor.z() >= output->getDepth())
                continue;

            // Check that voxel is not already segmented
            if(outputData[neighbor.x()+neighbor.y()*output->getWidth()+neighbor.z()*output->getWidth()*output->getHeight()] == 1)
                continue;

            // Check condition
            T value = input[neighbor.x()+neighbor.y()*output->getWidth()+neighbor.z()*output->getWidth()*output->getHeight()];
            if(value >= mMinimumIntensity && value <= mMaximumIntensity) {
                // add it to segmentation
                outputData[neighbor.x()+neighbor.y()*output->getWidth()+neighbor.z()*output->getWidth()*output->getHeight()] = 1;

                // Add to queue
                queue.push(neighbor);
            }
        }}}
    }
}

static igtl::ImageMessage::Pointer createIGTLImageMessage(Image::pointer image) {
    // size parameters
    int   size[3]     = {(int)image->getWidth(), (int)image->getHeight(), (int)image->getDepth()};       // image dimension
    float spacing[3]  = {image->getSpacing().x(), image->getSpacing().y(), image->getSpacing().z()};     // spacing (mm/pixel)
    int   svoffset[3] = {0, 0, 0};           // sub-volume offset
    int   scalarType;
    size_t totalSize = image->getWidth()*image->getHeight()*image->getDepth()*image->getNrOfChannels();
    switch(image->getDataType()) {
        case TYPE_UINT8:
            scalarType = igtl::ImageMessage::TYPE_UINT8;
            totalSize *= sizeof(unsigned char);
            break;
        case TYPE_INT8:
            scalarType = igtl::ImageMessage::TYPE_INT8;
            totalSize *= sizeof(char);
            break;
        case TYPE_UINT16:
            scalarType = igtl::ImageMessage::TYPE_UINT16;
            totalSize *= sizeof(unsigned short);
            break;
        case TYPE_INT16:
            scalarType = igtl::ImageMessage::TYPE_INT16;
            totalSize *= sizeof(short);
            break;
        case TYPE_FLOAT:
            scalarType = igtl::ImageMessage::TYPE_FLOAT32;
            totalSize *= sizeof(float);
            break;
    }

    //------------------------------------------------------------
    // Create a new IMAGE type message
    igtl::ImageMessage::Pointer imgMsg = igtl::ImageMessage::New();
    imgMsg->SetDimensions(size);
    imgMsg->SetSpacing(spacing);
    imgMsg->SetNumComponents(image->getNrOfChannels());
    imgMsg->SetScalarType(scalarType);
    imgMsg->SetDeviceName("DummyImage");
    imgMsg->SetSubVolume(size, svoffset);
    imgMsg->AllocateScalars();

    ImageAccess::pointer access = image->getImageAccess(ACCESS_READ);
    memcpy(imgMsg->GetScalarPointer(), access->get(), totalSize);

    return imgMsg;
}

void executeAlgorithmOnHost(Image::pointer input, Image::pointer output, unsigned char group, unsigned char window, float strength, unsigned char sigma) {
	throw Exception("This is on host, does not work atm");
    
    ImageAccess::pointer inputAccess = input->getImageAccess(ACCESS_READ);
    ImageAccess::pointer outputAccess = output->getImageAccess(ACCESS_READ_WRITE);

    T * inputData = (T*)inputAccess->get();
    T * outputData = (T*)outputAccess->get();

    unsigned int width = input->getWidth();
    unsigned int height = input->getHeight();
    //Window is window-1/2
    //group is group-1/2
    //strength is strength*strength
    //sigma is sigma*sigma
    //Not working atm with the T
    //Does not work with outofbounds atm
    //So this code is for all pixels inbound, meaning x + group + window < width / x - group - window > 0 //same for y
    for (int x = 0; x < width; x++){
        for (int y = 0; y < height; y++){
            
            double normSum = 0.0;
            double totSum = 0.0;
            double indi = 0.0;
            double groupTot = 0.0;
            double value = 0.0;

            for (int i = x - window; i <= x + window; i++){
                for (int j = y - window; j <= y + window; j++){
                    if (i != x && j != y){

                        int mX = x - group;
                        int mY = y - group;
                        for (int k = i - group; k <= i + group; k++, mX++){
                            for (int l = j - group; l <= j + group; l++, mY++){
                                //This is wrong, need to fix T
                                //indi = inputData[mX][mY] - inputData[k][l];
                                indi = abs(indi*indi);
                                indi = exp( - (indi/strength));
                                groupTot += indi;
                            }
                        }
                        //This is wrong, need to fix T
                        //value = inputData[i][j];
                        double pA[] = {i,j};
                        double pB[] = {x,y};
                        //double dist = i, j - x, y;
                        double dist = std::inner_product(std::begin(pA), std::end(pA), std::begin(pB), 0.0);
                        
                        double gaussWeight = exp(-(dist / (2.0 * sigma)));
                        gaussWeight = gaussWeight / (2.0 * sigma);
                        groupTot *= gaussWeight;

                        normSum += groupTot;
                        totSum += groupTot * value;
                        groupTot = 0.0;
                    }
                }
            }
            value = totSum / normSum;
            /*
            Not sure it needed
            if (value < 0){
                value = 0;
            }
            if (value > 1.0){
                value = 1.0f;
            }
            */
            //This is wrong, need to fix T
            //outputData[x][y] = (T)value;
            
        }
    }
}

void SeededRegionGrowing::execute() {
    if(mSeedPoints.size() == 0)
        throw Exception("No seed points supplied to SeededRegionGrowing");

    Image::pointer input = getStaticInputData<Image>();
    if(input->getNrOfComponents() != 1)
        throw Exception("Seeded region growing currently doesn't support images with several components.");

    Segmentation::pointer output = getStaticOutputData<Segmentation>();

    // Initialize output image
    output->createFromImage(input, getMainDevice());

    if(getMainDevice()->isHost()) {
        ImageAccess::pointer inputAccess = input->getImageAccess(ACCESS_READ);
        void* inputData = inputAccess->get();
        switch(input->getDataType()) {
            fastSwitchTypeMacro(executeOnHost<FAST_TYPE>((FAST_TYPE*)inputData, output));
        }
    } else {
        OpenCLDevice::pointer device = getMainDevice();

        recompileOpenCLCode(input);

        ImageAccess::pointer access = output->getImageAccess(ACCESS_READ_WRITE);
        uchar* outputData = (uchar*)access->get();
        // Initialize to all 0s
        memset(outputData,0,sizeof(uchar)*output->getWidth()*output->getHeight()*output->getDepth());

        // Add sedd points
        for(int i = 0; i < mSeedPoints.size(); i++) {
            Vector3ui pos = mSeedPoints[i];

            // Check if seed point is in bounds
            if(pos.x() < 0 || pos.y() < 0 || pos.z() < 0 ||
                pos.x() >= output->getWidth() || pos.y() >= output->getHeight() || pos.z() >= output->getDepth())
                throw Exception("One of the seed points given to SeededRegionGrowing was out of bounds.");

            outputData[pos.x() + pos.y()*output->getWidth() + pos.z()*output->getWidth()*output->getHeight()] = 2;
        }
        access->release();

        cl::NDRange globalSize;
        if(output->getDimensions() == 2) {
            globalSize = cl::NDRange(input->getWidth(),input->getHeight());
            OpenCLImageAccess2D::pointer inputAccess = input->getOpenCLImageAccess2D(ACCESS_READ, device);
            mKernel.setArg(0, *inputAccess->get());
        } else {
            globalSize = cl::NDRange(input->getWidth(),input->getHeight(), input->getDepth());
            OpenCLImageAccess3D::pointer inputAccess = input->getOpenCLImageAccess3D(ACCESS_READ, device);
            mKernel.setArg(0, *inputAccess->get());
        }

        OpenCLBufferAccess::pointer outputAccess = output->getOpenCLBufferAccess(ACCESS_READ_WRITE, device);
        cl::Buffer stopGrowingBuffer = cl::Buffer(
                device->getContext(),
                CL_MEM_READ_WRITE,
                sizeof(char));
        cl::CommandQueue queue = device->getCommandQueue();
        mKernel.setArg(1, *outputAccess->get());
        mKernel.setArg(2, stopGrowingBuffer);
        mKernel.setArg(3, mMinimumIntensity);
        mKernel.setArg(4, mMaximumIntensity);

        bool stopGrowing = false;
        char stopGrowingInit = 1;
        char * stopGrowingResult = new char;
        int iterations = 0;
        do {
            iterations++;
            queue.enqueueWriteBuffer(stopGrowingBuffer, CL_TRUE, 0, sizeof(char), &stopGrowingInit);

            queue.enqueueNDRangeKernel(
                    mKernel,
                    cl::NullRange,
                    globalSize,
                    cl::NullRange
            );

            queue.enqueueReadBuffer(stopGrowingBuffer, CL_TRUE, 0, sizeof(char), stopGrowingResult);
            if(*stopGrowingResult == 1)
                stopGrowing = true;
        } while(!stopGrowing);
    }

}