/* 
get positional bias
NOTE: INCOMPLETE
*/
#include <iostream>
#include <fstream>
#include "auxiliaryio.h"
#include "aligh.h"

using namespace std;


Annotation ALLANNO;
GeneRange GENE_RANGE;
string samfile="-";

void printmsg(){
  cerr<<"Usage: getposbias [REFERENCE BED FILE] [SAM FILE]\n";
}

int pareopt(int argc, char* argv[]){

  if(argc<3){
    printmsg();
    return -1;
  }
  //parsing annotation file
  string annobedfile(argv[argc-2]);
  if(parseannobed(annobedfile,ALLANNO)==-1) return -1;
  vector<string> allchr=ALLANNO.getChrom();
  for(int ai=0;ai<allchr.size();ai++){
    vector<range_t> allr=ALLANNO.getCluster(allchr[ai]);
    for(int ri=0;ri<allr.size();ri++)
      GENE_RANGE.push_back(allchr[ai],allr[ri]);

  samfile=string(argv[argc-1]);

  return 0;
}

int completelyoverlap(vi& sn, vi& en, vi& refs, vi& refe,int readlen){
  long totallen=0;
  for(int i=0;i<sn.size();i++){
    for(int j=0;j<refs.size();j++){
      long d=getoverlaplength(sn[i],en[i],refs[j],refe[j]);
      totallen+=d;
      if(totallen>=readlen) return totallen;
      if(refs[j]>en[i])break;
   
    }
  }
  return totallen;
}



int processreadgroups(ReadGroup& rd,string chr, range_t currentrange){
  if(rd.size()==0)return 0;
  //get the annotation
  ReadGroup ref=ALLANNO.getReadGroup(chr,currentrange);
  if(ref.size()==0)return 0;
  
  //go through the reads
  mi &rs=rd.s();
  mi &re=rd.e();
  
  mi &refs=ref.s();
  mi &refe=ref.e();
  mi &refd=ref.getDirection();
  
  for(int ii=0;ii<rs.size();ii++){
    vi& sn=rs[ii];
    vi& en=re[ii];
    int readlen=rd.getReadLen(ii);
    for(int jj=0;jj<refs.size();jj++){
      vi& crs=refs[jj];
      vi& cre=refe[jj];
      int matchlen=completelyoverlap(sn,en,crs,cre);
      if(matchlen==readlen){
        //completely match
      }
    }
  }
  
  return 0;
}

int gothroughsam(){
  ifstream ifs;
  if(samfile!="-"){
    ifs.open(samfile.c_str());
    if(!ifs.is_open()){
      cerr<<"Error opening sam file "<<samfile<<endl;
      return -1;
    }
  }
  istream & alli=(samfile=="-"?cin:ifs);
  
  ReadGroup rd;
  Align calign;
  range_t currentrange(0,0);
  string prevchr="";
  int linecount=0;
  //main loop
  while(true){
    string oneline;
    getline(alli,oneline);
    if(alli.eof())break;
    linecount++;
    //parse
    if(calign.parse(oneline)==-1)
      cerr<<"Error parsing SAM records at line "<<linecount<<endl;
   
    if(calign.rname!=prevchr || calign.pos>currentrange.second){
      processreadgroups(rd,currentrange);
      rd.clear();
      if(GENE_RANGE.getNext(calign.rname,prevchr,currentrange)==-1)
         break;
    }
    
    rd.add(calign);
    prevchr=calign.rname;
    
  }
  if(ifs.is_open()){
    ifs.close();
  }
  return 0;
}


int main(int argc, char* argv[]){
  if(parseopt(argc,argv)==-1)return 0;
  gothroughsam();
  return 0;
}