Page MenuHomec4science

CorrelationMatrixCreatorApplicationParallel.cpp
No OneTemporary

File Metadata

Created
Fri, May 3, 20:00

CorrelationMatrixCreatorApplicationParallel.cpp

#include <CorrelationMatrixCreatorApplication.hpp>
#include <MatrixCreator.hpp> // MatrixCreator::methods
#include <CorrelationMatrixCreator.hpp>
#include <boost/program_options.hpp>
#include <iostream>
#include <string>
#include <stdexcept> // std::invalid_argument
namespace po = boost::program_options ;
// the valid values for --method option
std::string method_read = "read" ;
std::string method_read_atac = "read_atac" ;
std::string method_fragment = "fragment" ;
std::string method_fragment_center = "fragment_center" ;
CorrelationMatrixCreatorApplication::CorrelationMatrixCreatorApplication(int argn, char** argv)
: file_bed(""), file_bam(""), from(0), to(0), bin_size(0),
method(MatrixCreator::FRAGMENT), runnable(true), n_threads(1)
{
// parse command line options and set the fields
this->parseOptions(argn, argv) ;
}
int CorrelationMatrixCreatorApplication::run()
{ if(this->runnable)
{ CorrelationMatrixCreator mc(this->file_bed,
this->file_bam,
this->file_bai,
this->from,
this->to,
this->bin_size,
this->method,
this->n_threads) ;
std::cout << mc.create_matrix() << std::endl ;
return EXIT_SUCCESS ;
}
else
{ return EXIT_FAILURE ; }
}
void CorrelationMatrixCreatorApplication::parseOptions(int argn, char** argv)
{
// no option to parse
if(argv == nullptr)
{ std::string message = "no options to parse!" ;
throw std::invalid_argument(message) ;
}
// help messages
std::string desc_msg = "\n"
"CorrelationMatrixCreator is an application that creates a "
"count matrix from a BED file and a BAM file and returnes it "
"through stdout.\n"
"The matrix contains one row per reference region present in the "
"BED file. The region center is computed and then a region covering the "
"interval [from,to] is build around the middle and divided into "
"equally sized bins. Finally, each bin is assigned the number of "
"target present in the BAM file that are mapped at that position.\n\n" ;
std::string opt_help_msg = "Produces this help message." ;
std::string opt_bed_msg = "The path to the BED file containing the references";
std::string opt_bam_msg = "The path to the BAM file containing the targets";
std::string opt_bai_msg = "The path to the BAM index file of the BAM file containing the targets";
std::string opt_from_msg = "The upstream limit - in relative coordinate - of the region to build "
"around each reference center." ;
std::string opt_to_msg = "The downstream limit - in relative coordinate - of the region to build "
"around each reference center." ;
std::string opt_thread_msg = "The number of threads to use." ;
std::string opt_binsize_msg = "The size of the bins." ;
char tmp[4096] ;
sprintf(tmp,
"How the data in the BAM file should be handled when computing "
"the number of counts in each bin.\n"
"\t\"%s\" uses each position within the reads (by default)\n"
"\t\"%s\" uses only the insertion site for ATAC-seq data\n"
"\t\"%s\" uses each position within the fragments\n"
"\t\"%s\" uses only the fragment central positions\n",
method_read.c_str(),
method_read_atac.c_str(),
method_fragment.c_str(),
method_fragment_center.c_str()) ;
std::string opt_method_msg = tmp ;
// option parser
boost::program_options::variables_map vm ;
boost::program_options::options_description desc(desc_msg) ;
std::string method(method_read) ;
desc.add_options()
("help,h", opt_help_msg.c_str())
("bed", po::value<std::string>(&(this->file_bed)), opt_bed_msg.c_str())
("bam", po::value<std::string>(&(this->file_bam)), opt_bam_msg.c_str())
("bai", po::value<std::string>(&(this->file_bai)), opt_bai_msg.c_str())
("from,f", po::value<int>(&(this->from)), opt_from_msg.c_str())
("to,t", po::value<int>(&(this->to)), opt_to_msg.c_str())
("binSize,b", po::value<int>(&(this->bin_size)), opt_binsize_msg.c_str())
("method,m", po::value<std::string>(&(method)), opt_method_msg.c_str())
("parallel,p", po::value<size_t>(&(this->n_threads)), opt_thread_msg.c_str()) ;
// parse
try
{ po::store(po::parse_command_line(argn, argv, desc), vm) ;
po::notify(vm) ;
}
catch(std::invalid_argument& e)
{ std::string msg = std::string("Error! Invalid option given!\n") + std::string(e.what()) ;
throw std::invalid_argument(msg) ;
}
catch(...)
{ throw std::invalid_argument("An unknown error occured while parsing the options") ; }
bool help = vm.count("help") ;
// checks unproper option settings
if(this->file_bed == "" and (not help))
{ std::string msg("Error! No BED file was given (--bed)!") ;
throw std::invalid_argument(msg) ;
}
else if(this->file_bam == "" and (not help))
{ std::string msg("Error! No BAM file was given (--bam)!") ;
throw std::invalid_argument(msg) ;
}
else if(this->file_bam == "" and (not help))
{ std::string msg("Error! No BAM index file was given (--bai)!") ;
throw std::invalid_argument(msg) ;
}
else if(this->from == 0 and this->to == 0 and (not help))
{ std::string msg("Error! No range given (--from and --to)!") ;
throw std::invalid_argument(msg) ;
}
else if(this->from >= this->to and (not help))
{ std::string msg("Error! from shoud be smaller than to (--from and --to)!") ;
throw std::invalid_argument(msg) ;
}
else if(this->bin_size <= 0 and (not help))
{ std::string msg("Error! bin size should be bigger than 0 (--binSize)!") ;
throw std::invalid_argument(msg) ;
}
else if(method != method_read and
method != method_read_atac and
method != method_fragment and
method != method_fragment_center)
{ char msg[4096] ;
sprintf(msg, "Error! method should be %s, %s, %s or %s (--method)",
method_read.c_str(),
method_read_atac.c_str(),
method_fragment.c_str(),
method_fragment_center.c_str()) ;
throw std::invalid_argument(msg) ;
}
else if(this->n_threads == 0)
{ std::string msg("Error! at least one thread should be used (--parallel)!") ;
throw std::invalid_argument(msg) ;
}
// set method
if(method == method_read)
{ this->method = MatrixCreator::READ ; }
else if(method == method_read_atac)
{ this->method = MatrixCreator::READ_ATAC ; }
else if(method == method_fragment)
{ this->method = MatrixCreator::FRAGMENT ; }
else if(method == method_fragment_center)
{ this->method = MatrixCreator::FRAGMENT_CENTER ; }
// help invoked, run() cannot be invoked
if(help)
{ std::cout << desc << std::endl ;
this->runnable = false ;
return ;
}
// everything fine, run() can be called
else
{ this->runnable = true ;
return ;
}
}
int main(int argn, char** argv)
{ CorrelationMatrixCreatorApplication app(argn, argv) ;
return app.run() ;
}

Event Timeline