websoft/bin/concise2blc

   1 #!/usr/bin/perl
   2
   3 #
   4 # Generates a blc file from the concise file, reads in the sequence
   5 # alignments those records in @data
   6 #
   7
   8 #
   9 #  CC 22/06/07 - bug fix in 'conversion' of coil-coiled predictions
  10 #  CC 22/05/06 - bug fix in final foreach loop
  11 #
  12
  13 use strict;
  14 use warnings;
  15
  16 if ($ARGV[0]) { open(IN, "<$ARGV[0]") or die($!); }
  17 else { open(IN, "<-") or die($!); }
  18
  19 # The names of the records from the concise file that we want
  20 my @data = ("JNETALIGN", "JNETHMM", "jnetpred", "JNETPSSM", "JNETCONF", "JNETSOL25", "JNETSOL5", "JNETSOL0", "Lupas_21", "Lupas_14", "Lupas_28");
  21
  22 my %convert = ("JNETALIGN" => "jalign",
  23         "JNETHMM" => "jhmm",
  24         "jnetpred" => "jnet",
  25         "JNETPSSM" => "jpssm",
  26         "JNETCONF" => "conf",
  27         "JNETSOL25" => "sol25",
  28         "JNETSOL5" => "sol5",
  29         "JNETSOL0" => "sol0",
  30         "Lupas_21" => "lupas_21",
  31         "Lupas_14" => "lupas_14",
  32         "Lupas_28" => "lupas_28",
  33         );
  34
  35 my (@seq, %seq, @pred, %pred);   # CC - these could be replaced with a tied hash, for future edits.
  36
  37 #
  38 # Read in the concise file and extract the data
  39 #
  40
  41 while (<IN>) {
  42         if (/^\n/) { next; }
  43         my ($id, $seq) = split(":", $_);
  44         if (!$id || !$seq) { next; }            # Check we have proper values
  45         chomp($seq);
  46         $seq =~ s/,//g;
  47         if ($id =~ /align\d*;/) {
  48                 @_ = split(";", $id);           # Then its an alignment
  49                 push @seq, $_[1];
  50                 $seq{$_[1]} = $seq;
  51                 }
  52         foreach (@data) {
  53                 if ($id eq $_) {
  54                         push @pred, $_;
  55                         $pred{$_} = $seq;
  56                         }
  57                 }
  58         }
  59 close(IN);
  60
  61 my @blc;
  62 my $gap;
  63 my $seq_len = length($seq{$seq[0]}) - 1;
  64 foreach (0..$seq_len) { $gap .= " "; }
  65
  66 push @blc, $gap;
  67
  68 foreach (@seq) {
  69         $seq{$_} =~ s/\.|-/ /g;
  70         push @blc, $seq{$_};
  71         }
  72 push @blc, $gap;
  73 push @blc, $gap;
  74 push @blc, $gap;
  75 foreach (@pred) {
  76         ## CC 22/06/07 - Convert 8-state secondary structure types to 3-state types unless
  77         ## it's a coil-coil prediction where 'C' is a valid prediction.
  78         ## I think this is a throw-back from when Jpred was a consensus server?
  79    if ($_ =~ /Lupas/i) {
  80       $pred{$_} =~ s/-/ /g;
  81    } else {
  82       $pred{$_} =~ s/[TCYWXZ_SI\?-]/ /g;
  83       $pred{$_} =~ s/G/H/g;
  84    }
  85    push @blc, $pred{$_};
  86 }
  87
  88 foreach (@seq) {
  89         print ">$_\n";
  90         }
  91 print ">\n>\n>\n";
  92 foreach (@pred) {
  93         print ">$convert{$_}\n";
  94         }
  95 print " * iteration 1\n";
  96 #print " $seq_len\n";
  97
  98 foreach my $i (0..$seq_len) {
  99         foreach (@blc) {
 100       ### CC 22/05/06 - below has been changed due to bug where if
 101       # statement return false when $b == 0. This should actually
 102       # be true when parsing the JNETCONF string of integers.
 103       #
 104       # Now we test whether the variable is defined rather than true.
 105       # NB: not sure what do if the variable is not defined?
 106
 107                 #if (my $b = substr($_, $i, 1)) { print $b; }
 108
 109       my $value = substr($_, $i, 1);
 110       print $value if defined($value);
 111                 }
 112         print "\n";
 113         }
 114 print " *";