view release on metacpan or  search on metacpan

lib/Bio/App/SELEX/Stockholm.pm  view on Meta::CPAN

=head2 to_file

    $stock->to_file ($filename);

Writes a Stockholm object to a file.

=cut

# to_file method
sub to_file {
    my ( $self, $filename, $maxcols ) = @_;

    local *FILE;
    open FILE, ">$filename"
      or croak "Couldn't open '$filename' for writing: $!";
    print FILE $self->to_string($maxcols);
    close FILE or croak "Couldn't close '$filename': $!";

    return $filename;
}

=head2 to_string

    print $stock->to_string ($maxcols)
    print $stock->to_string ($maxcols, ARG1=>VAL1, ARG2=>VAL2 ...)
    print $stock->to_string (MAXCOLS=>$maxcols, ARG1=>VAL1, ARG2=>VAL2 ...)

Returns the object as a Stockholm-formatted string.

ARGs can include...
        MAXCOLS    -- limit maximum number of columns (can also be specified as a single scalar arg)
        NOSEQDATA  -- don\'t print any sequence data (can be used this to compress output)

=cut

# to_string method
# Usage:
#        $stock->to_string ($maxcols)
#        $stock->to_string ($maxcols, ARG1=>VAL1, ARG2=>VAL2 ...)
#        $stock->to_string (MAXCOLS=>$maxcols, ARG1=>VAL1, ARG2=>VAL2 ...)
# Args:
#        MAXCOLS    -- limit maximum number of columns (can also be specified as a single scalar arg)
#        NOSEQDATA  -- don't print any sequence data ("windowlicker.pl" uses this to compress output)
sub to_string {
    my ( $self, @args ) = @_;
    my ( %args, $maxcols );
    if ( @args % 2 == 1 ) {
        $maxcols = shift @args;
        %args    = @args;
    }
    else {
        %args    = @args;
        $maxcols = $args{'MAXCOLS'};
    }
    $maxcols = 80 unless defined $maxcols;    # default 80-column output

    # init output array
    my @out;
    push @out, "# STOCKHOLM 1.0";

    # determine alignment columns, legend columns & effective columns per line
    my $acols   = $self->columns;
    my $lcols   = $self->lcols;
    my $colstep = $maxcols < 1 ? $acols : $maxcols - $lcols - 1;
    $colstep = $maxcols
      if $colstep < 1;    # protect against negative and 0 colstep...

    # GF lines
    # check for gfOrder (insane, fragile Stockholm line ordering strikes again)
    if ( @{ $self->gfOrder } == map { (@$_) } values %{ $self->gf } )
    {                     # gfOrder same number of lines as #=GF block?
        my %gfCursor = map ( ( $_ => 0 ), keys %{ $self->gf } );
        foreach my $feature ( @{ $self->gfOrder } ) {
            push @out,
              $self->prettify(
                $lcols,
                "#=GF $feature",
                $self->gf_($feature)->[ $gfCursor{$feature}++ ]
              );
        }
    }
    else {
        @{ $self->gfOrder } = ();    # gfOrder is useless, so flush it
        foreach my $feature ( sort { $a cmp $b } keys %{ $self->gf } ) {
            push @out,
              $self->prettify(
                $lcols,
                "#=GF $feature",
                @{ $self->gf_($feature) }
              );
        }
    }

    # GS lines
    my @gs_seqname = @{ $self->seqname };
    my %gs_seqname_hash = map ( ( $_ => 1 ),
        grep ( !exists $self->seqdata->{$_},
            map ( keys( %{ $self->gs_($_) } ), keys %{ $self->gs } ) ) );
    push @gs_seqname, keys %gs_seqname_hash;
    foreach my $feature ( sort { $a cmp $b } keys %{ $self->gs } ) {
        my $hash = $self->gs_($feature);
        foreach my $seqname ( grep ( exists( $hash->{$_} ), @gs_seqname ) ) {
            push @out,
              $self->prettify(
                $lcols,
                "#=GS $seqname $feature",
                @{ $hash->{$seqname} }
              );
        }
    }

    my @gcfeat          = sort { $a cmp $b } keys %{ $self->gc };
    my @gr_seqname      = @{ $self->seqname };
    my %gr_seqname_hash = map ( ( $_ => 1 ),
        grep ( !exists $self->seqdata->{$_},
            map ( keys( %{ $self->gr_($_) } ), keys %{ $self->gr } ) ) );
    push @gr_seqname, keys %gr_seqname_hash;

    # Loop over columns
    for ( my $col = 0 ; $col < $acols ; $col += $colstep ) {

lib/Bio/App/SELEX/Stockholm.pm  view on Meta::CPAN

=head2 gr

    my $gr = $stock->gr->{FEATURE}->{SEQNAME};
    my $gr = $stock->gr_FEATURE->{SEQNAME};

Returns the line beginning '#=GR SEQNAME FEATURE ...'

=cut


# catch all methods by default
sub AUTOLOAD {
    my ($self, @args) = @_;
    my $sub = our $AUTOLOAD;
    $sub =~ s/.*:://;

    # check for DESTROY
    return if $sub eq "DESTROY";

    # check for GF, GC, GS, GR tag_feature accessors
    # e.g. $self->GF_ID
    # can also use $self->GF_('ID')
    if ($sub =~ /^(gf|gc|gs|gr)_(\S*)$/i) {
	my ($tag, $feature) = ($1, $2);
	$tag = lc $tag;
	$feature = shift(@args) unless length $feature;
	my $hash = $self->{$tag};
	if (@args) {
	    $hash->{$feature} = shift(@args);   # TODO: check ref-type of arg
	}
	cluck "Warning: ignoring extra arguments to ${tag}_$feature"
	    if @args > 0;
	if (!defined $hash->{$feature}) {
	    $hash->{$feature} = 
		$tag eq "gf" ? [] :
		$tag eq "gc" ? "" :
		$tag eq "gs" ? {} :
		$tag eq "gr" ? {} :
		croak "Unreachable";
	}
	return $hash->{$feature};
    }

    # check for ordinary accessors
    if (exists $self->{$sub}) {
	croak "Usage: $sub() or $sub(newValue)" if @args > 1;
	return
	    @args
	    ? $self->{$sub} = $args[0]
	    : $self->{$sub};
    }

    # croak
    croak "Unsupported method: $sub";
}


# pretty print line(s)
sub prettify {
    my ($self, $lcols, $legend, @data) = @_;
    # This horribly inefficient/redundant series of transformations comes out with something I like (IH, 7/24/07)
    # Trim it down? pah! Like I have nothing better to do
    $legend = sprintf ("% ${lcols}s", $legend);
    $legend =~ s/^(\s+)(\#=\S\S)(\s+\S+)$/$2$1$3/;
    $legend =~ s/^(\s+)(\#=\S\S\s+\S+)(\s+\S+)$/$2$1$3/;
    $legend =~ s/^(\s\s\s\s\s)(\s+)([^\#]\S+)/$1$3$2/;
    return map ("$legend $_", @data);
    # (pukes into cold coffee cup)
}

# legend width (subtract this from maxcols-1 to get number of columns available for sequence display)
sub lcols {
    my ($self) = @_;
    my $lcols = max ($self->maxNameLen,
		     map(length("#=GF $_"), keys(%{$self->gf})),
		     map(length("#=GC $_"), keys(%{$self->gc})));
    while (my ($gr_key, $gr_hash) = each %{$self->gr}) {
	$lcols = max ($lcols, map(length("#=GR $gr_key $_"), keys(%$gr_hash)));
    }
    while (my ($gs_key, $gs_hash) = each %{$self->gs}) {
	$lcols = max ($lcols, map(length("#=GS $gs_key $_"), keys(%$gs_hash)));
    }
    return $lcols;
}

# safer version of substr (doesn't complain if startpoint outside of string)
sub safe_substr {
    my ($string, $start, $len) = @_;
    return "" if !defined($string) || $start > length $string;
    if (defined $len) {
      return substr ($string, $start, $len);
    } else {
      return substr ($string, $start);  # print from $start to end if no length specified
    }
}

=head2 add_gf

    $stock->add_gf (FEATURE, $line)
    $stock->add_gf (FEATURE, $line1, $line2, ...)
    $stock->add_gf (FEATURE, @lines)
    $stock->add_gf (FEATURE, "$line1\n$line2\n$...")

Add '#=GF FEATURE' annotation, preserving the order of the '#=GF' lines
(damn those crazy context-sensitive Stockholm semantics!)

Each list entry gets printed on its own line.

=cut

# Add #=GF annotation, preserving the order of the #=GF lines
# (damn those crazy context-sensitive Stockholm semantics!)
#
# Usage:
#        $stock->add_gf ($gfFeature, $line);
#        $stock->add_gf ($gfFeature, $line1, $line2, ...);
#        $stock->add_gf ($gfFeature, @lines);
#        $stock->add_gf ($gfFeature, "$line1\n$line2\n$...");
#
# Each list entry gets printed on its own line.
#
sub add_gf {
    my ($self, $feature, @data) = @_;
    unless (@data > 0) {
	carp "Missing parameters to 'add_gf'; nothing will be done";
	return;
    }

    foreach my $line (map {$_ eq '' ? '' : split ("\n", $_, -1)} @data) {
	push @{$self->gf_($feature)}, $line;
	push @{$self->gfOrder}, $feature;