#!/usr/bin/perl -w
#Filename:
#Author: Tian Dongmei
#Email: tiandm@big.ac.cn
#Date: 2009-05-06
#Modified:
#Description: ɾ³ýmatched reads 
my $version=1.00;

use strict;
use Getopt::Long;
use File::Basename;

my %opts;
GetOptions(\%opts,"i=s","outdir:s","column1:i","mark1=s","depth1:i","depth2:i","column2:i","mark2=s","h");
if (!(defined $opts{i} and defined $opts{mark1} and defined $opts{mark2}) || defined $opts{h}) { #necessary arguments
&usage;
}

my $filein=$opts{'i'};
my $outputdir=defined $opts{'outdir'} ? $opts{'outdir'} : "./";
unless ($outputdir=~/\/$/) {$outputdir .="/";}
my $column1=defined $opts{column1} ? $opts{column1} : 3;
my $column2=defined $opts{column2} ? $opts{column2} : 4;
my $mark1=$opts{mark1};
my $mark2=$opts{mark2};
my $fileout=$outputdir."degseq.R";

open OUT,">$fileout"; #output file  
#my ($name,$dir);
#$name=basename($filein);
print OUT "library(DEGseq)\n";
print OUT "geneExpFile <- system.file(package=\"DEGseq\")\n";
print OUT "geneExpFile<-file.path(\"$filein\")\n";
print OUT "layout(matrix(c(1,2,3,4,5,6), 3, 2, byrow=TRUE))\npar(mar=c(2, 2, 2,2))\n";
print OUT "outputdir<-file.path(\"$outputdir\")\n";
print OUT "geneExpMatrix1 <- readGeneExp(file=geneExpFile, geneCol=1, valCol=c($column1))\n";
print OUT "geneExpMatrix2 <- readGeneExp(file=geneExpFile, geneCol=1, valCol=c($column2))\n";
if(defined $opts{'depth1'} && defined $opts{'depth2'}){
print OUT "DEGexp(geneExpMatrix1=geneExpMatrix1, geneCol1=1, expCol1=c(2), groupLabel1=\"$mark1\",geneExpMatrix2=geneExpMatrix2, geneCol2=1, expCol2=c(2), groupLabel2=\"$mark2\",depth1=$opts{depth1},depth2=$opts{depth2},outputDir=outputdir,method=\"MARS\")\n";
}
else{
print OUT "DEGexp(geneExpMatrix1=geneExpMatrix1, geneCol1=1, expCol1=c(2), groupLabel1=\"$mark1\",geneExpMatrix2=geneExpMatrix2, geneCol2=1, expCol2=c(2), groupLabel2=\"$mark2\",outputDir=outputdir,method=\"MARS\")\n";
}
close OUT;


my $re=`which R`;
print $re;
system("R CMD BATCH $fileout");

wait;

sub usage{
print <<"USAGE";
Version $version
Usage:
$0 -i -outdir -column1 -mark1 -column2 -mark2 -depth1 -depth2
options:
-i input file
-outdir output file dir
-column1 the first column for DEGseq
-mark1 the name of the column1
-depth1 depth for the first file,use for normalize
-column2 the second column for DEGseq
-mark2 the name of the column2
-depth2 depth for the second file,use for normalize

-h help
USAGE
exit(1);
}