branch-1-6/Bio/Map/PositionI.pm

   1 # $Id$
   2 #
   3 # BioPerl module for Bio::Map::PositionI
   4 #
   5 # Please direct questions and support issues to <bioperl-l@bioperl.org>
   6 #
   7 # Cared for by Sendu Bala <bix@sendu.me.uk>
   8 #
   9 # Copyright Jason Stajich
  10 #
  11 # You may distribute this module under the same terms as perl itself
  12
  13 # POD documentation - main docs before the code
  14
  15 =head1 NAME
  16
  17 Bio::Map::PositionI - Abstracts the notion of a position having a value in the context of a marker and a Map
  18
  19 =head1 SYNOPSIS
  20
  21     # do not use this module directly
  22     # See Bio::Map::Position for an example of
  23     # implementation.
  24
  25 =head1 DESCRIPTION
  26
  27 This object stores one of the postions that a mappable object
  28 (e.g. Marker) may have in a map.
  29
  30 Positions can have non-numeric values or other methods to store the locations,
  31 so they have a method numeric() which does the conversion. numeric()
  32 returns the position in a form that can be compared between other positions of
  33 the same type. It is not necessarily a value suitable for sorting positions (it
  34 may be the distance from the previous position); for that purpose the result of
  35 sortable() should be used.
  36
  37 A 'position', in addition to being a single point, can also be an area and so
  38 can be imagined as a range and compared with other positions on the basis of
  39 overlap, intersection etc.
  40
  41 =head1 FEEDBACK
  42
  43 =head2 Mailing Lists
  44
  45 User feedback is an integral part of the evolution of this and other
  46 Bioperl modules. Send your comments and suggestions preferably to
  47 the Bioperl mailing list.  Your participation is much appreciated.
  48
  49   bioperl-l@bioperl.org                  - General discussion
  50   http://bioperl.org/wiki/Mailing_lists  - About the mailing lists
  51
  52 =head2 Support
  53
  54 Please direct usage questions or support issues to the mailing list:
  55
  56 I<bioperl-l@bioperl.org>
  57
  58 rather than to the module maintainer directly. Many experienced and
  59 reponsive experts will be able look at the problem and quickly
  60 address it. Please include a thorough description of the problem
  61 with code and data examples if at all possible.
  62
  63 =head2 Reporting Bugs
  64
  65 Report bugs to the Bioperl bug tracking system to help us keep track
  66 of the bugs and their resolution. Bug reports can be submitted via the
  67 web:
  68
  69   http://bugzilla.open-bio.org/
  70
  71 =head1 AUTHOR - Jason Stajich
  72
  73 Email jason-at-bioperl.org
  74
  75 =head1 CONTRIBUTORS
  76
  77 Lincoln Stein, lstein-at-cshl.org
  78 Heikki Lehvaslaiho, heikki-at-bioperl-dot-org
  79 Sendu Bala, bix@sendu.me.uk
  80
  81 =head1 APPENDIX
  82
  83 The rest of the documentation details each of the object methods.
  84 Internal methods are usually preceded with a _
  85
  86 =cut
  87
  88 # Let the code begin...
  89
  90 package Bio::Map::PositionI;
  91 use strict;
  92 use Bio::Map::PositionHandler;
  93 use Bio::Map::Mappable;
  94 use Scalar::Util qw(looks_like_number);
  95
  96 use base qw(Bio::Map::EntityI Bio::RangeI);
  97
  98 =head2 EntityI methods
  99
 100  These are fundamental to coordination of Positions and other entities, so are
 101  implemented at the interface level
 102
 103 =cut
 104
 105 =head2 get_position_handler
 106
 107  Title   : get_position_handler
 108  Usage   : my $position_handler = $entity->get_position_handler();
 109  Function: Gets a PositionHandlerI that $entity is registered with.
 110  Returns : Bio::Map::PositionHandlerI object
 111  Args    : none
 112
 113 =cut
 114
 115 sub get_position_handler {
 116     my $self = shift;
 117     unless (defined $self->{_eh}) {
 118         my $ph = Bio::Map::PositionHandler->new(-self => $self);
 119         $self->{_eh} = $ph;
 120         $ph->register;
 121     }
 122     return $self->{_eh};
 123 }
 124
 125 =head2 PositionHandlerI-related methods
 126
 127  These are fundamental to coordination of Positions and other entities, so are
 128  implemented at the interface level
 129
 130 =cut
 131
 132 =head2 map
 133
 134  Title   : map
 135  Usage   : my $map = $position->map();
 136            $position->map($map);
 137  Function: Get/Set the map the position is in.
 138  Returns : L<Bio::Map::MapI>
 139  Args    : none to get
 140            new L<Bio::Map::MapI> to set
 141
 142 =cut
 143
 144 sub map {
 145     my ($self, $map) = @_;
 146     return $self->get_position_handler->map($map);
 147 }
 148
 149 =head2 element
 150
 151  Title   : element
 152  Usage   : my $element = $position->element();
 153            $position->element($element);
 154  Function: Get/Set the element the position is for.
 155  Returns : L<Bio::Map::MappableI>
 156  Args    : none to get
 157            new L<Bio::Map::MappableI> to set
 158
 159 =cut
 160
 161 sub element {
 162     my ($self, $element) = @_;
 163     return $self->get_position_handler->element($element);
 164 }
 165
 166 =head2 marker
 167
 168  Title   : marker
 169  Function: This is a synonym of the element() method
 170  Status  : deprecated, will be removed in the next version
 171
 172 =cut
 173
 174 *marker = \&element;
 175
 176 =head2 PositionI-specific methods
 177
 178 =cut
 179
 180 =head2 value
 181
 182  Title   : value
 183  Usage   : my $pos = $position->value();
 184  Function: Get/Set the value for this position
 185  Returns : scalar, value
 186  Args    : [optional] new value to set
 187
 188 =cut
 189
 190 sub value {
 191     my $self = shift;
 192     $self->throw_not_implemented();
 193 }
 194
 195 =head2 numeric
 196
 197  Title   : numeric
 198  Usage   : my $num = $position->numeric;
 199  Function: Read-only method that is guaranteed to return a numeric
 200            representation of the start of this position.
 201  Returns : scalar numeric
 202  Args    : none to get the co-ordinate normally (see absolute() method), OR
 203            Bio::Map::RelativeI to get the co-ordinate converted to be
 204            relative to what this Relative describes.
 205
 206 =cut
 207
 208 sub numeric {
 209     my $self = shift;
 210     $self->throw_not_implemented();
 211 }
 212
 213 =head2 sortable
 214
 215  Title   : sortable
 216  Usage   : my $num = $position->sortable();
 217  Function: Read-only method that is guaranteed to return a value suitable
 218            for correctly sorting this kind of position amongst other positions
 219            of the same kind on the same map. Note that sorting different kinds
 220            of position together is unlikely to give sane results.
 221  Returns : numeric
 222  Args    : none
 223
 224 =cut
 225
 226 sub sortable {
 227     my $self = shift;
 228     $self->throw_not_implemented();
 229 }
 230
 231 =head2 relative
 232
 233   Title   : relative
 234   Usage   : my $relative = $position->relative();
 235             $position->relative($relative);
 236   Function: Get/set the thing this Position's coordinates (numerical(), start(),
 237             end()) are relative to, as described by a Relative object.
 238   Returns : Bio::Map::RelativeI (default is one describing "relative to the
 239             start of the Position's map")
 240   Args    : none to get, OR
 241             Bio::Map::RelativeI to set
 242
 243 =cut
 244
 245 sub relative {
 246     my $self = shift;
 247     $self->throw_not_implemented();
 248 }
 249
 250 =head2 absolute
 251
 252   Title   : absolute
 253   Usage   : my $absolute = $position->absolute();
 254             $position->absolute($absolute);
 255   Function: Get/set how this Position's co-ordinates (numerical(), start(),
 256             end()) are reported. When absolute is off, co-ordinates are
 257             relative to the thing described by relative(). Ie. the value
 258             returned by start() will be the same as the value you set start()
 259             to. When absolute is on, co-ordinates are converted to be relative
 260             to the start of the map.
 261
 262             So if relative() currently points to a Relative object describing
 263             "relative to another position which is 100 bp from the start of
 264             the map", this Position's start() had been set to 50 and absolute()
 265             returns 1, $position->start() will return 150. If absolute() returns
 266             0 in the same situation, $position->start() would return 50.
 267
 268   Returns : boolean (default 0)
 269   Args    : none to get, OR
 270             boolean to set
 271
 272 =cut
 273
 274 sub absolute {
 275     my $self = shift;
 276     $self->throw_not_implemented();
 277 }
 278
 279 =head2 RangeI-based methods
 280
 281 =cut
 282
 283 =head2 start
 284
 285   Title   : start
 286   Usage   : my $start = $position->start();
 287             $position->start($start);
 288   Function: Get/set the start co-ordinate of this position.
 289   Returns : the start of this position
 290   Args    : scalar numeric to set, OR
 291             none to get the co-ordinate normally (see absolute() method), OR
 292             Bio::Map::RelativeI to get the co-ordinate converted to be
 293             relative to what this Relative describes.
 294
 295 =cut
 296
 297 =head2 end
 298
 299   Title   : end
 300   Usage   : my $end = $position->end();
 301             $position->end($end);
 302   Function: Get/set the end co-ordinate of this position.
 303   Returns : the end of this position
 304   Args    : scalar numeric to set, OR
 305             none to get the co-ordinate normally (see absolute() method), OR
 306             Bio::Map::RelativeI to get the co-ordinate converted to be
 307             relative to what this Relative describes.
 308
 309 =cut
 310
 311 =head2 length
 312
 313   Title   : length
 314   Usage   : $length = $position->length();
 315   Function: Get the length of this position.
 316   Returns : the length of this position
 317   Args    : none
 318
 319 =cut
 320
 321 =head2 strand
 322
 323   Title   : strand
 324   Usage   : $strand = $position->strand();
 325   Function: Get the strand of this position; it is always 1 since maps to not
 326             have strands.
 327   Returns : 1
 328   Args    : none
 329
 330 =cut
 331
 332 sub strand {
 333     return 1;
 334 }
 335
 336 =head2 toString
 337
 338   Title   : toString
 339   Usage   : print $position->toString(), "\n";
 340   Function: stringifies this range
 341   Returns : a string representation of the range of this Position
 342   Args    : optional Bio::Map::RelativeI to have the co-ordinates reported
 343             relative to the thing described by that Relative
 344
 345 =cut
 346
 347 sub toString {
 348     my $self = shift;
 349     $self->throw_not_implemented();
 350 }
 351
 352 =head1 RangeI-related methods
 353
 354 These methods work by considering only the values of start() and end(), as
 355 modified by considering every such co-ordinate relative to the start of the map
 356 (ie. absolute(1) is set temporarily during the calculation), or any supplied
 357 Relative. For the boolean methods, when the comparison Position is on the same
 358 map as the calling Position, there is no point supplying a Relative since the
 359 answer will be the same as without. Relative is most useful when comparing
 360 Positions on different maps and you have a Relative that describes some special
 361 place on each map like 'the start of the gene', where the actual start of the
 362 gene relative to the start of the map is different for each map.
 363
 364 The methods do not consider maps during their calculations - things on different
 365 maps can overlap/contain/intersect/etc. each other.
 366
 367 The geometrical methods (intersect, union etc.) do things to the geometry of
 368 ranges, and return Bio::Map::PositionI compliant objects or triplets (start,
 369 stop, strand) from which new positions could be built. When a PositionI is made
 370 it will have a map transferred to it if all the arguments shared the same map.
 371 If a Relative was supplied the result will have that same Relative.
 372
 373 Note that the strand-testing args are there for compatability with the RangeI
 374 interface. They have no meaning when only using PositionI objects since maps do
 375 not have strands. Typically you will just set the argument to undef if you want
 376 to supply the argument after it.
 377
 378 =head2 equals
 379
 380   Title   : equals
 381   Usage   : if ($p1->equals($p2)) {...}
 382   Function: Test whether $p1 has the same start, end, length as $p2.
 383   Returns : true if they are describing the same position (regardless of map)
 384   Args    : arg #1 = a Bio::RangeI (eg. a Bio::Map::Position) to compare this
 385                      one to (mandatory)
 386             arg #2 = optional strand-testing arg ('strong', 'weak', 'ignore')
 387             arg #3 = optional Bio::Map::RelativeI to ask if the Positions
 388                      equal in terms of their relative position to the thing
 389                      described by that Relative
 390
 391 =cut
 392
 393 sub equals {
 394     # overriding the RangeI implementation so we can handle Relative
 395     my ($self, $other, $so, $rel) = @_;
 396
 397     my ($own_start, $own_end) = $self->_pre_rangei($self, $rel);
 398     my ($other_start, $other_end) = $self->_pre_rangei($other, $rel);
 399
 400     return ($self->_testStrand($other, $so) and
 401             $own_start == $other_start and $own_end == $other_end);
 402 }
 403
 404
 405 =head2 less_than
 406
 407  Title   : less_than
 408  Usage   : if ($position->less_than($other_position)) {...}
 409  Function: Ask if this Position ends before another starts.
 410  Returns : boolean
 411  Args    : arg #1 = a Bio::RangeI (eg. a Bio::Map::Position) to compare this
 412                     one to (mandatory)
 413            arg #2 = optional Bio::Map::RelativeI to ask if the Position is less
 414                     in terms of their relative position to the thing described
 415                     by that Relative
 416
 417 =cut
 418
 419 sub less_than {
 420     my ($self, $other, $rel) = @_;
 421
 422     my ($own_start, $own_end) = $self->_pre_rangei($self, $rel);
 423     my ($other_start, $other_end) = $self->_pre_rangei($other, $rel);
 424
 425     return $own_end < $other_start;
 426 }
 427
 428 =head2 greater_than
 429
 430  Title   : greater_than
 431  Usage   : if ($position->greater_than($other_position)) {...}
 432  Function: Ask if this Position starts after another ends.
 433  Returns : boolean
 434  Args    : arg #1 = a Bio::RangeI (eg. a Bio::Map::Position) to compare this
 435                     one to (mandatory)
 436            arg #2 = optional Bio::Map::RelativeI to ask if the Position is
 437                     greater in terms of their relative position to the thing
 438                     described by that Relative
 439
 440 =cut
 441
 442 sub greater_than {
 443     my ($self, $other, $rel) = @_;
 444
 445     my ($own_start, $own_end) = $self->_pre_rangei($self, $rel);
 446     my ($other_start, $other_end) = $self->_pre_rangei($other, $rel);
 447
 448     return $own_start > $other_end;
 449 }
 450
 451 =head2 overlaps
 452
 453   Title   : overlaps
 454   Usage   : if ($p1->overlaps($p2)) {...}
 455   Function: Tests if $p1 overlaps $p2.
 456   Returns : True if the positions overlap (regardless of map), false otherwise
 457   Args    : arg #1 = a Bio::RangeI (eg. a Bio::Map::Position) to compare this
 458                      one to (mandatory)
 459             arg #2 = optional strand-testing arg ('strong', 'weak', 'ignore')
 460             arg #3 = optional Bio::Map::RelativeI to ask if the Positions
 461                      overlap in terms of their relative position to the thing
 462                      described by that Relative
 463             arg #4 = optional minimum percentage length of the overlap before
 464                      reporting an overlap exists (default 0)
 465
 466 =cut
 467
 468 sub overlaps {
 469     # overriding the RangeI implementation so we can handle Relative
 470     my ($self, $other, $so, $rel, $min_percent) = @_;
 471     $min_percent ||= 0;
 472
 473     my ($own_min, $other_min) = (0, 0);
 474     if ($min_percent > 0) {
 475         $own_min = (($self->length / 100) * $min_percent) - 1;
 476         $other_min = (($other->length / 100) * $min_percent) - 1;
 477     }
 478
 479     my ($own_start, $own_end) = $self->_pre_rangei($self, $rel);
 480     my ($other_start, $other_end) = $self->_pre_rangei($other, $rel);
 481
 482     return ($self->_testStrand($other, $so) and not
 483             (($own_start + $own_min > $other_end or $own_end - $own_min < $other_start) ||
 484              ($own_start > $other_end - $other_min or $own_end < $other_start + $other_min)));
 485 }
 486
 487 =head2 contains
 488
 489   Title   : contains
 490   Usage   : if ($p1->contains($p2)) {...}
 491   Function: Tests whether $p1 totally contains $p2.
 492   Returns : true if the argument is totally contained within this position
 493             (regardless of map), false otherwise
 494   Args    : arg #1 = a Bio::RangeI (eg. a Bio::Map::Position) to compare this
 495                      one to, or scalar number (mandatory)
 496             arg #2 = optional strand-testing arg ('strong', 'weak', 'ignore')
 497             arg #3 = optional Bio::Map::RelativeI to ask if the Position
 498                      is contained in terms of their relative position to the
 499                      thing described by that Relative
 500
 501 =cut
 502
 503 sub contains {
 504     # overriding the RangeI implementation so we can handle Relative
 505     my ($self, $other, $so, $rel) = @_;
 506
 507     my ($own_start, $own_end) = $self->_pre_rangei($self, $rel);
 508     my ($other_start, $other_end) = $self->_pre_rangei($other, $rel);
 509
 510     return ($self->_testStrand($other, $so) and
 511                         $other_start >= $own_start and $other_end <= $own_end);
 512 }
 513
 514 =head2 intersection
 515
 516  Title   : intersection
 517  Usage   : ($start, $stop, $strand) = $p1->intersection($p2)
 518            ($start, $stop, $strand) = Bio::Map::Position->intersection(\@positions);
 519            $mappable = $p1->intersection($p2, undef, $relative);
 520            $mappable = Bio::Map::Position->intersection(\@positions);
 521  Function: gives the range that is contained by all ranges
 522  Returns : undef if they do not overlap, OR
 523            Bio::Map::Mappable object who's positions are the
 524            cross-map-calculated intersection of the input positions on all the
 525            maps that the input positions belong to, OR, in list context, a three
 526            element array (start, end, strand)
 527  Args    : arg #1 = [REQUIRED] a Bio::RangeI (eg. a Bio::Map::Position) to
 528                     compare this one to, or an array ref of Bio::RangeI
 529            arg #2 = optional strand-testing arg ('strong', 'weak', 'ignore')
 530            arg #3 = optional Bio::Map::RelativeI to ask how the Positions
 531                     intersect in terms of their relative position to the thing
 532                     described by that Relative
 533
 534 =cut
 535
 536 sub intersection {
 537     # overriding the RangeI implementation so we can transfer map and handle
 538     # Relative
 539     my ($self, $given, $so, $rel) = @_;
 540         $self->throw("missing arg: you need to pass in another argument") unless $given;
 541
 542     my @positions;
 543     if ($self eq "Bio::Map::PositionI") {
 544                 $self = "Bio::Map::Position";
 545                 $self->warn("calling static methods of an interface is deprecated; use $self instead");
 546         }
 547         if (ref $self) {
 548                 push(@positions, $self);
 549         }
 550     ref($given) eq 'ARRAY' ? push(@positions, @{$given}) : push(@positions, $given);
 551     $self->throw("Need at least 2 Positions") unless @positions >= 2;
 552
 553     my ($intersect, $i_start, $i_end, $c_start, $c_end, %known_maps);
 554     while (@positions > 0) {
 555         unless ($intersect) {
 556             $intersect = shift(@positions);
 557             ($i_start, $i_end) = $self->_pre_rangei($intersect, $rel);
 558             my $map = $intersect->map;
 559             $known_maps{$map->unique_id} = $map;
 560         }
 561
 562         my $compare = shift(@positions);
 563         ($c_start, $c_end) = $self->_pre_rangei($compare, $rel);
 564         return unless $compare->_testStrand($intersect, $so);
 565         if ($compare->isa('Bio::Map::PositionI')) {
 566             my $this_map = $compare->map;
 567             if ($this_map) {
 568                 $known_maps{$this_map->unique_id} = $this_map;
 569             }
 570         }
 571         else {
 572             $self->throw("Only Bio::Map::PositionI objects are supported, not [$compare]");
 573         }
 574
 575         my @starts = sort {$a <=> $b} ($i_start, $c_start);
 576         my @ends   = sort {$a <=> $b} ($i_end, $c_end);
 577
 578         my $start = pop @starts; # larger of the 2 starts
 579         my $end = shift @ends;   # smaller of the 2 ends
 580
 581         my $intersect_strand;    # strand for the intersection
 582         if (defined($intersect->strand) && defined($compare->strand) && $intersect->strand == $compare->strand) {
 583             $intersect_strand = $compare->strand;
 584         }
 585         else {
 586             $intersect_strand = 0;
 587         }
 588
 589         if ($start > $end) {
 590             return;
 591         }
 592         else {
 593             $intersect = $self->new(-start  => $start,
 594                                     -end    => $end,
 595                                     -strand => $intersect_strand);
 596         }
 597     }
 598
 599     $intersect || return;
 600     my ($start, $end, $strand) = ($intersect->start, $intersect->end, $intersect->strand);
 601
 602     my @intersects;
 603     foreach my $known_map (values %known_maps) {
 604         my $new_intersect = $intersect->new(-start => $start,
 605                                             -end => $end,
 606                                             -strand => $strand,
 607                                             -map => $known_map);
 608         $new_intersect->relative($rel) if $rel;
 609         push(@intersects, $new_intersect);
 610     }
 611     unless (@intersects) {
 612         $intersect->relative($rel) if $rel;
 613         @intersects = ($intersect);
 614     }
 615
 616     my $result = Bio::Map::Mappable->new();
 617     $result->add_position(@intersects); # sneaky, add_position can take a list of positions
 618     return $result;
 619 }
 620
 621 =head2 union
 622
 623  Title   : union
 624  Usage   : ($start, $stop, $strand) = $p1->union($p2);
 625            ($start, $stop, $strand) = Bio::Map::Position->union(@positions);
 626            my $mappable = $p1->union($p2);
 627            my $mappable = Bio::Map::Position->union(@positions);
 628  Function: finds the minimal position/range that contains all of the positions
 629  Returns : Bio::Map::Mappable object who's positions are the
 630            cross-map-calculated union of the input positions on all the maps
 631            that the input positions belong to, OR, in list context, a three
 632            element array (start, end, strand)
 633  Args    : a Bio::Map::PositionI to compare this one to, or a list of such
 634            OR
 635            a single Bio::Map::PositionI or array ref of such AND a
 636            Bio::Map::RelativeI to ask for the Position's union in terms of their
 637            relative position to the thing described by that Relative
 638
 639 =cut
 640
 641 sub union {
 642     # overriding the RangeI implementation so we can transfer map and handle
 643     # Relative
 644     my ($self, @args) = @_;
 645     $self->throw("Not enough arguments") unless @args >= 1;
 646
 647     my @positions;
 648     my $rel;
 649     if ($self eq "Bio::Map::PositionI") {
 650                 $self = "Bio::Map::Position";
 651                 $self->warn("calling static methods of an interface is deprecated; use $self instead");
 652         }
 653         if (ref $self) {
 654                 push(@positions, $self);
 655         }
 656     if (ref $args[0] eq 'ARRAY') {
 657         push(@positions, @{shift(@args)});
 658     }
 659     else {
 660         push(@positions, shift(@args));
 661     }
 662     if ($args[0] && $args[0]->isa('Bio::Map::RelativeI')) {
 663         $rel = shift(@args);
 664     }
 665     foreach my $arg (@args) {
 666         # avoid pushing undefined values into @positions
 667         push(@positions, $arg) if $arg;
 668     }
 669     $self->throw("Need at least 2 Positions") unless @positions >= 2;
 670
 671     my (@starts, @ends, %known_maps, $union_strand);
 672     foreach my $compare (@positions) {
 673         # RangeI union allows start or end to be undefined; however _pre_rangei
 674         # will throw
 675         my ($start, $end) = $self->_pre_rangei($compare, $rel);
 676
 677         if ($compare->isa('Bio::Map::PositionI')) {
 678             my $this_map = $compare->map;
 679             if ($this_map) {
 680                 $known_maps{$this_map->unique_id} = $this_map;
 681             }
 682         }
 683         else {
 684             $self->throw("Only Bio::Map::PositionI objects are supported, not [$compare]");
 685         }
 686
 687         if (! defined $union_strand) {
 688                         $union_strand = $compare->strand;
 689                 }
 690         else {
 691                         if (! defined $compare->strand or $union_strand ne $compare->strand) {
 692                                 $union_strand = 0;
 693                         }
 694                 }
 695
 696         push(@starts, $start);
 697         push(@ends, $end);
 698     }
 699
 700         @starts = sort { $a <=> $b } @starts;
 701         @ends   = sort { $a <=> $b } @ends;
 702         my $start = shift @starts;
 703         my $end = pop @ends;
 704
 705     my @unions;
 706     foreach my $known_map (values %known_maps) {
 707         my $new_union = $self->new(-start => $start,
 708                                    -end => $end,
 709                                    -strand => $union_strand,
 710                                    -map => $known_map);
 711         $new_union->relative($rel) if $rel;
 712         push(@unions, $new_union);
 713     }
 714     unless (@unions) {
 715         @unions = ($self->new(-start => $start,
 716                          -end => $end,
 717                          -strand => $union_strand));
 718         $unions[0]->relative($rel) if $rel;
 719     }
 720
 721     my $result = Bio::Map::Mappable->new();
 722     $result->add_position(@unions); # sneaky, add_position can take a list of positions
 723     return $result;
 724 }
 725
 726 =head2 overlap_extent
 727
 728  Title   : overlap_extent
 729  Usage   : ($a_unique,$common,$b_unique) = $a->overlap_extent($b)
 730  Function: Provides actual amount of overlap between two different
 731            positions
 732  Example :
 733  Returns : array of values containing the length unique to the calling
 734            position, the length common to both, and the length unique to
 735            the argument position
 736  Args    : a position
 737
 738 =cut
 739
 740 #*** should this be overridden from RangeI?
 741
 742 =head2 disconnected_ranges
 743
 744  Title   : disconnected_ranges
 745  Usage   : my @disc_ranges = Bio::Map::Position->disconnected_ranges(@ranges);
 746  Function: Creates the minimal set of positions such that each input position is
 747            fully contained by at least one output position, and none of the
 748            output positions overlap.
 749  Returns : Bio::Map::Mappable with the calculated disconnected ranges
 750  Args    : a Bio::Map::PositionI to compare this one to, or a list of such,
 751            OR
 752            a single Bio::Map::PositionI or array ref of such AND a
 753            Bio::Map::RelativeI to consider all Position's co-ordinates in terms
 754            of their relative position to the thing described by that Relative,
 755            AND, optionally, an int for the minimum percentage of overlap that
 756            must be present before considering two ranges to be overlapping
 757            (default 0)
 758
 759 =cut
 760
 761 sub disconnected_ranges {
 762     # overriding the RangeI implementation so we can transfer map and handle
 763     # Relative
 764     my ($self, @args) = @_;
 765     $self->throw("Not enough arguments") unless @args >= 1;
 766
 767     my @positions;
 768     my $rel;
 769     my $overlap = 0;
 770     if ($self eq "Bio::Map::PositionI") {
 771                 $self = "Bio::Map::Position";
 772                 $self->warn("calling static methods of an interface is deprecated; use $self instead");
 773         }
 774         if (ref $self) {
 775                 push(@positions, $self);
 776         }
 777     if (ref $args[0] eq 'ARRAY') {
 778         push(@positions, @{shift(@args)});
 779     }
 780     else {
 781         push(@positions, shift(@args));
 782     }
 783     if ($args[0] && $args[0]->isa('Bio::Map::RelativeI')) {
 784         $rel = shift(@args);
 785         $overlap = shift(@args);
 786     }
 787     foreach my $arg (@args) {
 788         push(@positions, $arg) if $arg;
 789     }
 790     $self->throw("Need at least 2 Positions") unless @positions >= 2;
 791
 792     my %known_maps;
 793     foreach my $pos (@positions) {
 794         $pos->isa('Bio::Map::PositionI') || $self->throw("Must supply only Bio::Map::PositionI objects, not [$pos]");
 795         my $map = $pos->map || next;
 796         $known_maps{$map->unique_id} = $map;
 797     }
 798     my %prior_positions;
 799     foreach my $map (values %known_maps) {
 800         foreach my $pos ($map->get_positions) {
 801             $prior_positions{$pos} = 1;
 802         }
 803     }
 804
 805     my @outranges = ();
 806     foreach my $inrange (@positions) {
 807         my @outranges_new = ();
 808         my %overlapping_ranges = ();
 809
 810         for (my $i=0; $i<@outranges; $i++) {
 811             my $outrange = $outranges[$i];
 812             if ($inrange->overlaps($outrange, undef, $rel, $overlap)) {
 813                 my $union_able = $inrange->union($outrange, $rel); # using $inrange->union($outrange, $rel); gives >6x speedup,
 814                                                                    # but different answer, not necessarily incorrect...
 815                 foreach my $pos ($union_able->get_positions) {
 816                     $overlapping_ranges{$pos->toString} = $pos; # we flatten down to a result on a single map
 817                                                                 # to avoid creating 10s of thousands of positions during this process;
 818                                                                 # we then apply the final answer to all maps at the very end
 819                     last;
 820                 }
 821             }
 822             else {
 823                 push(@outranges_new, $outrange);
 824             }
 825         }
 826
 827         @outranges = @outranges_new;
 828
 829         my @overlappers = values %overlapping_ranges;
 830         if (@overlappers) {
 831             if (@overlappers > 1) {
 832                 my $merged_range_able = shift(@overlappers)->union(\@overlappers, $rel);
 833                 push(@outranges, $merged_range_able->get_positions);
 834             }
 835             else {
 836                 push(@outranges, @overlappers);
 837             }
 838         }
 839         else {
 840             push(@outranges, $self->new(-start => $inrange->start($rel), -end => $inrange->end($rel), -strand => $inrange->strand, -map => $inrange->map, -relative => $rel));
 841         }
 842     }
 843
 844     # purge positions that were created whilst calculating the answer, but
 845     # aren't the final answer and weren't there previously
 846     my %answers = map { $_ => 1 } @outranges;
 847     foreach my $map (values %known_maps) {
 848         foreach my $pos ($map->get_positions) {
 849             if (! exists $prior_positions{$pos} && ! exists $answers{$pos}) {
 850                 $map->purge_positions($pos);
 851             }
 852         }
 853     }
 854
 855     my %post_positions;
 856     foreach my $map (values %known_maps) {
 857         foreach my $pos ($map->get_positions) {
 858             $post_positions{$pos} = 1;
 859         }
 860     }
 861
 862     @outranges || return;
 863
 864     # make an outrange on all known maps
 865     my @final_positions;
 866     foreach my $map (values %known_maps) {
 867         foreach my $pos (@outranges) {
 868             if ($pos->map eq $map) {
 869                 push(@final_positions, $pos);
 870             }
 871             else {
 872                 push(@final_positions, $pos->new(-start => $pos->start,
 873                                                  -end => $pos->end,
 874                                                  -relative => $pos->relative,
 875                                                  -map => $map));
 876             }
 877         }
 878     }
 879
 880     # assign the positions to a result mappable
 881     my $result = Bio::Map::Mappable->new();
 882     $result->add_position(@final_positions); # sneaky, add_position can take a list of positions
 883     return $result;
 884 }
 885
 886 # get start & end suitable for rangeI methods, taking relative into account
 887 sub _pre_rangei {
 888     my ($self, $other, $rel) = @_;
 889     $self->throw("Must supply an object") unless $other;
 890     if ($rel) {
 891         $self->throw("Must supply an object for the Relative argument") unless ref($rel);
 892         $self->throw("This is [$rel], not a Bio::Map::RelativeI") unless $rel->isa('Bio::Map::RelativeI');
 893     }
 894
 895     my ($other_start, $other_end);
 896     if (ref($other)) {
 897         if (ref($other) eq 'ARRAY') {
 898             $self->throw("_pre_rangei got an array");
 899         }
 900         $self->throw("This is [$other], not a Bio::RangeI object") unless defined $other && $other->isa('Bio::RangeI');
 901
 902         if ($other->isa('Bio::Map::PositionI')) {
 903             # to get the desired start/end we need the position to be on a map;
 904             # if it isn't on one temporarily place it on self's map
 905             # - this lets us have 'generic' positions that aren't on any map
 906             # but have a relative defined and can thus be usefully compared to
 907             # positions that /are/ on maps
 908             my $other_map = $other->map;
 909             unless ($other_map) {
 910                 my $self_map = $self->map || $self->throw("Trying to compare two positions but neither had been placed on a map");
 911                 $other->map($self_map);
 912             }
 913
 914             # want start and end positions relative to the supplied rel or map start
 915             $rel ||= $other->absolute_relative;
 916             $other_start = $other->start($rel);
 917             $other_end = $other->end($rel);
 918
 919             unless ($other_map) {
 920                 $self->map->purge_positions($other);
 921             }
 922         }
 923         else {
 924             $other_start = $other->start;
 925             $other_end = $other->end;
 926         }
 927     }
 928     else {
 929         $self->throw("not a number") unless looks_like_number($other);
 930         $other_start = $other_end = $other;
 931     }
 932
 933         $other->throw("start is undefined") unless defined $other_start;
 934         $other->throw("end is undefined") unless defined $other_end;
 935
 936     return ($other_start, $other_end);
 937 }
 938
 939 1;