diff --git a/src/cuffToTrans/Makefile b/src/cuffToTrans/Makefile deleted file mode 100644 index 6e417369a1cd5f4a15ec1c39c99b6ed405cfe11c..0000000000000000000000000000000000000000 --- a/src/cuffToTrans/Makefile +++ /dev/null @@ -1,44 +0,0 @@ -UTILITIES_DIR = ../utils/ -OBJ_DIR = ../../obj/ -BIN_DIR = ../../bin/ - -# ------------------- -# define our includes -# ------------------- -INCLUDES = -I$(UTILITIES_DIR)/bedFile/ -I$(UTILITIES_DIR)/sequenceUtilities/ -I$(UTILITIES_DIR)/lineFileUtilities/ -I$(UTILITIES_DIR)/version/ -I$(UTILITIES_DIR)/gzstream/ -I$(UTILITIES_DIR)/fileType/ - -# ---------------------------------- -# define our source and object files -# ---------------------------------- -SOURCES= cuffToTransMain.cpp cuffToTrans.cpp Fasta.cpp split.cpp -OBJECTS= $(SOURCES:.cpp=.o) -_EXT_OBJECTS=bedFile.o sequenceUtils.o lineFileUtilities.o gzstream.o fileType.o -EXT_OBJECTS=$(patsubst %,$(OBJ_DIR)/%,$(_EXT_OBJECTS)) -BUILT_OBJECTS= $(patsubst %,$(OBJ_DIR)/%,$(OBJECTS)) -PROGRAM= cuffToTrans - - -all: $(PROGRAM) - -.PHONY: all - -$(PROGRAM): $(BUILT_OBJECTS) $(EXT_OBJECTS) - @echo " * linking $(PROGRAM)" - @$(CXX) $(LDFLAGS) $(CXXFLAGS) -o $(BIN_DIR)/$@ $^ $(LIBS) - -$(BUILT_OBJECTS): $(SOURCES) - @echo " * compiling" $(*F).cpp - @$(CXX) -c -o $@ $(*F).cpp $(LDFLAGS) $(CXXFLAGS) $(INCLUDES) - -$(EXT_OBJECTS): - @$(MAKE) --no-print-directory -C $(UTILITIES_DIR)/sequenceUtilities/ - @$(MAKE) --no-print-directory -C $(UTILITIES_DIR)/lineFileUtilities/ - @$(MAKE) --no-print-directory -C $(UTILITIES_DIR)/bedFile/ - @$(MAKE) --no-print-directory -C $(UTILITIES_DIR)/gzstream/ - @$(MAKE) --no-print-directory -C $(UTILITIES_DIR)/fileType/ - -clean: - @echo "Cleaning up." - @rm -f $(OBJ_DIR)/* $(BIN_DIR)/* - -.PHONY: clean \ No newline at end of file diff --git a/src/fjoin/Makefile b/src/fjoin/Makefile deleted file mode 100644 index c98c6696430d7099b3de5cb3a43e17729b1eb8e0..0000000000000000000000000000000000000000 --- a/src/fjoin/Makefile +++ /dev/null @@ -1,42 +0,0 @@ -UTILITIES_DIR = ../utils/ -OBJ_DIR = ../../obj/ -BIN_DIR = ../../bin/ - -# ------------------- -# define our includes -# ------------------- -INCLUDES = -I$(UTILITIES_DIR)/bedFile/ -I$(UTILITIES_DIR)/lineFileUtilities/ -I$(UTILITIES_DIR)/version/ -I$(UTILITIES_DIR)/gzstream/ -I$(UTILITIES_DIR)/fileType/ - -# ---------------------------------- -# define our source and object files -# ---------------------------------- -SOURCES= fjoinMain.cpp fjoin.cpp -OBJECTS= $(SOURCES:.cpp=.o) -_EXT_OBJECTS=bedFile.o lineFileUtilities.o gzstream.o fileType.o -EXT_OBJECTS=$(patsubst %,$(OBJ_DIR)/%,$(_EXT_OBJECTS)) -BUILT_OBJECTS= $(patsubst %,$(OBJ_DIR)/%,$(OBJECTS)) -PROGRAM= fjoin - -all: $(PROGRAM) - -.PHONY: all - -$(PROGRAM): $(BUILT_OBJECTS) $(EXT_OBJECTS) - @echo " * linking $(PROGRAM)" - @$(CXX) $(LDFLAGS) $(CXXFLAGS) -o $(BIN_DIR)/$@ $^ $(LIBS) - -$(BUILT_OBJECTS): $(SOURCES) - @echo " * compiling" $(*F).cpp - @$(CXX) -c -o $@ $(*F).cpp $(LDFLAGS) $(CXXFLAGS) $(INCLUDES) - -$(EXT_OBJECTS): - @$(MAKE) --no-print-directory -C $(UTILITIES_DIR)/bedFile/ - @$(MAKE) --no-print-directory -C $(UTILITIES_DIR)/lineFileUtilities/ - @$(MAKE) --no-print-directory -C $(UTILITIES_DIR)/gzstream/ - @$(MAKE) --no-print-directory -C $(UTILITIES_DIR)/fileType/ - -clean: - @echo "Cleaning up." - @rm -f $(OBJ_DIR)/* $(BIN_DIR)/* - -.PHONY: clean diff --git a/src/fjoin/fjoin.cpp b/src/fjoin/fjoin.cpp deleted file mode 100644 index b1de416ae202304b92cb208cb267213c9f3ba41b..0000000000000000000000000000000000000000 --- a/src/fjoin/fjoin.cpp +++ /dev/null @@ -1,350 +0,0 @@ -/***************************************************************************** - intersectBed.cpp - - (c) 2009 - Aaron Quinlan - Hall Laboratory - Department of Biochemistry and Molecular Genetics - University of Virginia - aaronquinlan@gmail.com - - Licenced under the GNU General Public License 2.0 license. -******************************************************************************/ -#include "lineFileUtilities.h" -#include "fjoin.h" -#include <queue> -#include <set> - -bool leftOf(const BED &a, const BED &b); - - -bool BedIntersect::processHits(BED &a, vector<BED> &hits) { - // how many overlaps are there b/w the bed and the set of hits? - int s, e, overlapBases; - int numOverlaps = 0; - bool hitsFound = false; - int aLength = (a.end - a.start); // the length of a in b.p. - - // loop through the hits and report those that meet the user's criteria - vector<BED>::const_iterator h = hits.begin(); - vector<BED>::const_iterator hitsEnd = hits.end(); - for (; h != hitsEnd; ++h) { - s = max(a.start, h->start); - e = min(a.end, h->end); - overlapBases = (e - s); // the number of overlapping bases b/w a and b - - // is there enough overlap relative to the user's request? (default ~ 1bp) - if ( ( (float) overlapBases / (float) aLength ) >= _overlapFraction ) { - // Report the hit if the user doesn't care about reciprocal overlap between A and B. - if (_reciprocal == false) { - hitsFound = true; - numOverlaps++; - if (_printable == true) - ReportOverlapDetail(overlapBases, a, *h, s, e); - } - // we require there to be sufficient __reciprocal__ overlap - else { - int bLength = (h->end - h->start); - float bOverlap = ( (float) overlapBases / (float) bLength ); - if (bOverlap >= _overlapFraction) { - hitsFound = true; - numOverlaps++; - if (_printable == true) - ReportOverlapDetail(overlapBases, a, *h, s, e); - } - } - } - } - // report the summary of the overlaps if requested. - ReportOverlapSummary(a, numOverlaps); - // were hits found for this BED feature? - return hitsFound; -} - -/* - Constructor -*/ -BedIntersect::BedIntersect(string bedAFile, string bedBFile, bool anyHit, - bool writeA, bool writeB, bool writeOverlap, bool writeAllOverlap, - float overlapFraction, bool noHit, bool writeCount, bool forceStrand, - bool reciprocal, bool obeySplits, bool bamInput, bool bamOutput) { - - _bedAFile = bedAFile; - _bedBFile = bedBFile; - _anyHit = anyHit; - _noHit = noHit; - _writeA = writeA; - _writeB = writeB; - _writeOverlap = writeOverlap; - _writeAllOverlap = writeAllOverlap; - _writeCount = writeCount; - _overlapFraction = overlapFraction; - _forceStrand = forceStrand; - _reciprocal = reciprocal; - _obeySplits = obeySplits; - _bamInput = bamInput; - _bamOutput = bamOutput; - - if (_anyHit || _noHit || _writeCount) - _printable = false; - else - _printable = true; - - // create new BED file objects for A and B - _bedA = new BedFile(bedAFile); - _bedB = new BedFile(bedBFile); - - IntersectBed(); -} - - -/* - Destructor -*/ -BedIntersect::~BedIntersect(void) { -} - - -bool leftOf(const BED &a, const BED &b) { - return (a.end <= b.start); -} - - -void BedIntersect::ReportOverlapDetail(const int &overlapBases, const BED &a, const BED &b, - const CHRPOS &s, const CHRPOS &e) { - // default. simple intersection only - if (_writeA == false && _writeB == false && _writeOverlap == false) { - _bedA->reportBedRangeNewLine(a,s,e); - } - // -wa -wbwrite the original A and B - else if (_writeA == true && _writeB == true) { - _bedA->reportBedTab(a); - _bedB->reportBedNewLine(b); - } - // -wa write just the original A - else if (_writeA == true) { - _bedA->reportBedNewLine(a); - } - // -wb write the intersected portion of A and the original B - else if (_writeB == true) { - _bedA->reportBedRangeTab(a,s,e); - _bedB->reportBedNewLine(b); - } - // -wo write the original A and B plus the no. of overlapping bases. - else if (_writeOverlap == true) { - _bedA->reportBedTab(a); - _bedB->reportBedTab(b); - printf("%d\n", overlapBases); - } -} - - -void BedIntersect::ReportOverlapSummary(const BED &a, const int &numOverlapsFound) { - // -u just report the fact that there was >= 1 overlaps - if (_anyHit && (numOverlapsFound >= 1)) { - _bedA->reportBedNewLine(a); - } - // -c report the total number of features overlapped in B - else if (_writeCount) { - _bedA->reportBedTab(a); - printf("%d\n", numOverlapsFound); - } - // -v report iff there were no overlaps - else if (_noHit && (numOverlapsFound == 0)) { - _bedA->reportBedNewLine(a); - } - // -wao the user wants to force the reporting of 0 overlap - else if (_writeAllOverlap && (numOverlapsFound == 0)) { - _bedA->reportBedTab(a); - _bedB->reportNullBedTab(); - printf("0\n"); - } -} - - - -void BedIntersect::Scan(BED *x, vector<BED *> *windowX, BedLineStatus xStatus, - const BED &y, vector<BED *> *windowY, BedLineStatus yStatus) { - - if (xStatus != BED_VALID) { - return; - } - - std::vector<BED *>::iterator wYIter = windowY->begin(); - while (wYIter != windowY->end()) { - if (leftOf(*(*wYIter), *x) == true) { - (*wYIter)->finished = true; - wYIter = windowY->erase(wYIter); // erase auto-increments to the next position - } - else if (overlaps((*wYIter)->start, (*wYIter)->end, x->start, x->end) > 0) { - if (_lastPick == 0) { - AddHits(x, *(*wYIter)); - } - else { - AddHits(*wYIter, *x); - } - ++wYIter; // force incrementing - } - } - if (leftOf(*x,y) == false) - windowX->push_back(x); - else { - x->finished = true; - } - // dump the buffered results (if any) - FlushOutputBuffer(); -} - - -void BedIntersect::AddHits(BED *x, const BED &y) { - if (x->overlaps.empty() == true) - _outputBuffer.push(x); - x->overlaps.push_back(y); -} - - -void BedIntersect::FlushOutputBuffer(bool final) { - while (_outputBuffer.empty() == false) - { - if (final == false && _outputBuffer.front()->finished == false) - break; - - processHits(*_outputBuffer.front(), _outputBuffer.front()->overlaps); - // remove the finished BED entry from the heap - delete _outputBuffer.front(); - _outputBuffer.pop(); - } -} - - -vector<BED*>* BedIntersect::GetWindow(const string &chrom, bool isA) { - - // iterator to test if a window for a given chrom exists. - map<string, vector<BED*> >::iterator it; - - // grab the current window for A or B, depending on - // the request. if a window hasn't yet been created - // for the requested chrom, create one. - - if (isA) { - it = _windowA.find(chrom); - if (it != _windowA.end()) { - return & _windowA[chrom]; - } - else { - _windowA.insert(pair<string, vector<BED *> >(chrom, vector<BED *>())); - return & _windowA[chrom]; - } - } - else { - it = _windowB.find(chrom); - if (it != _windowB.end()) { - return & _windowB[chrom]; - } - else { - _windowB.insert(pair<string, vector<BED *> >(chrom, vector<BED *>())); - return & _windowB[chrom]; - } - } -} - - -void BedIntersect::ChromSwitch(const string &chrom) { - - vector<BED*>::iterator windowAIter = _windowA[chrom].begin(); - vector<BED*>::iterator windowAEnd = _windowA[chrom].end(); - for (; windowAIter != windowAEnd; ++windowAIter) - (*windowAIter)->finished = true; - - vector<BED*>::iterator windowBIter = _windowB[chrom].begin(); - vector<BED*>::iterator windowBEnd = _windowB[chrom].end(); - for (; windowBIter != windowBEnd; ++windowBIter) - (*windowBIter)->finished = true; - - FlushOutputBuffer(); -} - - -void BedIntersect::IntersectBed() { - - int aLineNum = 0; - int bLineNum = 0; - - // current feature from each file - BED *a, *b, *prevA, *prevB; - - // status of the current lines - BedLineStatus aStatus, bStatus; - - // open the files; get the first line from each - _bedA->Open(); - _bedB->Open(); - - prevA = NULL; - prevB = NULL; - a = new BED(); - b = new BED(); - aStatus = _bedA->GetNextBed(*a, aLineNum); - bStatus = _bedB->GetNextBed(*b, bLineNum); - - cout << a->chrom << " " << a->start << " " << a->chrom << " " << b->start << endl; - while (aStatus != BED_INVALID || bStatus != BED_INVALID) { - - if ((a->start <= b->start) && (a->chrom == b->chrom)) { - prevA = a; - _lastPick = 0; - Scan(a, GetWindow(a->chrom, true), aStatus, - *b, GetWindow(a->chrom, false), bStatus); - - a = new BED(); - aStatus = _bedA->GetNextBed(*a, aLineNum); - } - else if ((a->start > b->start) && (a->chrom == b->chrom)) { - prevB = b; - _lastPick = 1; - Scan(b, GetWindow(b->chrom, false), bStatus, - *a, GetWindow(b->chrom, true), aStatus); - - b = new BED(); - bStatus = _bedB->GetNextBed(*b, bLineNum); - } - else if (a->chrom != b->chrom) { - // A was most recently read - if (_lastPick == 0) { - prevB = b; - while (b->chrom == prevA->chrom){ - _windowB[prevA->chrom].push_back(b); - b = new BED(); - bStatus = _bedB->GetNextBed(*b, bLineNum); - } - Scan(prevA, GetWindow(prevA->chrom, true), aStatus, - *prevB, GetWindow(prevA->chrom, false), bStatus); - } - // B was most recently read - else { - prevA = a; - while (a->chrom == prevB->chrom) { - _windowA[prevB->chrom].push_back(a); - a = new BED(); - aStatus = _bedA->GetNextBed(*a, aLineNum); - } - Scan(prevB, GetWindow(prevB->chrom, false), bStatus, - *prevA, GetWindow(prevB->chrom, true), aStatus); - } - FlushOutputBuffer(true); - } - if (prevA!=NULL&&prevB!=NULL) - //cout << prevA->chrom << " " << a->chrom << " " << a->start << " " - // << prevB->chrom << " " << b->chrom << " " << b->start << "\n"; - if (aStatus == BED_INVALID) a->start = INT_MAX; - if (bStatus == BED_INVALID) b->start = INT_MAX; - } - - // clear out the final bit of staged output - FlushOutputBuffer(true); - - // close the files - _bedA->Close(); - _bedB->Close(); -} - - diff --git a/src/fjoin/fjoin.h b/src/fjoin/fjoin.h deleted file mode 100644 index c7aabd46444f78dd37e632096b543fa630ff4b53..0000000000000000000000000000000000000000 --- a/src/fjoin/fjoin.h +++ /dev/null @@ -1,114 +0,0 @@ -/***************************************************************************** - intersectBed.h - - (c) 2009 - Aaron Quinlan - Hall Laboratory - Department of Biochemistry and Molecular Genetics - University of Virginia - aaronquinlan@gmail.com - - Licenced under the GNU General Public License 2.0 license. -******************************************************************************/ -#ifndef INTERSECTBED_H -#define INTERSECTBED_H - -#include "bedFile.h" -// #include "BamReader.h" -// #include "BamWriter.h" -// #include "BamAncillary.h" -// #include "BamAux.h" -// using namespace BamTools; - - -#include <vector> -#include <queue> -#include <iostream> -#include <fstream> -#include <stdlib.h> -using namespace std; - - - -class BedIntersect { - -public: - - // constructor - BedIntersect(string bedAFile, string bedBFile, bool anyHit, - bool writeA, bool writeB, bool writeOverlap, bool writeAllOverlap, - float overlapFraction, bool noHit, bool writeCount, bool forceStrand, - bool reciprocal, bool obeySplits, bool bamInput, bool bamOutput); - - // destructor - ~BedIntersect(void); - -private: - - //------------------------------------------------ - // private attributes - //------------------------------------------------ - string _bedAFile; - string _bedBFile; - - bool _writeA; // should the original A feature be reported? - bool _writeB; // should the original B feature be reported? - bool _writeOverlap; - bool _writeAllOverlap; - - bool _forceStrand; - bool _reciprocal; - float _overlapFraction; - - bool _anyHit; - bool _noHit; - bool _writeCount; // do we want a count of the number of overlaps in B? - bool _obeySplits; - bool _bamInput; - bool _bamOutput; - - bool _printable; - - queue<BED*> _outputBuffer; - bool _lastPick; - - map<string, vector<BED*> > _windowA; - map<string, vector<BED*> > _windowB; - - // instance of a bed file class. - BedFile *_bedA, *_bedB; - - //------------------------------------------------ - // private methods - //------------------------------------------------ - void IntersectBed(istream &bedInput); - - void Scan(BED *x, vector<BED *> *windowX, BedLineStatus xStatus, - const BED &y, vector<BED *> *windowY, BedLineStatus yStatus); - - void AddHits(BED *x, const BED &y); - - void FlushOutputBuffer(bool final = false); - - vector<BED*>* GetWindow(const string &chrom, bool isA); - - void ChromSwitch(const string &chrom); - - void IntersectBed(); - - void IntersectBam(string bamFile); - - bool processHits(BED &a, vector<BED> &hits); - - bool FindOverlaps(const BED &a, vector<BED> &hits); - - bool FindOneOrMoreOverlap(const BED &a); - - void ReportOverlapDetail(const int &overlapBases, const BED &a, const BED &b, - const CHRPOS &s, const CHRPOS &e); - void ReportOverlapSummary(const BED &a, const int &numOverlapsFound); - - void ReportHits(set<BED> &A, set<BED> &B); - -}; - -#endif /* INTERSECTBED_H */ diff --git a/src/fjoin/fjoinMain.cpp b/src/fjoin/fjoinMain.cpp deleted file mode 100644 index 48a2ad982a1c6c6074ce53879beb1287cdfb2ef9..0000000000000000000000000000000000000000 --- a/src/fjoin/fjoinMain.cpp +++ /dev/null @@ -1,271 +0,0 @@ -/***************************************************************************** - intersectMain.cpp - - (c) 2009 - Aaron Quinlan - Hall Laboratory - Department of Biochemistry and Molecular Genetics - University of Virginia - aaronquinlan@gmail.com - - Licenced under the GNU General Public License 2.0 license. -******************************************************************************/ -#include "fjoin.h" -#include "version.h" - -using namespace std; - -// define our program name -#define PROGRAM_NAME "fjoin" - - -// define our parameter checking macro -#define PARAMETER_CHECK(param, paramLen, actualLen) (strncmp(argv[i], param, min(actualLen, paramLen))== 0) && (actualLen == paramLen) - -// function declarations -void ShowHelp(void); - -int main(int argc, char* argv[]) { - - // our configuration variables - bool showHelp = false; - - // input files - string bedAFile; - string bedBFile; - - // input arguments - float overlapFraction = 1E-9; - - bool haveBedA = false; - bool haveBedB = false; - bool noHit = false; - bool anyHit = false; - bool writeA = false; - bool writeB = false; - bool writeCount = false; - bool writeOverlap = false; - bool writeAllOverlap = false; - bool haveFraction = false; - bool reciprocalFraction = false; - bool forceStrand = false; - bool obeySplits = false; - bool inputIsBam = false; - bool outputIsBam = true; - - // check to see if we should print out some help - if(argc <= 1) showHelp = true; - - for(int i = 1; i < argc; i++) { - int parameterLength = (int)strlen(argv[i]); - - if((PARAMETER_CHECK("-h", 2, parameterLength)) || - (PARAMETER_CHECK("--help", 5, parameterLength))) { - showHelp = true; - } - } - - if(showHelp) ShowHelp(); - - // do some parsing (all of these parameters require 2 strings) - for(int i = 1; i < argc; i++) { - - int parameterLength = (int)strlen(argv[i]); - - if(PARAMETER_CHECK("-a", 2, parameterLength)) { - if ((i+1) < argc) { - haveBedA = true; - outputIsBam = false; - bedAFile = argv[i + 1]; - i++; - } - } - else if(PARAMETER_CHECK("-abam", 5, parameterLength)) { - if ((i+1) < argc) { - haveBedA = true; - inputIsBam = true; - bedAFile = argv[i + 1]; - i++; - } - } - else if(PARAMETER_CHECK("-b", 2, parameterLength)) { - if ((i+1) < argc) { - haveBedB = true; - bedBFile = argv[i + 1]; - i++; - } - } - else if(PARAMETER_CHECK("-bed", 4, parameterLength)) { - outputIsBam = false; - } - else if(PARAMETER_CHECK("-u", 2, parameterLength)) { - anyHit = true; - } - else if(PARAMETER_CHECK("-f", 2, parameterLength)) { - if ((i+1) < argc) { - haveFraction = true; - overlapFraction = atof(argv[i + 1]); - i++; - } - } - else if(PARAMETER_CHECK("-wa", 3, parameterLength)) { - writeA = true; - } - else if(PARAMETER_CHECK("-wb", 3, parameterLength)) { - writeB = true; - } - else if(PARAMETER_CHECK("-wo", 3, parameterLength)) { - writeOverlap = true; - } - else if(PARAMETER_CHECK("-wao", 4, parameterLength)) { - writeAllOverlap = true; - writeOverlap = true; - } - else if(PARAMETER_CHECK("-c", 2, parameterLength)) { - writeCount = true; - } - else if(PARAMETER_CHECK("-r", 2, parameterLength)) { - reciprocalFraction = true; - } - else if (PARAMETER_CHECK("-v", 2, parameterLength)) { - noHit = true; - } - else if (PARAMETER_CHECK("-s", 2, parameterLength)) { - forceStrand = true; - } - else if (PARAMETER_CHECK("-split", 6, parameterLength)) { - obeySplits = true; - } - else { - cerr << endl << "*****ERROR: Unrecognized parameter: " << argv[i] << " *****" << endl << endl; - showHelp = true; - } - } - - // make sure we have both input files - if (!haveBedA || !haveBedB) { - cerr << endl << "*****" << endl << "*****ERROR: Need -a and -b files. " << endl << "*****" << endl; - showHelp = true; - } - - if (anyHit && noHit) { - cerr << endl << "*****" << endl << "*****ERROR: Request either -u OR -v, not both." << endl << "*****" << endl; - showHelp = true; - } - - if (writeB && writeCount) { - cerr << endl << "*****" << endl << "*****ERROR: Request either -wb OR -c, not both." << endl << "*****" << endl; - showHelp = true; - } - - if (writeCount && writeOverlap) { - cerr << endl << "*****" << endl << "*****ERROR: Request either -wb OR -wo, not both." << endl << "*****" << endl; - showHelp = true; - } - - if (writeA && writeOverlap) { - cerr << endl << "*****" << endl << "*****ERROR: Request either -wa OR -wo, not both." << endl << "*****" << endl; - showHelp = true; - } - - if (writeB && writeOverlap) { - cerr << endl << "*****" << endl << "*****ERROR: Request either -wb OR -wo, not both." << endl << "*****" << endl; - showHelp = true; - } - - if (reciprocalFraction && !haveFraction) { - cerr << endl << "*****" << endl << "*****ERROR: If using -r, you need to define -f." << endl << "*****" << endl; - showHelp = true; - } - - if (anyHit && writeCount) { - cerr << endl << "*****" << endl << "*****ERROR: Request either -u OR -c, not both." << endl << "*****" << endl; - showHelp = true; - } - - if (anyHit && writeB) { - cerr << endl << "*****" << endl << "*****ERROR: Request either -u OR -wb, not both." << endl << "*****" << endl; - showHelp = true; - } - - if (anyHit && writeOverlap) { - cerr << endl << "*****" << endl << "*****ERROR: Request either -u OR -wo, not both." << endl << "*****" << endl; - showHelp = true; - } - - - if (!showHelp) { - - BedIntersect *bi = new BedIntersect(bedAFile, bedBFile, anyHit, writeA, writeB, writeOverlap, - writeAllOverlap, overlapFraction, noHit, writeCount, forceStrand, - reciprocalFraction, obeySplits, inputIsBam, outputIsBam); - delete bi; - return 0; - } - else { - ShowHelp(); - } -} - -void ShowHelp(void) { - - cerr << endl << "Program: " << PROGRAM_NAME << " (v" << VERSION << ")" << endl; - - cerr << "Author: Aaron Quinlan (aaronquinlan@gmail.com)" << endl; - - cerr << "Summary: Report overlaps between two feature files." << endl << endl; - - cerr << "Usage: " << PROGRAM_NAME << " [OPTIONS] -a <bed/gff/vcf> -b <bed/gff/vcf>" << endl << endl; - - cerr << "Options: " << endl; - - cerr << "\t-abam\t" << "The A input file is in BAM format. Output will be BAM as well." << endl << endl; - - cerr << "\t-bed\t" << "When using BAM input (-abam), write output as BED. The default" << endl; - cerr << "\t\tis to write output in BAM when using -abam." << endl << endl; - - cerr << "\t-wa\t" << "Write the original entry in A for each overlap." << endl << endl; - - cerr << "\t-wb\t" << "Write the original entry in B for each overlap." << endl; - cerr << "\t\t- Useful for knowing _what_ A overlaps. Restricted by -f and -r." << endl << endl; - - cerr << "\t-wo\t" << "Write the original A and B entries plus the number of base" << endl; - cerr << "\t\tpairs of overlap between the two features." << endl; - cerr << "\t\t- Overlaps restricted by -f and -r." << endl; - cerr << "\t\t Only A features with overlap are reported." << endl << endl; - - cerr << "\t-wao\t" << "Write the original A and B entries plus the number of base" << endl; - cerr << "\t\tpairs of overlap between the two features." << endl; - cerr << "\t\t- Overlapping features restricted by -f and -r." << endl; - cerr << "\t\t However, A features w/o overlap are also reported" << endl; - cerr << "\t\t with a NULL B feature and overlap = 0." << endl << endl; - - cerr << "\t-u\t" << "Write the original A entry _once_ if _any_ overlaps found in B." << endl; - cerr << "\t\t- In other words, just report the fact >=1 hit was found." << endl; - cerr << "\t\t- Overlaps restricted by -f and -r." << endl << endl; - - cerr << "\t-c\t" << "For each entry in A, report the number of overlaps with B." << endl; - cerr << "\t\t- Reports 0 for A entries that have no overlap with B." << endl; - cerr << "\t\t- Overlaps restricted by -f and -r." << endl << endl; - - cerr << "\t-v\t" << "Only report those entries in A that have _no overlaps_ with B." << endl; - cerr << "\t\t- Similar to \"grep -v\" (an homage)." << endl << endl; - - cerr << "\t-f\t" << "Minimum overlap required as a fraction of A." << endl; - cerr << "\t\t- Default is 1E-9 (i.e., 1bp)." << endl; - cerr << "\t\t- FLOAT (e.g. 0.50)" << endl << endl; - - cerr << "\t-r\t" << "Require that the fraction overlap be reciprocal for A and B." << endl; - cerr << "\t\t- In other words, if -f is 0.90 and -r is used, this requires" << endl; - cerr << "\t\t that B overlap 90% of A and A _also_ overlaps 90% of B." << endl << endl; - - cerr << "\t-s\t" << "Force strandedness. That is, only report hits in B that" << endl; - cerr << "\t\toverlap A on the same strand." << endl; - cerr << "\t\t- By default, overlaps are reported without respect to strand." << endl << endl; - - cerr << "\t-split\t" << "Treat \"split\" BAM or BED12 entries as distinct BED intervals." << endl << endl; - - - // end the program here - exit(1); - -}