view release on metacpan or  search on metacpan

Bio/Assembly/IO/ace.pm  view on Meta::CPAN

        # Ignore it
        #(/^AS\s+(\d+)\s+(\d+)/) && do {
        #    my $nof_contigs = $1;
        #    my $nof_seq_in_contigs = $2;
        #};

        # Loading Whole Assembly tags
        /^WA\s*\{/ && do {
            my ($type,$source,$date) = split(' ',$self->_readline);
            my $extra_info = undef;
            while ($_ = $self->_readline) {
                last if (/\}/);
                $extra_info .= $_;

            }
            my $assembly_tags = join(" ","TYPE:",$type,"PROGRAM:",$source,
                "DATE:",$date,"DATA:",$extra_info);
            $assembly_tags = Bio::Annotation::SimpleValue->new(-value=>$assembly_tags);
            $assembly->annotation->add_Annotation('whole assembly',$assembly_tags);
        };

        # Loading Contig Tags (a.k.a. Bioperl features)
        /^CT\s*\{/ && do {
            my ($contigID,$type,$source,$start,$end,$date) = split(' ',$self->_readline);
            my %tags = ('source' => $source, 'creation_date' => $date);
            my $tag_type = 'extra_info';
            while ($_ = $self->_readline) {
                if (/COMMENT\s*\{/) {
                    $tag_type = 'comment';
                } elsif (/C\}/) {
                    $tag_type = 'extra_info';
                } elsif (/\}/) {
                    last;
                } else {
                    $tags{$tag_type} .= "$_";
                }
            }
            my $contig_tag = Bio::SeqFeature::Generic->new( -start   => $start,
                                                            -end     => $end,
                                                            -primary => $type,
                                                            -source  => 'ace',
                                                            -tag     => \%tags );
            my $contig = $assembly->get_contig_by_id($contigID) ||
                         $assembly->get_singlet_by_id($contigID);
            $self->throw("Cannot add feature to unknown contig '$contigID'")
              unless defined $contig;

            $contig->add_features([ $contig_tag ],1);
        };

    }
    return 1;
}


=head2 write_assembly

    Title   : write_assembly
    Usage   : $ass_io->write_assembly($assembly)
    Function: Write the assembly object in ACE compatible format. The contig IDs
              are sorted naturally if the Sort::Naturally module is present, or
              lexically otherwise. Internally, write_assembly use the
              write_contig, write_footer and write_header methods. Use these
              methods if you want more control on the writing process.
    Returns : 1 on success, 0 for error
    Args    : A Bio::Assembly::Scaffold object

=cut


=head2 write_contig

    Title   : write_contig
    Usage   : $ass_io->write_contig($contig)
    Function: Write a contig or singlet object in ACE compatible format. Quality
              scores are automatically generated if the contig does not contain
              any
    Returns : 1 on success, 0 for error
    Args    : A Bio::Assembly::Contig or Singlet object

=cut

sub write_contig {
    my ($self, @args) = @_;
    my ($contig) = $self->_rearrange([qw(CONTIG)], @args);

    # Sanity check
    if ( !$contig || !$contig->isa('Bio::Assembly::Contig') ) {
        $self->throw("Must provide a Bio::Assembly::Contig or Singlet object when calling write_contig");
    }

    # Contig consensus sequence
    my $contig_id        =  $contig->id;
    my $cons             =  $contig->get_consensus_sequence;
    my $cons_seq         =  $cons->seq;
    my $cons_len         =  $cons->length;
    my $contig_num_reads =  $contig->num_sequences;
    my $cons_strand      = ($contig->strand == -1) ? 'C' : 'U';
    my @bs_feats         = $contig->get_features_collection->get_features_by_type('_base_segments');
    my $nof_segments     = scalar @bs_feats;

    $self->_print(
        "CO $contig_id $cons_len $contig_num_reads $nof_segments $cons_strand\n".
        _formatted_seq($cons_seq, $line_width).
        "\n"
    );

    # Consensus quality scores
    $cons = $contig->get_consensus_quality;
    my $cons_qual = $cons->qual if defined $cons;
    $self->_print(
        "BQ\n".
        _formatted_qual($cons_qual, $cons_seq, $line_width, $qual_value).
        "\n"
    );
        
    # Read entries
    my @reads  = $contig->each_seq;
    for my $read (@reads) {
        my $read_id     =  $read->id;
        my $read_strand = ($read->strand == -1) ? 'C' : 'U';