[Bioperl-l] Fastq - is it flush, or not?

[Joel Martin]

i'd use fastx toolkit( fastx_reverse_complement),  or at least you
could validate the bug with that.  and your regex can exclude valid
sequences.  this find a problematic line.

#!/usr/bin/perl -w
use strict;
my $slen;
while ( <> ) {
  $_=<>;
  $slen=length;
  <>;$_=<>;
  if( length != $slen ) {
    print "quality length != ", $slen - 1, " at line $.\nbad->$_";
  }
}


2011/7/11 Nestor Zaburannyi <nestor at linuxmail.org>:
> Chris,
>
> cat myfile.fastq | grep -v ^@ | grep -v ^+$ | awk '{ print $0 " = " length($0) }' | grep -v 50
>
> gives empty string and that means every sequence AND every quality string has exactly 50 characters. So, no clipping present, that's for sure. Any other ideas anyone?
>
> Sincerely yours
> Nestor
>
>
>
> Monday, July 11, 2011, 2:17:07 PM, you wrote:
>
>> The error appears to be in revcom(), so my guess is something with clipping of the quality data?  Hard so say for sure, though, w/o seeing the data and testing it.
>
>> chris
>
>> On Jul 11, 2011, at 3:47 AM, Nestor Zaburannyi wrote:
>
>>> Dear All.
>
>>> I need to reverse-complement a huge *.fastq file. I try to
>
>>> while (my $seq = $in->next_seq)
>>>    {
>>>    $out->write_seq($seq->revcom);
>>>    }
>
>>> However, it stops with:
>
>>> ------------- EXCEPTION -------------
>>> MSG: Can not get a reverse complement. The object is not flush.
>>> STACK Bio::Seq::Meta::Array::revcom /usr/share/perl5/Bio/Seq/Meta/Array.pm:648
>
>
>>> Adding
>
>>> $seq->force_flush('1');
>
>>> helps, but is it necessary and proper way? Sequences seem to be "flush".
>
>>> Sincerely yours
>>> Nestor
>
>>> _______________________________________________
>>> Bioperl-l mailing list
>>> Bioperl-l at lists.open-bio.org
>>> http://lists.open-bio.org/mailman/listinfo/bioperl-l
>
>
>
>
> --
> З повагою,
>  Nestor                            mailto:nestor at linuxmail.org
>
>
> _______________________________________________
> Bioperl-l mailing list
> Bioperl-l at lists.open-bio.org
> http://lists.open-bio.org/mailman/listinfo/bioperl-l