use Data::Dumper;


if(-d "org") {
   system("rm -fr org/");   
}
system("mkdir org");

my %rf2type;
open GENOMESEQ, "family.txt" || die "$!";
while(<GENOMESEQ>) {
   chomp;
   my @temp = split(/\t/);
   $temp[18] =~ s/Gene; //;
   $temp[18] =~ s/Cis-reg; //;
   if($temp[18] =~ /snoRNA/) {
	  $temp[18] = "snoRNA";    
   } elsif($temp[18] =~ /snRNA/) {
	  $temp[18] = "snRNA";    
   } elsif($temp[18] =~ /sRNA/) {
	  $temp[18] = "sRNA";    
   } 
   $temp[18] =~ s/;//;
   
   $rf2type{$temp[0]} = $temp[18];
}
close GENOMESEQ;

my %acc2name;
my %acc2annotation;
open GENOMESEQ, "acc2name.txt" || die "$!";
while(<GENOMESEQ>) {
   chomp;
   my @temp = split(/\t/);
   $acc2name{$temp[1]} = $temp[0];
   $acc2annotation{$temp[1]} = $temp[2];
}
close GENOMESEQ;



my %id2acc;
open GENOMESEQ, "Rfam.full_region" || die "$!";
while(<GENOMESEQ>) {
   chomp;
   my @temp = split(/\t/);
   my ($start,$end);
   $start = $temp[2]; 
      $end = $temp[3]; 	  
   
   $id2acc{$temp[1]."/".$start."-".$end} = $temp[0];
}
close GENOMESEQ;


my %txid2org_mir;
open GENOMESEQ, "/var/www/html/sRNAtools/db/mirbase/v22.1/organisms.txt" || die "$!";
while(<GENOMESEQ>) {
   chomp;
   next if (/^#/);
   my @temp = split(/\t/);
   $txid2org_mir{$temp[4]} = $temp[0]."_".$temp[4];
}
close GENOMESEQ;


my %id2org;
my %ge2txid;
open GENOMESEQ, "rfamseq.txt" || die "$!";
while(<GENOMESEQ>) {
   chomp;
   my @temp = split(/\t/);
   my ($gename,$txonid);
   $gename = $temp[0]; 
   $txonid = $temp[3];
   my $orgabb = "";
   if($temp[6] !~ /uncultured/ and $temp[6] !~ / sp. /) {
	if($temp[6] =~ /^(\w)\w+\s+(\w\w)\w+/) {
	  $orgabb = lc($1).$2;
	  $id2org{$gename} = $orgabb."_".$txonid;
	}
   }
   $ge2txid{$gename} = $txonid;
}
close GENOMESEQ;

open GENOMESEQ, "all.rfam.fa" || die "$!";
$/=">";<GENOMESEQ>;$/="\n";
while(<GENOMESEQ>) {
   chomp;
   my $seqname = $_;  
   $/=">";
   my $seq=<GENOMESEQ>;
   chomp($seq);
   my @temp = split(/ /,$seqname);
   my @temp2 = split(/\//,$temp[0]);
   my $org="";
   my $txid = $ge2txid{$temp2[0]};
   if(exists $txid2org_mir{$txid}) {
	   $org = $txid2org_mir{$txid};
   } elsif(exists $id2org{$temp2[0]}) {
	   $org = $id2org{$temp2[0]};
   }
   if($org ne "") {
	open OUT, ">>org/$org.fa" || die "$!";
	my $pfid = $id2acc{$temp[0]};
	if(exists $rf2type{$pfid} and $rf2type{$pfid} !~ /miRNA|lncRNA|tRNA/) {
		print OUT ">".$temp[0]."|".$acc2name{$pfid}."|".$acc2annotation{$pfid}."|".$rf2type{$pfid}."\n".$seq;
	}
	close OUT;
   }
   $/="\n";
}
close GENOMESEQ;
