<commit_before>/*
 *  qimask.cpp
 *
 *  Created by <NAME> on 2015/09/22.
 *  Copyright (c) 2015 <NAME>.
 *
 *  This Source Code Form is subject to the terms of the Mozilla Public
 *  License, v. 2.0. If a copy of the MPL was not distributed with this
 *  file, You can obtain one at http://mozilla.org/MPL/2.0/.
 *
 */

#include <iostream>
#include <string>
#include <complex>

#include "itkComplexToModulusImageFilter.h"
#include "itkOtsuThresholdImageFilter.h"
#include "itkVotingBinaryIterativeHoleFillingImageFilter.h"
#include "itkBinaryCrossStructuringElement.h"
#include "itkBinaryBallStructuringElement.h"
#include "itkBinaryErodeImageFilter.h"
#include "itkBinaryDilateImageFilter.h"

#include "Args.h"
#include "Types.h"
#include "Util.h"
#include "IO.h"
#include "Masking.h"

int main(int argc, char **argv) {
    args::ArgumentParser parser(
    "Generates masks in stages.\n"
    "Stage 1 - Otsu thresholding to generate binary mask\n"
    "Stage 2 - RATs (optional)\n"
    "Stage 3 - Hole filling (optional)\n"
    "http://github.com/spinicist/QUIT"
    );

    args::Positional<std::string> input_path(parser, "INPUT_FILE", "Input file");

    args::HelpFlag help(parser, "HELP", "Show this help menu", {'h', "help"});
    args::Flag     verbose(parser, "VERBOSE", "Print more information", {'v', "verbose"});
    args::ValueFlag<std::string> outarg(parser, "OUTPUT FILE", "Set output filename, default is input + _mask", {'o', "out"});
    args::ValueFlag<int> volume(parser, "VOLUME", "Choose volume to mask in multi-volume file. Default 1, -1 selects last volume", {'v', "volume"}, 0);
    args::Flag     is_complex(parser, "COMPLEX", "Input data is complex, take magnitude first", {'x', "complex"});
    args::ValueFlag<float> lower_threshold(parser, "LOWER THRESHOLD", "Specify lower intensity threshold for 1st stage, otherwise Otsu's method is used", {'l', "lower"}, 0.);
    args::ValueFlag<float> upper_threshold(parser, "UPPER THRESHOLD", "Specify upper intensity threshold for 1st stage, otherwise Otsu's method is used", {'u', "upper"}, std::numeric_limits<float>::infinity());
    args::ValueFlag<float> rats(parser, "RATS", "Perform the RATS step, argument is size threshold for connected component", {'r', "rats"}, 0.);
    args::ValueFlag<int> fillh_radius(parser, "FILL HOLES", "Fill holes in thresholded mask with radius N", {'F', "fillh"}, 0);
    args::Flag     run_bet(parser, "RUN BET", "Run the Brain Extraction Tool stage", {'B', "bet"});

    QI::ParseArgs(parser, argc, argv);

    if (verbose) std::cout << "Reading input image: " << QI::CheckPos(input_path) << std::endl;
    QI::SeriesF::Pointer vols = ITK_NULLPTR;
    if (is_complex) {
        vols = QI::ReadMagnitudeImage<QI::SeriesF>(QI::CheckPos(input_path));
    } else {
        vols = QI::ReadImage<QI::SeriesF>(QI::CheckPos(input_path));
    }

    std::string out_path;
    if (outarg.Get() == "") {
        out_path = QI::StripExt(input_path.Get()) + "_mask" + QI::OutExt();
    } else {
        out_path = outarg.Get();
    }

    // Extract one volume to process
    auto region = vols->GetLargestPossibleRegion();
    if (std::abs(volume.Get()) < region.GetSize()[3]) {
        int volume_to_get = (region.GetSize()[3] + volume.Get()) % region.GetSize()[3];
        if (verbose) std::cout << "Using volume " << volume_to_get << std::endl;
        region.GetModifiableIndex()[3] = volume_to_get;
    } else {
        std::cerr << "Specified volume was invalid: " << volume.Get() << std::endl;
        return EXIT_FAILURE;
    }
    region.GetModifiableSize()[3] = 0;
    auto vol = itk::ExtractImageFilter<QI::SeriesF, QI::VolumeF>::New();
    vol->SetExtractionRegion(region);
    vol->SetInput(vols);
    vol->SetDirectionCollapseToSubmatrix();
    vol->Update();
    QI::VolumeF::Pointer intensity_image = vol->GetOutput();
    intensity_image->DisconnectPipeline();

    /*
     *  Stage 1 - Otsu or Threshold
     */
    QI::VolumeI::Pointer mask_image = ITK_NULLPTR;
    if (lower_threshold || upper_threshold) {
        if (verbose) std::cout << "Thresholding range: " << lower_threshold.Get() << "-" << upper_threshold.Get() << std::endl;
        mask_image = QI::ThresholdMask(intensity_image, lower_threshold.Get(), upper_threshold.Get());
    } else {
        if (verbose) std::cout << "Generating Otsu mask" << std::endl;
        mask_image = QI::OtsuMask(intensity_image);
    }

    /*
     *  Stage 2 - RATS
     */
     if (rats) {
        typedef itk::BinaryBallStructuringElement<int, 3> TBall;
        typedef itk::BinaryErodeImageFilter<QI::VolumeI, QI::VolumeI, TBall> TErode;
        typedef itk::BinaryDilateImageFilter<QI::VolumeI, QI::VolumeI, TBall> TDilate;
        float mask_volume = std::numeric_limits<float>::infinity();
        float voxel_volume = QI::VoxelVolume(mask_image);
        std::cout << "Voxel volume: " << voxel_volume << std::endl;
        int radius = 0;
        QI::VolumeI::Pointer mask_rats;
        while (mask_volume > rats.Get()) {
            radius++;
            TBall ball;
            TBall::SizeType radii; radii.Fill(radius);
            ball.SetRadius(radii);<commit_msg>BUG: Missing a verbose check.<commit_after>/*
 *  qimask.cpp
 *
 *  Created by <NAME> on 2015/09/22.
 *  Copyright (c) 2015 <NAME>.
 *
 *  This Source Code Form is subject to the terms of the Mozilla Public
 *  License, v. 2.0. If a copy of the MPL was not distributed with this
 *  file, You can obtain one at http://mozilla.org/MPL/2.0/.
 *
 */

#include <iostream>
#include <string>
#include <complex>

#include "itkComplexToModulusImageFilter.h"
#include "itkOtsuThresholdImageFilter.h"
#include "itkVotingBinaryIterativeHoleFillingImageFilter.h"
#include "itkBinaryCrossStructuringElement.h"
#include "itkBinaryBallStructuringElement.h"
#include "itkBinaryErodeImageFilter.h"
#include "itkBinaryDilateImageFilter.h"

#include "Args.h"
#include "Types.h"
#include "Util.h"
#include "IO.h"
#include "Masking.h"

int main(int argc, char **argv) {
    args::ArgumentParser parser(
    "Generates masks in stages.\n"
    "Stage 1 - Otsu thresholding to generate binary mask\n"
    "Stage 2 - RATs (optional)\n"
    "Stage 3 - Hole filling (optional)\n"
    "http://github.com/spinicist/QUIT"
    );

    args::Positional<std::string> input_path(parser, "INPUT_FILE", "Input file");

    args::HelpFlag help(parser, "HELP", "Show this help menu", {'h', "help"});
    args::Flag     verbose(parser, "VERBOSE", "Print more information", {'v', "verbose"});
    args::ValueFlag<std::string> outarg(parser, "OUTPUT FILE", "Set output filename, default is input + _mask", {'o', "out"});
    args::ValueFlag<int> volume(parser, "VOLUME", "Choose volume to mask in multi-volume file. Default 1, -1 selects last volume", {'v', "volume"}, 0);
    args::Flag     is_complex(parser, "COMPLEX", "Input data is complex, take magnitude first", {'x', "complex"});
    args::ValueFlag<float> lower_threshold(parser, "LOWER THRESHOLD", "Specify lower intensity threshold for 1st stage, otherwise Otsu's method is used", {'l', "lower"}, 0.);
    args::ValueFlag<float> upper_threshold(parser, "UPPER THRESHOLD", "Specify upper intensity threshold for 1st stage, otherwise Otsu's method is used", {'u', "upper"}, std::numeric_limits<float>::infinity());
    args::ValueFlag<float> rats(parser, "RATS", "Perform the RATS step, argument is size threshold for connected component", {'r', "rats"}, 0.);
    args::ValueFlag<int> fillh_radius(parser, "FILL HOLES", "Fill holes in thresholded mask with radius N", {'F', "fillh"}, 0);
    args::Flag     run_bet(parser, "RUN BET", "Run the Brain Extraction Tool stage", {'B', "bet"});

    QI::ParseArgs(parser, argc, argv);

    if (verbose) std::cout << "Reading input image: " << QI::CheckPos(input_path) << std::endl;
    QI::SeriesF::Pointer vols = ITK_NULLPTR;
    if (is_complex) {
        vols = QI::ReadMagnitudeImage<QI::SeriesF>(QI::CheckPos(input_path));
    } else {
        vols = QI::ReadImage<QI::SeriesF>(QI::CheckPos(input_path));
    }

    std::string out_path;
    if (outarg.Get() == "") {
        out_path = QI::StripExt(input_path.Get()) + "_mask" + QI::OutExt();
    } else {
        out_path = outarg.Get();
    }

    // Extract one volume to process
    auto region = vols->GetLargestPossibleRegion();
    if (std::abs(volume.Get()) < region.GetSize()[3]) {
        int volume_to_get = (region.GetSize()[3] + volume.Get()) % region.GetSize()[3];
        if (verbose) std::cout << "Using volume " << volume_to_get << std::endl;
        region.GetModifiableIndex()[3] = volume_to_get;
    } else {
        std::cerr << "Specified volume was invalid: " << volume.Get() << std::endl;
        return EXIT_FAILURE;
    }
    region.GetModifiableSize()[3] = 0;
    auto vol = itk::ExtractImageFilter<QI::SeriesF, QI::VolumeF>::New();
    vol->SetExtractionRegion(region);
    vol->SetInput(vols);
    vol->SetDirectionCollapseToSubmatrix();
    vol->Update();
    QI::VolumeF::Pointer intensity_image = vol->GetOutput();
    intensity_image->DisconnectPipeline();

    /*
     *  Stage 1 - Otsu or Threshold
     */
    QI::VolumeI::Pointer mask_image = ITK_NULLPTR;
    if (lower_threshold || upper_threshold) {
        if (verbose) std::cout << "Thresholding range: " << lower_threshold.Get() << "-" << upper_threshold.Get() << std::endl;
        mask_image = QI::ThresholdMask(intensity_image, lower_threshold.Get(), upper_threshold.Get());
    } else {
        if (verbose) std::cout << "Generating Otsu mask" << std::endl;
        mask_image = QI::OtsuMask(intensity_image);
    }

    /*
     *  Stage 2 - RATS
     */
     if (rats) {
        typedef itk::BinaryBallStructuringElement<int, 3> TBall;
        typedef itk::BinaryErodeImageFilter<QI::VolumeI, QI::VolumeI, TBall> TErode;
        typedef itk::BinaryDilateImageFilter<QI::VolumeI, QI::VolumeI, TBall> TDilate;
        float mask_volume = std::numeric_limits<float>::infinity();
        float voxel_volume = QI::VoxelVolume(mask_image);
        if (verbose) std::cout << "Voxel volume: " << voxel_volume << std::endl;
        int radius = 0;
        QI::VolumeI::Pointer mask_rats;
        while (mask_volume > rats.Get()) {
            radius++;
            TBall ball;
            TBall::SizeType radii; radii.Fill(radius);
            ball.SetRadius(radii);
{'fraction_non_alphanumeric': 0.15730772960915526, 'fraction_numerical': 0.007176801474153816, 'mean_word_length': 3.823199251637044, 'pattern_counts': {'":': 0, '<': 121, '<?xml version=': 0, '>': 59, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '6576300', 'n_tokens_mistral': 3512, 'n_tokens_neox': 3325, 'n_words': 827}