[med-svn] [gmap] 01/02: Imported Upstream version 2014-12-25
Alex Mestiashvili
malex-guest at moszumanska.debian.org
Fri Mar 20 11:07:24 UTC 2015
This is an automated email from the git hooks/post-receive script.
malex-guest pushed a commit to branch master
in repository gmap.
commit 877e0f46dd264974d19fcefa6a5bf4e31e7aeab1
Author: Alexandre Mestiashvili <alex at biotec.tu-dresden.de>
Date: Fri Mar 20 11:55:02 2015 +0100
Imported Upstream version 2014-12-25
---
ChangeLog | 11 ++
VERSION | 2 +-
configure | 24 +--
src/samprint.c | 86 ++++++-----
src/stage3hr.c | 472 +++++++++++++++++++++++++++++++++++++-------------------
src/stage3hr.h | 7 +-
src/substring.c | 10 +-
src/substring.h | 6 +-
8 files changed, 395 insertions(+), 223 deletions(-)
diff --git a/ChangeLog b/ChangeLog
index a1b9de3..822b46c 100644
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,14 @@
+2015-03-18 twu
+
+ * public-2014-12-17, src, stage3hr.c: Merged revision 161197 from trunk to
+ fix a fatal bug in adjust_hardclips
+
+ * VERSION, public-2014-12-17, src: Updated version number
+
+ * samprint.c, stage3hr.c, stage3hr.h, substring.c, substring.h: Merged
+ revisions 160877 through 161118 from trunk to undo use of
+ substring_hardclipped, and to use substring_LtoH instead
+
2015-03-13 twu
* outbuffer.c, samprint.c, samprint.h: Applied revision 160876 from trunk to
diff --git a/VERSION b/VERSION
index 83a8c64..4568d14 100644
--- a/VERSION
+++ b/VERSION
@@ -1 +1 @@
-2014-12-24
\ No newline at end of file
+2014-12-25
\ No newline at end of file
diff --git a/configure b/configure
index 4a3fa72..1631267 100755
--- a/configure
+++ b/configure
@@ -1,6 +1,6 @@
#! /bin/sh
# Guess values for system-dependent variables and create Makefiles.
-# Generated by GNU Autoconf 2.63 for gmap 2014-12-24.
+# Generated by GNU Autoconf 2.63 for gmap 2014-12-25.
#
# Report bugs to <Thomas Wu <twu at gene.com>>.
#
@@ -745,8 +745,8 @@ SHELL=${CONFIG_SHELL-/bin/sh}
# Identity of this package.
PACKAGE_NAME='gmap'
PACKAGE_TARNAME='gmap'
-PACKAGE_VERSION='2014-12-24'
-PACKAGE_STRING='gmap 2014-12-24'
+PACKAGE_VERSION='2014-12-25'
+PACKAGE_STRING='gmap 2014-12-25'
PACKAGE_BUGREPORT='Thomas Wu <twu at gene.com>'
ac_unique_file="src/gmap.c"
@@ -1512,7 +1512,7 @@ if test "$ac_init_help" = "long"; then
# Omit some internal or obsolete options to make the list less imposing.
# This message is too long to be a string in the A/UX 3.1 sh.
cat <<_ACEOF
-\`configure' configures gmap 2014-12-24 to adapt to many kinds of systems.
+\`configure' configures gmap 2014-12-25 to adapt to many kinds of systems.
Usage: $0 [OPTION]... [VAR=VALUE]...
@@ -1583,7 +1583,7 @@ fi
if test -n "$ac_init_help"; then
case $ac_init_help in
- short | recursive ) echo "Configuration of gmap 2014-12-24:";;
+ short | recursive ) echo "Configuration of gmap 2014-12-25:";;
esac
cat <<\_ACEOF
@@ -1718,7 +1718,7 @@ fi
test -n "$ac_init_help" && exit $ac_status
if $ac_init_version; then
cat <<\_ACEOF
-gmap configure 2014-12-24
+gmap configure 2014-12-25
generated by GNU Autoconf 2.63
Copyright (C) 1992, 1993, 1994, 1995, 1996, 1998, 1999, 2000, 2001,
@@ -1732,7 +1732,7 @@ cat >config.log <<_ACEOF
This file contains any messages produced by compilers while
running configure, to aid debugging if configure makes a mistake.
-It was created by gmap $as_me 2014-12-24, which was
+It was created by gmap $as_me 2014-12-25, which was
generated by GNU Autoconf 2.63. Invocation command line was
$ $0 $@
@@ -2102,8 +2102,8 @@ ac_compiler_gnu=$ac_cv_c_compiler_gnu
{ $as_echo "$as_me:$LINENO: checking package version" >&5
$as_echo_n "checking package version... " >&6; }
-{ $as_echo "$as_me:$LINENO: result: 2014-12-24" >&5
-$as_echo "2014-12-24" >&6; }
+{ $as_echo "$as_me:$LINENO: result: 2014-12-25" >&5
+$as_echo "2014-12-25" >&6; }
### Read defaults
@@ -4162,7 +4162,7 @@ fi
# Define the identity of the package.
PACKAGE=gmap
- VERSION=2014-12-24
+ VERSION=2014-12-25
cat >>confdefs.h <<_ACEOF
@@ -26480,7 +26480,7 @@ exec 6>&1
# report actual input values of CONFIG_FILES etc. instead of their
# values after options handling.
ac_log="
-This file was extended by gmap $as_me 2014-12-24, which was
+This file was extended by gmap $as_me 2014-12-25, which was
generated by GNU Autoconf 2.63. Invocation command line was
CONFIG_FILES = $CONFIG_FILES
@@ -26543,7 +26543,7 @@ Report bugs to <bug-autoconf at gnu.org>."
_ACEOF
cat >>$CONFIG_STATUS <<_ACEOF || ac_write_fail=1
ac_cs_version="\\
-gmap config.status 2014-12-24
+gmap config.status 2014-12-25
configured by $0, generated by GNU Autoconf 2.63,
with options \\"`$as_echo "$ac_configure_args" | sed 's/^ //; s/[\\""\`\$]/\\\\&/g'`\\"
diff --git a/src/samprint.c b/src/samprint.c
index c8cac55..eda8d82 100644
--- a/src/samprint.c
+++ b/src/samprint.c
@@ -1,4 +1,4 @@
-static char rcsid[] = "$Id: samprint.c 160877 2015-03-13 00:31:23Z twu $";
+static char rcsid[] = "$Id: samprint.c 161183 2015-03-18 17:04:33Z twu $";
#ifdef HAVE_CONFIG_H
#include <config.h>
#endif
@@ -699,9 +699,10 @@ adjust_hardclips (int *hardclip_low, Stage3end_T hit_low, int low_querylength,
Chrpos_T
SAM_compute_chrpos (int hardclip_low, int hardclip_high, Stage3end_T this, int querylength) {
- Substring_T substring_hardclipped;
+ Substring_T substring_low;
Chrpos_T chrpos;
int querystart, queryend;
+ int trim_low;
if (this == NULL) {
return 0U;
@@ -710,46 +711,57 @@ SAM_compute_chrpos (int hardclip_low, int hardclip_high, Stage3end_T this, int q
chrpos = Pair_genomicpos_low(hardclip_low,hardclip_high,Stage3end_pairarray(this),Stage3end_npairs(this),
querylength,/*watsonp*/Stage3end_plusp(this),hide_soft_clips_p);
- } else if (Stage3end_plusp(this) == true) {
- substring_hardclipped = Stage3end_substring_containing(this,hardclip_low);
- debug4(printf("Plus: Substring containing hardclip_low %d is %d..%d => %u..%u\n",
- hardclip_low,Substring_querystart(substring_hardclipped),Substring_queryend(substring_hardclipped),
- Substring_chrstart(substring_hardclipped),Substring_chrend(substring_hardclipped)));
+ } else if (hide_soft_clips_p == true) {
+ substring_low = Stage3end_substring_low(this,hardclip_low);
+ if (Stage3end_plusp(this) == true) {
+ /* Add 1 to report in 1-based coordinates */
+ chrpos = Substring_alignstart_chr(substring_low) + 1U;
+ querystart = Substring_querystart_orig(substring_low);
+ chrpos += hardclip_low - querystart;
- /* Add 1 to report in 1-based coordinates */
- if (hide_soft_clips_p == true) {
- chrpos = Substring_alignstart(substring_hardclipped) - Substring_chroffset(substring_hardclipped) + 1U;
- querystart = Substring_querystart_orig(substring_hardclipped);
- /* queryend = Substring_queryend_orig(Stage3end_substring_high(this)); */
} else {
- chrpos = Substring_alignstart_trim(substring_hardclipped) - Substring_chroffset(substring_hardclipped) + 1U;
- querystart = Substring_querystart(substring_hardclipped);
- /* queryend = Substring_queryend(Stage3end_substring_high(this)); */
+ /* Add 1 to report in 1-based coordinates */
+ chrpos = Substring_alignend_chr(substring_low) + 1U;
+ queryend = Substring_queryend_orig(substring_low);
+ chrpos += hardclip_low - (querylength - queryend);
}
- debug4(printf("Incrementing chrpos by %d\n",hardclip_low - querystart));
- chrpos += hardclip_low - querystart;
-
} else {
- substring_hardclipped = Stage3end_substring_containing(this,querylength - 1 - hardclip_low);
- debug4(printf("Minus: Substring containing querylength - 1 - hardclip_low %d is %d..%d => %u..%u\n",
- querylength - 1 - hardclip_low,
- Substring_querystart(substring_hardclipped),Substring_queryend(substring_hardclipped),
- Substring_chrstart(substring_hardclipped),Substring_chrend(substring_hardclipped)));
+ if (Stage3end_plusp(this) == true) {
+ if ((trim_low = Stage3end_trim_left_raw(this)) < hardclip_low) {
+ trim_low = hardclip_low;
+ }
+ substring_low = Stage3end_substring_low(this,trim_low);
+ debug4(printf("Plus: Substring containing trim_low %d is %d..%d => %u..%u\n",
+ trim_low,Substring_querystart(substring_low),Substring_queryend(substring_low),
+ Substring_alignstart_chr(substring_low),Substring_alignend_chr(substring_low)));
+
+ /* Add 1 to report in 1-based coordinates */
+ chrpos = Substring_alignstart_chr(substring_low) + 1U;
+ querystart = Substring_querystart_orig(substring_low);
+
+ debug4(printf("Incrementing chrpos by %d = %d - %d\n",trim_low - querystart,trim_low,querystart));
+ chrpos += trim_low - querystart;
- /* Add 1 to report in 1-based coordinates */
- if (hide_soft_clips_p == true) {
- chrpos = Substring_alignend(substring_hardclipped) - Substring_chroffset(substring_hardclipped) + 1U;
- /* querystart = Substring_querystart_orig(Stage3end_substring_high(this)); */
- queryend = Substring_queryend_orig(substring_hardclipped);
} else {
- chrpos = Substring_alignend_trim(substring_hardclipped) - Substring_chroffset(substring_hardclipped) + 1U;
- /* querystart = Substring_querystart(Stage3end_substring_high(this)); */
- queryend = Substring_queryend(substring_hardclipped);
- }
+ if ((trim_low = Stage3end_trim_right_raw(this)) < hardclip_low) {
+ trim_low = hardclip_low;
+ }
+
+ substring_low = Stage3end_substring_low(this,trim_low);
+ debug4(printf("Minus: Substring containing trim_low %d is %d..%d => %u..%u\n",
+ trim_low,Substring_querystart(substring_low),Substring_queryend(substring_low),
+ Substring_alignstart_chr(substring_low),Substring_alignend_chr(substring_low)));
+
+ /* Add 1 to report in 1-based coordinates */
+ chrpos = Substring_alignend_chr(substring_low) + 1U;
+ queryend = Substring_queryend_orig(substring_low);
+ debug4(printf("queryend is %d, chrpos is %u\n",queryend,chrpos));
- debug4(printf("Incrementing chrpos by %d\n",queryend - (querylength - hardclip_low)));
- chrpos += queryend - (querylength - hardclip_low);
+ debug4(printf("Incrementing chrpos by %d = %d - (querylength %d - %d)\n",
+ trim_low - (querylength - queryend),trim_low,querylength,queryend));
+ chrpos += trim_low - (querylength - queryend);
+ }
}
return chrpos;
@@ -5521,13 +5533,15 @@ print_exon_exon (FILE *fp, char *abbrev, Stage3end_T this, Stage3end_T mate,
querylength = Shortread_fulllength(queryseq);
donor_chrpos = SAM_compute_chrpos(hardclip_low,hardclip_high,this,querylength);
acceptor_chrpos = SAM_compute_chrpos(hardclip_low,hardclip_high,this,querylength);
- if (Stage3end_substring_low(this) == donor) {
+ if (Stage3end_substring_low(this,hardclip_low) == donor) {
concordant_chrpos = donor_chrpos;
- } else if (Stage3end_substring_low(this) == acceptor) {
+ } else if (Stage3end_substring_low(this,hardclip_low) == acceptor) {
concordant_chrpos = acceptor_chrpos;
} else {
+#if 0
fprintf(stderr,"Stage3end_substring_low %p is neither donor %p or acceptor %p\n",
Stage3end_substring_low(this),donor,acceptor);
+#endif
concordant_chrpos = 0U;
}
diff --git a/src/stage3hr.c b/src/stage3hr.c
index 8002396..8267333 100644
--- a/src/stage3hr.c
+++ b/src/stage3hr.c
@@ -1,4 +1,4 @@
-static char rcsid[] = "$Id: stage3hr.c 160873 2015-03-13 00:22:08Z twu $";
+static char rcsid[] = "$Id: stage3hr.c 161199 2015-03-18 18:26:53Z twu $";
#ifdef HAVE_CONFIG_H
#include <config.h>
#endif
@@ -568,8 +568,7 @@ struct T {
int npairs;
int nsegments;
- Substring_T substring_low; /* For SAM output */
- Substring_T substring_high; /* For SAM output */
+ List_T substring_LtoH; /* Chromosomal low-to-high, for computing chrpos */
bool paired_usedp;
bool paired_seenp; /* for paired-end. set to true by Stage3_pair_up(). */
@@ -771,10 +770,13 @@ Stage3end_genomic_alignment_length (T this) {
Chrpos_T
Stage3end_chrpos_low_trim (T this) {
+ Substring_T substring_low;
+
+ substring_low = (Substring_T) List_head(this->substring_LtoH);
if (this->plusp == true) {
- return Substring_alignstart_trim(this->substring_low) - Substring_chroffset(this->substring_low);
+ return Substring_alignstart_trim(substring_low) - Substring_chroffset(substring_low);
} else {
- return Substring_alignend_trim(this->substring_low) - Substring_chroffset(this->substring_low);
+ return Substring_alignend_trim(substring_low) - Substring_chroffset(substring_low);
}
}
@@ -987,20 +989,81 @@ Stage3end_substringA (T this) {
}
Substring_T
-Stage3end_substring_low (T this) {
+Stage3end_substring_low (T this, int trim_low) {
+ List_T p;
+ Substring_T substring;
+
if (this == NULL) {
return (Substring_T) NULL;
+
+ } else if (this->plusp == true) {
+ p = this->substring_LtoH;
+#ifdef DEBUG13
+ if (p != NULL) {
+ printf("Substring is %d..%d against trim_low %d\n",
+ Substring_querystart((Substring_T) List_head(p)),Substring_queryend((Substring_T) List_head(p)),
+ trim_low);
+ }
+#endif
+ while (p != NULL && Substring_queryend((Substring_T) List_head(p)) < trim_low) {
+ p = List_next(p);
+#ifdef DEBUG13
+ if (p != NULL) {
+ printf("Substring is %d..%d against trim_low %d\n",
+ Substring_querystart((Substring_T) List_head(p)),Substring_queryend((Substring_T) List_head(p)),
+ trim_low);
+ }
+#endif
+ }
+ assert(p != NULL);
+ if (p == NULL) {
+ return (Substring_T) NULL;
+ } else {
+ return (Substring_T) List_head(p);
+ }
+
} else {
- return this->substring_low;
+#ifdef DEBUG13
+ for (p = this->substring_LtoH; p != NULL; p = List_next(p)) {
+ printf("LtoH: %d..%d\n",
+ Substring_querystart((Substring_T) List_head(p)),Substring_queryend((Substring_T) List_head(p)));
+ }
+#endif
+
+ p = this->substring_LtoH;
+#ifdef DEBUG13
+ if (p != NULL) {
+ printf("Substring is %d..%d against %d = querylength %d - trim_low %d\n",
+ Substring_querystart((Substring_T) List_head(p)),Substring_queryend((Substring_T) List_head(p)),
+ this->querylength_adj - trim_low,this->querylength_adj,trim_low);
+ }
+#endif
+ while (p != NULL && Substring_querystart((Substring_T) List_head(p)) >= this->querylength_adj - trim_low) {
+ p = List_next(p);
+#ifdef DEBUG13
+ if (p != NULL) {
+ printf("Substring is %d..%d against %d = querylength %d - trim_low %d\n",
+ Substring_querystart((Substring_T) List_head(p)),Substring_queryend((Substring_T) List_head(p)),
+ this->querylength_adj - trim_low,this->querylength_adj,trim_low);
+ }
+#endif
+ }
+ assert(p != NULL);
+ if (p == NULL) {
+ return (Substring_T) NULL;
+ } else {
+ return (Substring_T) List_head(p);
+ }
}
}
+
Substring_T
Stage3end_substring_high (T this) {
if (this == NULL) {
return (Substring_T) NULL;
} else {
- return this->substring_high;
+ return (Substring_T) List_last_value(this->substring_LtoH);
}
}
@@ -1643,6 +1706,8 @@ Stage3end_free (T *old) {
Substring_free(&(*old)->substring0);
}
+ List_free(&(*old)->substring_LtoH);
+
FREE_OUT(*old);
return;
@@ -1850,6 +1915,7 @@ find_ilengths (int *ilength_low, int *ilength_high, Stage3end_T hit, Univcoord_T
*ilength_low = hit->pairarray[hit->npairs - 1].querypos - hit->pairarray[i].querypos + 1;
*ilength_high = hit->pairarray[i].querypos - hit->pairarray[0].querypos + 1;
}
+ debug13(printf("GMAP: Have ilength_low %d and ilength_high %d\n",*ilength_low,*ilength_high));
} else if (hit->plusp == true) {
debug13(printf("plus. Checking common genomicpos %llu against substring0 %p, substring1 %p, substring2 %p\n",
@@ -1885,6 +1951,7 @@ find_ilengths (int *ilength_low, int *ilength_high, Stage3end_T hit, Univcoord_T
} else {
abort();
}
+ debug13(printf("Plus: Have ilength_low %d and ilength_high %d\n",*ilength_low,*ilength_high));
} else {
debug13(printf("minus. Checking common genomicpos %llu against substring0 %p, substring1 %p, substring2 %p\n",
@@ -1919,9 +1986,9 @@ find_ilengths (int *ilength_low, int *ilength_high, Stage3end_T hit, Univcoord_T
} else {
abort();
}
+ debug13(printf("Minus: Have ilength_low %d and ilength_high %d\n",*ilength_low,*ilength_high));
}
- debug13(printf("Have ilength_low %d and ilength_high %d\n",*ilength_low,*ilength_high));
return;
}
@@ -2545,6 +2612,16 @@ adjust_hardclips (int *hardclip_low, Stage3end_T hit_low, int low_querylength,
high_querypos = high_querylength - 1 - (*hardclip_high);
while (low_querypos < low_querylength && high_querypos < high_querylength &&
+ (Pairarray_contains_p(low_pairarray,low_npairs,low_querypos) == false ||
+ Pairarray_contains_p(high_pairarray,high_npairs,high_querypos) == false)) {
+ (*hardclip_low)++;
+ (*hardclip_high)--;
+ low_querypos++;
+ high_querypos++;
+ debug13(printf("Trying hardclip_low %d and hardclip_high %d\n",*hardclip_low,*hardclip_high));
+ }
+
+ while (low_querypos < low_querylength && high_querypos < high_querylength &&
(Pairarray_contains_p(low_pairarray,low_npairs,low_querypos-1) == false ||
Pairarray_contains_p(high_pairarray,high_npairs,high_querypos-1) == false)) {
/* Shift clip querypos upward (rightward on chromosome) */
@@ -2581,6 +2658,16 @@ adjust_hardclips (int *hardclip_low, Stage3end_T hit_low, int low_querylength,
high_querypos = *hardclip_high;
while (low_querypos >= 0 && high_querypos >= 0 &&
+ (Pairarray_contains_p(low_pairarray,low_npairs,low_querypos) == false ||
+ Pairarray_contains_p(high_pairarray,high_npairs,high_querypos) == false)) {
+ (*hardclip_low)++;
+ (*hardclip_high)--;
+ low_querypos--;
+ high_querypos--;
+ debug13(printf("Trying hardclip_low %d and hardclip_high %d\n",*hardclip_low,*hardclip_high));
+ }
+
+ while (low_querypos >= 0 && high_querypos >= 0 &&
(Pairarray_contains_p(low_pairarray,low_npairs,low_querypos+1) == false ||
Pairarray_contains_p(high_pairarray,high_npairs,high_querypos+1) == false)) {
/* Shift clip querypos downward (rightward on chromosome) */
@@ -2622,6 +2709,18 @@ adjust_hardclips (int *hardclip_low, Stage3end_T hit_low, int low_querylength,
low_querypos = *hardclip_low;
high_querypos = high_querylength - 1 - (*hardclip_high);
high_substring = Stage3end_substring_containing(hit_high,high_querypos);
+
+ while (low_querypos < low_querylength && high_querypos < high_querylength &&
+ (Pairarray_contains_p(low_pairarray,low_npairs,low_querypos) == false ||
+ high_substring == NULL)) {
+ (*hardclip_low)++;
+ (*hardclip_high)--;
+ low_querypos++;
+ high_querypos++;
+ high_substring = Stage3end_substring_containing(hit_high,high_querypos);
+ debug13(printf("Trying hardclip_low %d and hardclip_high %d\n",*hardclip_low,*hardclip_high));
+ }
+
while (low_querypos < low_querylength && high_substring != NULL &&
(Pairarray_contains_p(low_pairarray,low_npairs,low_querypos-1) == false ||
Substring_contains_p(high_substring,high_querypos-1) == false)) {
@@ -2660,6 +2759,18 @@ adjust_hardclips (int *hardclip_low, Stage3end_T hit_low, int low_querylength,
high_querypos = *hardclip_high;
high_substring = Stage3end_substring_containing(hit_high,high_querypos);
+ while (low_querypos >= 0 && high_querypos >= 0 &&
+ (Pairarray_contains_p(low_pairarray,low_npairs,low_querypos) == false ||
+ high_substring == NULL)) {
+ /* Shift clip querypos downward (rightward on chromosome) */
+ (*hardclip_low)++;
+ (*hardclip_high)--;
+ low_querypos--;
+ high_querypos--;
+ high_substring = Stage3end_substring_containing(hit_high,high_querypos);
+ debug13(printf("Trying hardclip_low %d and hardclip_high %d\n",*hardclip_low,*hardclip_high));
+ }
+
while (low_querypos >= 0 && high_substring != NULL &&
(Pairarray_contains_p(low_pairarray,low_npairs,low_querypos+1) == false ||
Substring_contains_p(high_substring,high_querypos+1) == false)) {
@@ -2704,6 +2815,17 @@ adjust_hardclips (int *hardclip_low, Stage3end_T hit_low, int low_querylength,
high_querypos = high_querylength - 1 - (*hardclip_high);
low_substring = Stage3end_substring_containing(hit_low,low_querypos);
+ while (low_querypos < low_querylength && high_querypos < high_querylength &&
+ (low_substring == NULL ||
+ Pairarray_contains_p(high_pairarray,high_npairs,high_querypos) == false)) {
+ (*hardclip_low)++;
+ (*hardclip_high)--;
+ low_querypos++;
+ high_querypos++;
+ low_substring = Stage3end_substring_containing(hit_low,low_querypos);
+ debug13(printf("Trying hardclip_low %d and hardclip_high %d\n",*hardclip_low,*hardclip_high));
+ }
+
while (low_substring != NULL && high_querypos < high_querylength &&
(Substring_contains_p(low_substring,low_querypos-1) == false ||
Pairarray_contains_p(high_pairarray,high_npairs,high_querypos-1) == false)) {
@@ -2742,6 +2864,17 @@ adjust_hardclips (int *hardclip_low, Stage3end_T hit_low, int low_querylength,
high_querypos = *hardclip_high;
low_substring = Stage3end_substring_containing(hit_low,low_querypos);
+ while (low_querypos >= 0 && high_querypos >= 0 &&
+ (low_substring == NULL ||
+ Pairarray_contains_p(high_pairarray,high_npairs,high_querypos) == false)) {
+ (*hardclip_low)++;
+ (*hardclip_high)--;
+ low_querypos--;
+ high_querypos--;
+ low_substring = Stage3end_substring_containing(hit_low,low_querypos);
+ debug13(printf("Trying hardclip_low %d and hardclip_high %d\n",*hardclip_low,*hardclip_high));
+ }
+
while (low_substring != NULL && high_querypos >= 0 &&
(Substring_contains_p(low_substring,low_querypos+1) == false ||
Pairarray_contains_p(high_pairarray,high_npairs,high_querypos+1) == false)) {
@@ -2785,6 +2918,19 @@ adjust_hardclips (int *hardclip_low, Stage3end_T hit_low, int low_querylength,
low_substring = Stage3end_substring_containing(hit_low,low_querypos);
high_substring = Stage3end_substring_containing(hit_high,high_querypos);
+ debug13(printf("Shifting upward with low_substring %p and high_substring %p\n",low_substring,high_substring));
+ while (low_querypos < low_querylength && high_querypos < high_querylength &&
+ (low_substring == NULL || high_substring == NULL)) {
+ (*hardclip_low)++;
+ (*hardclip_high)--;
+ low_querypos++;
+ high_querypos++;
+ low_substring = Stage3end_substring_containing(hit_low,low_querypos);
+ high_substring = Stage3end_substring_containing(hit_high,high_querypos);
+ debug13(printf("Trying hardclip_low %d and hardclip_high %d\n",*hardclip_low,*hardclip_high));
+ }
+
+ debug13(printf("Shifting upward with low_substring %p and high_substring %p\n",low_substring,high_substring));
while (low_substring != NULL && high_substring != NULL &&
(Substring_contains_p(low_substring,low_querypos-1) == false ||
Substring_contains_p(high_substring,high_querypos-1) == false)) {
@@ -2801,6 +2947,7 @@ adjust_hardclips (int *hardclip_low, Stage3end_T hit_low, int low_querylength,
debug13(printf("Trying hardclip_low %d and hardclip_high %d\n",*hardclip_low,*hardclip_high));
}
+ debug13(printf("Shifting downward with low_substring %p and high_substring %p\n",low_substring,high_substring));
while (low_substring != NULL && high_substring != NULL &&
(Substring_contains_p(low_substring,low_querypos+1) == false ||
Substring_contains_p(high_substring,high_querypos+1) == false)) {
@@ -2830,6 +2977,19 @@ adjust_hardclips (int *hardclip_low, Stage3end_T hit_low, int low_querylength,
low_substring = Stage3end_substring_containing(hit_low,low_querypos);
high_substring = Stage3end_substring_containing(hit_high,high_querypos);
+ debug13(printf("Shifting downward with low_substring %p and high_substring %p\n",low_substring,high_substring));
+ while (low_querypos >= 0 && high_querypos >= 0 &&
+ (low_substring == NULL || high_substring == NULL)) {
+ (*hardclip_low)++;
+ (*hardclip_high)--;
+ low_querypos--;
+ high_querypos--;
+ low_substring = Stage3end_substring_containing(hit_low,low_querypos);
+ high_substring = Stage3end_substring_containing(hit_high,high_querypos);
+ debug13(printf("Trying hardclip_low %d and hardclip_high %d\n",*hardclip_low,*hardclip_high));
+ }
+
+ debug13(printf("Shifting downward with low_substring %p and high_substring %p\n",low_substring,high_substring));
while (low_substring != NULL && high_substring != NULL &&
(Substring_contains_p(low_substring,low_querypos+1) == false ||
Substring_contains_p(high_substring,high_querypos+1) == false)) {
@@ -2846,6 +3006,7 @@ adjust_hardclips (int *hardclip_low, Stage3end_T hit_low, int low_querylength,
debug13(printf("Trying hardclip_low %d and hardclip_high %d\n",*hardclip_low,*hardclip_high));
}
+ debug13(printf("Shifting upward with low_substring %p and high_substring %p\n",low_substring,high_substring));
while (low_substring != NULL && high_substring != NULL &&
(Substring_contains_p(low_substring,low_querypos-1) == false ||
Substring_contains_p(high_substring,high_querypos-1) == false)) {
@@ -2883,10 +3044,8 @@ adjust_hardclips (int *hardclip_low, Stage3end_T hit_low, int low_querylength,
int
Stage3pair_overlap (int *hardclip5_low, int *hardclip5_high, int *hardclip3_low, int *hardclip3_high, Stage3pair_T this) {
Stage3end_T hit5, hit3;
- int totallength, insertlength;
int overlap;
int clipdir;
- int hit5_trimmed_length, hit3_trimmed_length;
int ilength53, ilength35, ilength5_low, ilength5_high, ilength3_low, ilength3_high;
int common_shift, common_left, common_right;
Univcoord_T common_genomicpos;
@@ -2912,16 +3071,17 @@ Stage3pair_overlap (int *hardclip5_low, int *hardclip5_high, int *hardclip3_low,
hit3->trim_left,hit3->start_amb_length,hit3->trim_right,hit3->end_amb_length));
if (hit5->plusp == true) {
/* plus */
+#if 0
hit5_trimmed_length = hit5->querylength - hit5->trim_left - hit5->trim_right - hit5->start_amb_length - hit5->end_amb_length;
hit3_trimmed_length = hit3->querylength - hit3->trim_left - hit3->trim_right - hit3->start_amb_length - hit3->end_amb_length;
totallength = hit5_trimmed_length + hit3_trimmed_length;
debug13(printf("totallength = %d, hit5 trimmed length = %d, hit3 trimmed length = %d\n",
totallength,hit5_trimmed_length,hit3_trimmed_length));
-
debug13(printf("original insertlength: %d, trim+amb5: %d..%d, trim+amb3: %d..%d\n",
this->insertlength,hit5->trim_left + hit5->start_amb_length,
hit5->trim_right + hit5->end_amb_length,hit3->trim_left + hit3->start_amb_length,
hit3->trim_right + hit3->end_amb_length));
+#endif
if ((common_genomicpos = pair_common_genomicpos(hit5,hit3)) == 0) {
debug13(printf("Cannot determine a common point, so returning 0\n"));
@@ -2947,13 +3107,10 @@ Stage3pair_overlap (int *hardclip5_low, int *hardclip5_high, int *hardclip3_low,
if ((ilength53 = ilength5_low + ilength3_high - 1) >= (ilength35 = ilength3_low + ilength5_high - 1)) {
/* Use >=, not >, so we favor clipping heads over clipping tails in case of a tie */
debug13(printf("plus, ilength53 is longer. Clipping heads.\n"));
- if ((overlap = totallength - ilength53) < 0) {
- debug13(printf("Overlap %d is negative, so returning 0\n",overlap));
- return 0;
- } else {
- debug13(printf("Overlap is %d\n",overlap));
- clipdir = +1;
- }
+ overlap = common_left + common_right - 1;
+ debug13(printf("Overlap is %d = common_left %d + common_right %d - 1\n",
+ overlap,common_left,common_right));
+ clipdir = +1;
if (hit5->hittype == GMAP || hit3->hittype == GMAP) {
/* Revise only for paired-ends involving GMAP and when successful. Observed to be the correct action. */
@@ -2963,7 +3120,7 @@ Stage3pair_overlap (int *hardclip5_low, int *hardclip5_high, int *hardclip3_low,
/* Want to clip 5 high and 3 low */
*hardclip5_high = ilength5_high - common_shift;
*hardclip3_low = overlap - (*hardclip5_high);
- debug13(printf("Clip for ilength53 plus is hardclip5 %d..%d and hardclip3 %d..%d\n",
+ debug13(printf("Overlap clip for ilength53 plus is hardclip5 %d..%d and hardclip3 %d..%d\n",
*hardclip5_low,*hardclip5_high,*hardclip3_low,*hardclip3_high));
*hardclip5_high += hit5->trim_right + hit5->end_amb_length;
*hardclip3_low += hit3->trim_left + hit3->start_amb_length;
@@ -2985,13 +3142,10 @@ Stage3pair_overlap (int *hardclip5_low, int *hardclip5_high, int *hardclip3_low,
} else {
debug13(printf("plus, ilength35 is longer. Clipping tails.\n"));
- if ((overlap = totallength - ilength35) < 0) {
- debug13(printf("Overlap %d is negative, so returning 0\n",overlap));
- return 0;
- } else {
- debug13(printf("Overlap is %d\n",overlap));
- clipdir = -1;
- }
+ overlap = common_left + common_right - 1;
+ debug13(printf("Overlap is %d = common_left %d + common_right %d - 1\n",
+ overlap,common_left,common_right));
+ clipdir = -1;
if (hit5->hittype == GMAP || hit3->hittype == GMAP) {
/* Revise only for paired-ends involving GMAP and when successful. Observed to be the correct action. */
@@ -3001,7 +3155,7 @@ Stage3pair_overlap (int *hardclip5_low, int *hardclip5_high, int *hardclip3_low,
/* Want to clip 5 low and 3 high */
*hardclip5_low = ilength5_low + common_shift;
*hardclip3_high = overlap - (*hardclip5_low);
- debug13(printf("Clip for ilength35 plus is hardclip5 %d..%d and hardclip3 %d..%d\n",
+ debug13(printf("Overlap clip for ilength35 plus is hardclip5 %d..%d and hardclip3 %d..%d\n",
*hardclip5_low,*hardclip5_high,*hardclip3_low,*hardclip3_high));
*hardclip5_low += hit5->trim_left + hit5->start_amb_length;
*hardclip3_high += hit3->trim_right + hit3->end_amb_length;
@@ -3028,16 +3182,17 @@ Stage3pair_overlap (int *hardclip5_low, int *hardclip5_high, int *hardclip3_low,
} else {
/* minus */
+#if 0
hit5_trimmed_length = hit5->querylength - hit5->trim_left - hit5->trim_right - hit5->start_amb_length - hit5->end_amb_length;
hit3_trimmed_length = hit3->querylength - hit3->trim_left - hit3->trim_right - hit3->start_amb_length - hit3->end_amb_length;
totallength = hit5_trimmed_length + hit3_trimmed_length;
debug13(printf("totallength = %d, hit5 trimmed length = %d, hit3 trimmed length = %d\n",
totallength,hit5_trimmed_length,hit3_trimmed_length));
-
debug13(printf("original insertlength: %d, trim+amb5: %d..%d, trim+amb3: %d..%d\n",
this->insertlength,hit5->trim_left + hit5->start_amb_length,
hit5->trim_right + hit5->end_amb_length,hit3->trim_left + hit3->start_amb_length,
hit3->trim_right + hit3->end_amb_length));
+#endif
if ((common_genomicpos = pair_common_genomicpos(hit5,hit3)) == 0) {
debug13(printf("Cannot determine a common point, so returning 0\n"));
@@ -3063,13 +3218,10 @@ Stage3pair_overlap (int *hardclip5_low, int *hardclip5_high, int *hardclip3_low,
if ((ilength53 = ilength5_low + ilength3_high - 1) > (ilength35 = ilength3_low + ilength5_high - 1)) {
/* Use >, not >=, so we favor clipping heads over clipping tails in case of a tie */
debug13(printf("minus, ilength53 is longer. Clipping tails.\n"));
- if ((overlap = totallength - ilength53) < 0) {
- debug13(printf("Overlap %d is negative, so returning 0\n",overlap));
- return 0;
- } else {
- debug13(printf("Overlap is %d\n",overlap));
- clipdir = +1;
- }
+ overlap = common_left + common_right - 1;
+ debug13(printf("Overlap is %d = common_left %d + common_right %d - 1\n",
+ overlap,common_left,common_right));
+ clipdir = +1;
if (hit5->hittype == GMAP || hit3->hittype == GMAP) {
/* Revise only for paired-ends involving GMAP and when successful. Observed to be the correct action. */
@@ -3079,7 +3231,7 @@ Stage3pair_overlap (int *hardclip5_low, int *hardclip5_high, int *hardclip3_low,
/* Want to clip 5 high and 3 low */
*hardclip5_high = ilength5_high - common_shift;
*hardclip3_low = overlap - (*hardclip5_high);
- debug13(printf("Clip for ilength53 minus is hardclip5 %d..%d and hardclip3 %d..%d\n",
+ debug13(printf("Overlap clip for ilength53 minus is hardclip5 %d..%d and hardclip3 %d..%d\n",
*hardclip5_low,*hardclip5_high,*hardclip3_low,*hardclip3_high));
*hardclip5_high += hit5->trim_left + hit5->start_amb_length;
*hardclip3_low += hit3->trim_right + hit3->end_amb_length;
@@ -3100,13 +3252,10 @@ Stage3pair_overlap (int *hardclip5_low, int *hardclip5_high, int *hardclip3_low,
} else {
debug13(printf("minus, ilength35 is longer. Clipping heads.\n"));
- if ((overlap = totallength - ilength35) < 0) {
- debug13(printf("Overlap %d is negative, so returning 0\n",overlap));
- return 0;
- } else {
- debug13(printf("Overlap is %d\n",overlap));
- clipdir = -1;
- }
+ overlap = common_left + common_right - 1;
+ debug13(printf("Overlap is %d = common_left %d + common_right %d - 1\n",
+ overlap,common_left,common_right));
+ clipdir = -1;
if (hit5->hittype == GMAP || hit3->hittype == GMAP) {
/* Revise only for paired-ends involving GMAP and when successful. Observed to be the correct action. */
@@ -3116,7 +3265,7 @@ Stage3pair_overlap (int *hardclip5_low, int *hardclip5_high, int *hardclip3_low,
/* Want to clip 5 low and 3 high. Verified. */
*hardclip5_low = ilength5_low + common_shift;
*hardclip3_high = overlap - (*hardclip5_low);
- debug13(printf("Clip for ilength35 minus is hardclip5 %d..%d and hardclip3 %d..%d\n",
+ debug13(printf("Overlap clip for ilength35 minus is hardclip5 %d..%d and hardclip3 %d..%d\n",
*hardclip5_low,*hardclip5_high,*hardclip3_low,*hardclip3_high));
*hardclip5_low += hit5->trim_right + hit5->end_amb_length;
*hardclip3_high += hit3->trim_left + hit3->start_amb_length;
@@ -3416,18 +3565,43 @@ Stage3end_copy (T old) {
new->nchimera_known = old->nchimera_known;
new->nchimera_novel = old->nchimera_novel;
- new->substring1 = Substring_copy(old->substring1);
- new->substring2 = Substring_copy(old->substring2);
- new->substring0 = Substring_copy(old->substring0);
+ new->substring_LtoH = (List_T) NULL;
if (old->hittype == GMAP) {
new->pairarray = Pairpool_copy_array(old->pairarray,old->npairs);
new->npairs = old->npairs;
new->nsegments = old->nsegments;
+
+ new->substring1 = (Substring_T) NULL;
+ new->substring2 = (Substring_T) NULL;
+ new->substring0 = (Substring_T) NULL;
+
} else {
new->pairarray = (struct Pair_T *) NULL;
new->npairs = 0;
new->nsegments = 0;
+
+ new->substring1 = Substring_copy(old->substring1);
+ new->substring2 = Substring_copy(old->substring2);
+ new->substring0 = Substring_copy(old->substring0);
+
+ if (new->plusp == true) {
+ if (new->substring2 != NULL) {
+ new->substring_LtoH = List_push(new->substring_LtoH,(void *) new->substring2);
+ }
+ new->substring_LtoH = List_push(new->substring_LtoH,(void *) new->substring1);
+ if (new->substring0 != NULL) {
+ new->substring_LtoH = List_push(new->substring_LtoH,(void *) new->substring0);
+ }
+ } else {
+ if (new->substring0 != NULL) {
+ new->substring_LtoH = List_push(new->substring_LtoH,(void *) new->substring0);
+ }
+ new->substring_LtoH = List_push(new->substring_LtoH,(void *) new->substring1);
+ if (new->substring2 != NULL) {
+ new->substring_LtoH = List_push(new->substring_LtoH,(void *) new->substring2);
+ }
+ }
}
if (old->substring_donor == NULL) {
@@ -3490,43 +3664,6 @@ Stage3end_copy (T old) {
Except_raise(&Copy_Substring, __FILE__, __LINE__);
}
-
- if (old->substring_low == NULL) {
- new->substring_low = NULL;
- } else if (old->substring_low == old->substring1) {
- new->substring_low = new->substring1;
- } else if (old->substring_low == old->substring2) {
- new->substring_low = new->substring2;
- } else if (old->substring_low == old->substring0) {
- new->substring_low = new->substring0;
- } else {
- fprintf(stderr,"substring_low for type %s is not NULL, substring1, or substring2, or substring0\n",
- hittype_string(old->hittype));
- fprintf(stderr,"substring_low %p\n",old->substring_low);
- fprintf(stderr,"substring1 %p\n",old->substring1);
- fprintf(stderr,"substring2 %p\n",old->substring2);
- fprintf(stderr,"substring0 %p\n",old->substring0);
- Except_raise(&Copy_Substring, __FILE__, __LINE__);
- }
-
- if (old->substring_high == NULL) {
- new->substring_high = NULL;
- } else if (old->substring_high == old->substring1) {
- new->substring_high = new->substring1;
- } else if (old->substring_high == old->substring2) {
- new->substring_high = new->substring2;
- } else if (old->substring_high == old->substring0) {
- new->substring_high = new->substring0;
- } else {
- fprintf(stderr,"substring_high for type %s is not NULL, substring1, or substring2, or substring0\n",
- hittype_string(old->hittype));
- fprintf(stderr,"substring_high %p\n",old->substring_high);
- fprintf(stderr,"substring1 %p\n",old->substring1);
- fprintf(stderr,"substring2 %p\n",old->substring2);
- fprintf(stderr,"substring0 %p\n",old->substring0);
- Except_raise(&Copy_Substring, __FILE__, __LINE__);
- }
-
new->paired_usedp = old->paired_usedp;
new->paired_seenp = old->paired_seenp;
new->concordantp = old->concordantp;
@@ -3683,7 +3820,7 @@ Stage3end_new_exact (int *found_score, Univcoord_T left, int genomiclength, Comp
new->substring0 = (Substring_T) NULL;
new->substring_donor = new->substring_acceptor = (Substring_T) NULL;
new->substringD = new->substringA = (Substring_T) NULL;
- new->substring_low = new->substring_high = new->substring1;
+ new->substring_LtoH = List_push(NULL,(void *) new->substring1);
new->pairarray = (struct Pair_T *) NULL;
@@ -3821,7 +3958,7 @@ Stage3end_new_substitution (int *found_score, int nmismatches_whole, Univcoord_T
new->substring0 = (Substring_T) NULL;
new->substring_donor = new->substring_acceptor = (Substring_T) NULL;
new->substringD = new->substringA = (Substring_T) NULL;
- new->substring_low = new->substring_high = new->substring1;
+ new->substring_LtoH = List_push(NULL,(void *) new->substring1);
new->pairarray = (struct Pair_T *) NULL;
@@ -4013,12 +4150,10 @@ Stage3end_new_insertion (int *found_score, int nindels, int indel_pos, int nmism
new->indel_pos = indel_pos;
if (plusp == true) {
- new->substring_low = new->substring1;
- new->substring_high = new->substring2;
+ new->substring_LtoH = List_push(List_push(NULL,new->substring2),new->substring1);
new->indel_low = indel_pos;
} else {
- new->substring_low = new->substring2;
- new->substring_high = new->substring1;
+ new->substring_LtoH = List_push(List_push(NULL,new->substring1),new->substring2);
new->indel_low = querylength - indel_pos;
}
@@ -4248,12 +4383,10 @@ Stage3end_new_deletion (int *found_score, int nindels, int indel_pos, int nmisma
new->deletion = (char *) NULL;
new->indel_pos = indel_pos;
if (plusp == true) {
- new->substring_low = new->substring1;
- new->substring_high = new->substring2;
+ new->substring_LtoH = List_push(List_push(NULL,new->substring1),new->substring2);
new->indel_low = indel_pos;
} else {
- new->substring_low = new->substring2;
- new->substring_high = new->substring1;
+ new->substring_LtoH = List_push(List_push(NULL,new->substring2),new->substring1);
new->indel_low = querylength - indel_pos;
}
#endif
@@ -4395,38 +4528,18 @@ Stage3end_new_splice (int *found_score, int nmismatches_donor, int nmismatches_a
new->querylength_adj = new->querylength = querylength;
if (donor == NULL) {
-#ifdef DONORIS1
new->substring1 = copy_acceptor_p ? Substring_copy(acceptor) : acceptor;
new->substring2 = (Substring_T) NULL;
new->substring_donor = (Substring_T) NULL;
new->substring_acceptor = new->substring1;
-#else
- new->substring1 = copy_acceptor_p ? Substring_copy(acceptor) : acceptor;
- new->substring2 = (Substring_T) NULL;
- new->substring_donor = (Substring_T) NULL;
- new->substring_acceptor = new->substring1;
-#endif
} else if (acceptor == NULL) {
-#ifdef DONORIS1
new->substring1 = copy_donor_p ? Substring_copy(donor) : donor;
new->substring2 = (Substring_T) NULL;
new->substring_donor = new->substring1;
new->substring_acceptor = (Substring_T) NULL;
-#else
- new->substring1 = copy_donor_p ? Substring_copy(donor) : donor;
- new->substring2 = (Substring_T) NULL;
- new->substring_donor = new->substring1;
- new->substring_acceptor = (Substring_T) NULL;
-#endif
} else {
-#ifdef DONORIS1
- new->substring1 = copy_donor_p ? Substring_copy(donor) : donor;
- new->substring2 = copy_acceptor_p ? Substring_copy(acceptor) : acceptor;
- new->substring_donor = new->substring1;
- new->substring_acceptor = new->substring2;
-#else
if (sensedir == SENSE_FORWARD) {
new->substring1 = copy_donor_p ? Substring_copy(donor) : donor;
new->substring2 = copy_acceptor_p ? Substring_copy(acceptor) : acceptor;
@@ -4440,7 +4553,6 @@ Stage3end_new_splice (int *found_score, int nmismatches_donor, int nmismatches_a
} else {
abort();
}
-#endif
}
new->substring0 = (Substring_T) NULL;
new->substringD = new->substringA = (Substring_T) NULL;
@@ -4745,24 +4857,24 @@ Stage3end_new_splice (int *found_score, int nmismatches_donor, int nmismatches_a
if (invert_first_p == false) {
if (Substring_queryend(acceptor) > Substring_queryend(donor)) {
substring_for_concordance = acceptor;
- new->substring_low = new->substring_high = new->substring_acceptor;
+ new->substring_LtoH = List_push(NULL,(void *) new->substring_acceptor);
new->effective_chrnum = Substring_chrnum(acceptor);
new->other_chrnum = Substring_chrnum(donor);
} else {
substring_for_concordance = donor;
- new->substring_low = new->substring_high = new->substring_donor;
+ new->substring_LtoH = List_push(NULL,(void *) new->substring_donor);
new->effective_chrnum = Substring_chrnum(donor);
new->other_chrnum = Substring_chrnum(acceptor);
}
} else {
if (Substring_querystart(acceptor) < Substring_querystart(donor)) {
substring_for_concordance = acceptor;
- new->substring_low = new->substring_high = new->substring_acceptor;
+ new->substring_LtoH = List_push(NULL,(void *) new->substring_acceptor);
new->effective_chrnum = Substring_chrnum(acceptor);
new->other_chrnum = Substring_chrnum(donor);
} else {
substring_for_concordance = donor;
- new->substring_low = new->substring_high = new->substring_donor;
+ new->substring_LtoH = List_push(NULL,(void *) new->substring_donor);
new->effective_chrnum = Substring_chrnum(donor);
new->other_chrnum = Substring_chrnum(acceptor);
}
@@ -4772,24 +4884,24 @@ Stage3end_new_splice (int *found_score, int nmismatches_donor, int nmismatches_a
if (invert_second_p == false) {
if (Substring_queryend(acceptor) > Substring_queryend(donor)) {
substring_for_concordance = acceptor;
- new->substring_low = new->substring_high = new->substring_acceptor;
+ new->substring_LtoH = List_push(NULL,(void *) new->substring_acceptor);
new->effective_chrnum = Substring_chrnum(acceptor);
new->other_chrnum = Substring_chrnum(donor);
} else {
substring_for_concordance = donor;
- new->substring_low = new->substring_high = new->substring_donor;
+ new->substring_LtoH = List_push(NULL,(void *) new->substring_donor);
new->effective_chrnum = Substring_chrnum(donor);
new->other_chrnum = Substring_chrnum(acceptor);
}
} else {
if (Substring_querystart(acceptor) < Substring_querystart(donor)) {
substring_for_concordance = acceptor;
- new->substring_low = new->substring_high = new->substring_acceptor;
+ new->substring_LtoH = List_push(NULL,(void *) new->substring_acceptor);
new->effective_chrnum = Substring_chrnum(acceptor);
new->other_chrnum = Substring_chrnum(donor);
} else {
substring_for_concordance = donor;
- new->substring_low = new->substring_high = new->substring_donor;
+ new->substring_LtoH = List_push(NULL,(void *) new->substring_donor);
new->effective_chrnum = Substring_chrnum(donor);
new->other_chrnum = Substring_chrnum(acceptor);
}
@@ -4807,19 +4919,11 @@ Stage3end_new_splice (int *found_score, int nmismatches_donor, int nmismatches_a
new->other_chrnum = 0;
if (donor == NULL) {
- new->substring_low = new->substring_high = new->substring1;
+ new->substring_LtoH = List_push(NULL,(void *) new->substring1);
} else if (acceptor == NULL) {
- new->substring_low = new->substring_high = new->substring1;
+ new->substring_LtoH = List_push(NULL,(void *) new->substring1);
} else if (sensedir == SENSE_FORWARD) {
-#ifdef DONORIS1
- if (new->plusp == true) {
- new->substring_low = new->substring1; /* donor */
- new->substring_high = new->substring2; /* acceptor */
- } else {
- new->substring_low = new->substring2; /* acceptor */
- new->substring_high = new->substring1; /* donor */
- }
-#else
+#if 0
if (new->plusp == true) {
new->substring_low = new->substring_donor;
new->substring_high = new->substring_acceptor;
@@ -4827,18 +4931,18 @@ Stage3end_new_splice (int *found_score, int nmismatches_donor, int nmismatches_a
new->substring_low = new->substring_acceptor;
new->substring_high = new->substring_donor;
}
-#endif
-
- } else if (sensedir == SENSE_ANTI) {
-#ifdef DONORIS1
+#else
if (new->plusp == true) {
- new->substring_low = new->substring2; /* acceptor */
- new->substring_high = new->substring1; /* donor */
+ /* Order is donor, acceptor. Same as substring1, substring2, as expected */
+ new->substring_LtoH = List_push(List_push(NULL,(void *) new->substring_acceptor),new->substring_donor);
} else {
- new->substring_low = new->substring1; /* donor */
- new->substring_high = new->substring2; /* acceptor */
+ /* Order is acceptor, donor. Same as substring2, substring1, as expected */
+ new->substring_LtoH = List_push(List_push(NULL,(void *) new->substring_donor),new->substring_acceptor);
}
-#else
+#endif
+
+ } else if (sensedir == SENSE_ANTI) {
+#if 0
if (new->plusp == true) {
new->substring_low = new->substring_acceptor;
new->substring_high = new->substring_donor;
@@ -4846,6 +4950,14 @@ Stage3end_new_splice (int *found_score, int nmismatches_donor, int nmismatches_a
new->substring_low = new->substring_donor;
new->substring_high = new->substring_acceptor;
}
+#else
+ if (new->plusp == true) {
+ /* Order is acceptor, donor. Same as substring1, substring2, as expected */
+ new->substring_LtoH = List_push(List_push(NULL,(void *) new->substring_donor),new->substring_acceptor);
+ } else {
+ /* Order is donor, acceptor. Same as substring2, substring1, as expected */
+ new->substring_LtoH = List_push(List_push(NULL,(void *) new->substring_acceptor),new->substring_donor);
+ }
#endif
} else {
@@ -4967,8 +5079,8 @@ Stage3end_new_shortexon (int *found_score, Substring_T donor, Substring_T accept
int ignore;
new = (T) MALLOC_OUT(sizeof(*new));
- debug0(printf("Stage3end_new_shortexon %p, amb_donor %d, amb_acceptor %d\n",
- new,amb_length_donor,amb_length_acceptor));
+ debug0(printf("Stage3end_new_shortexon %p, amb_donor %d, amb_acceptor %d, sensedir %d\n",
+ new,amb_length_donor,amb_length_acceptor,sensedir));
assert(Substring_match_length_orig(donor) + Substring_match_length_orig(shortexon) + Substring_match_length_orig(acceptor) +
amb_length_donor + amb_length_acceptor == querylength);
@@ -5187,7 +5299,8 @@ Stage3end_new_shortexon (int *found_score, Substring_T donor, Substring_T accept
} else {
abort();
}
-#else
+
+#elif 0
if (new->plusp == true) {
new->substring_low = (new->substring0 != NULL ? new->substring0 : new->substring1);
new->substring_high = (new->substring2 != NULL ? new->substring2 : new->substring1);
@@ -5195,6 +5308,26 @@ Stage3end_new_shortexon (int *found_score, Substring_T donor, Substring_T accept
new->substring_low = (new->substring2 != NULL ? new->substring2 : new->substring1);
new->substring_high = (new->substring0 != NULL ? new->substring0 : new->substring1);
}
+
+#else
+ new->substring_LtoH = (List_T) NULL;
+ if (new->plusp == true) {
+ if (new->substring2 != NULL) {
+ new->substring_LtoH = List_push(new->substring_LtoH,(void *) new->substring2);
+ }
+ new->substring_LtoH = List_push(new->substring_LtoH,(void *) new->substring1);
+ if (new->substring0 != NULL) {
+ new->substring_LtoH = List_push(new->substring_LtoH,(void *) new->substring0);
+ }
+ } else {
+ if (new->substring0 != NULL) {
+ new->substring_LtoH = List_push(new->substring_LtoH,(void *) new->substring0);
+ }
+ new->substring_LtoH = List_push(new->substring_LtoH,(void *) new->substring1);
+ if (new->substring2 != NULL) {
+ new->substring_LtoH = List_push(new->substring_LtoH,(void *) new->substring2);
+ }
+ }
#endif
@@ -5407,7 +5540,7 @@ Stage3end_new_terminal (int querystart, int queryend, Univcoord_T left, Compress
new->substring0 = (Substring_T) NULL;
new->substring_donor = new->substring_acceptor = (Substring_T) NULL;
new->substringD = new->substringA = (Substring_T) NULL;
- new->substring_low = new->substring_high = new->substring1;
+ new->substring_LtoH = List_push(NULL,(void *) new->substring1);
new->pairarray = (struct Pair_T *) NULL;
@@ -5569,7 +5702,7 @@ ATAGCCCACACGTTCCCCTTAAATAAGACATCACGATGGATCACAGGTCTATCACCCTATTAACCACTCACGGGAG
new->substring0 = (Substring_T) NULL;
new->substring_donor = new->substring_acceptor = (Substring_T) NULL;
new->substringD = new->substringA = (Substring_T) NULL;
- new->substring_low = new->substring_high = (Substring_T) NULL;
+ new->substring_LtoH = (List_T) NULL;
new->pairarray = pairarray;
new->npairs = npairs;
@@ -9707,20 +9840,27 @@ Stage3end_filter_bymatch (List_T hitlist) {
static Chrpos_T
overlap5_gmap_plus (int *querypos, Chrpos_T *genomicstart, Chrpos_T *genomicend,
Stage3end_T hit5, Stage3end_T gmap) {
- debug10(printf("Entered overlap5_gmap_plus\n"));
- *genomicstart = Substring_chrstart(hit5->substring_high);
- *genomicend = Substring_chrend(hit5->substring_high);
+ Substring_T substring_high;
+
+ debug10(printf("Entered overlap5_gmap_plus with gmap %d..%d\n",
+ gmap->pairarray[0].querypos,gmap->pairarray[gmap->npairs - 1].querypos));
+ substring_high = (Substring_T) List_last_value(hit5->substring_LtoH);
+ *genomicstart = Substring_alignstart_chr(substring_high);
+ *genomicend = Substring_alignend_chr(substring_high);
return Pair_binary_search_ascending(&(*querypos),/*lowi*/0,/*highi*/gmap->npairs,gmap->pairarray,
*genomicstart,*genomicend);
}
-
static Chrpos_T
overlap3_gmap_plus (int *querypos, Chrpos_T *genomicstart, Chrpos_T *genomicend,
Stage3end_T hit3, Stage3end_T gmap) {
- debug10(printf("Entered overlap3_gmap_plus\n"));
- *genomicstart = Substring_chrstart(hit3->substring_low);
- *genomicend = Substring_chrend(hit3->substring_low);
+ Substring_T substring_low;
+
+ debug10(printf("Entered overlap3_gmap_plus with gmap %d..%d\n",
+ gmap->pairarray[0].querypos,gmap->pairarray[gmap->npairs - 1].querypos));
+ substring_low = (Substring_T) List_head(hit3->substring_LtoH);
+ *genomicstart = Substring_alignstart_chr(substring_low);
+ *genomicend = Substring_alignend_chr(substring_low);
return Pair_binary_search_ascending(&(*querypos),/*lowi*/0,/*highi*/gmap->npairs,gmap->pairarray,
*genomicstart,*genomicend);
}
@@ -9729,9 +9869,13 @@ overlap3_gmap_plus (int *querypos, Chrpos_T *genomicstart, Chrpos_T *genomicend,
static Chrpos_T
overlap5_gmap_minus (int *querypos, Chrpos_T *genomicstart, Chrpos_T *genomicend,
Stage3end_T hit5, Stage3end_T gmap) {
- debug10(printf("Entered overlap5_gmap_minus\n"));
- *genomicstart = Substring_chrstart(hit5->substring_low);
- *genomicend = Substring_chrend(hit5->substring_low);
+ Substring_T substring_low;
+
+ debug10(printf("Entered overlap5_gmap_minus with gmap %d..%d\n",
+ gmap->pairarray[0].querypos,gmap->pairarray[gmap->npairs - 1].querypos));
+ substring_low = (Substring_T) List_head(hit5->substring_LtoH);
+ *genomicstart = Substring_alignstart_chr(substring_low);
+ *genomicend = Substring_alignend_chr(substring_low);
return Pair_binary_search_descending(&(*querypos),/*lowi*/0,/*highi*/gmap->npairs,gmap->pairarray,
*genomicstart,*genomicend);
}
@@ -9739,9 +9883,13 @@ overlap5_gmap_minus (int *querypos, Chrpos_T *genomicstart, Chrpos_T *genomicend
static Chrpos_T
overlap3_gmap_minus (int *querypos, Chrpos_T *genomicstart, Chrpos_T *genomicend,
Stage3end_T hit3, Stage3end_T gmap) {
- debug10(printf("Entered overlap3_gmap_minus\n"));
- *genomicstart = Substring_chrstart(hit3->substring_high);
- *genomicend = Substring_chrend(hit3->substring_high);
+ Substring_T substring_high;
+
+ debug10(printf("Entered overlap3_gmap_minus with gmap %d..%d\n",
+ gmap->pairarray[0].querypos,gmap->pairarray[gmap->npairs - 1].querypos));
+ substring_high = (Substring_T) List_last_value(hit3->substring_LtoH);
+ *genomicstart = Substring_alignstart_chr(substring_high);
+ *genomicend = Substring_alignend_chr(substring_high);
return Pair_binary_search_descending(&(*querypos),/*lowi*/0,/*highi*/gmap->npairs,gmap->pairarray,
*genomicstart,*genomicend);
}
diff --git a/src/stage3hr.h b/src/stage3hr.h
index b1cb8e9..f1debe6 100644
--- a/src/stage3hr.h
+++ b/src/stage3hr.h
@@ -1,4 +1,4 @@
-/* $Id: stage3hr.h 154778 2014-12-06 03:32:33Z twu $ */
+/* $Id: stage3hr.h 161183 2015-03-18 17:04:33Z twu $ */
#ifndef STAGE3HR_INCLUDED
#define STAGE3HR_INCLUDED
@@ -172,15 +172,14 @@ Stage3end_substringD (T this);
extern Substring_T
Stage3end_substringA (T this);
extern Substring_T
-Stage3end_substring_low (T this);
-extern Substring_T
-Stage3end_substring_high (T this);
+Stage3end_substring_low (T this, int hardclip_low);
extern Substring_T
Stage3end_substring_containing (T this, int querypos);
extern struct Pair_T *
Stage3end_pairarray (T this);
extern int
Stage3end_npairs (T this);
+
extern Chrpos_T
Stage3end_distance (T this);
extern Chrpos_T
diff --git a/src/substring.c b/src/substring.c
index 3de8a72..2e5f828 100644
--- a/src/substring.c
+++ b/src/substring.c
@@ -1,4 +1,4 @@
-static char rcsid[] = "$Id: substring.c 154591 2014-12-04 02:00:32Z twu $";
+static char rcsid[] = "$Id: substring.c 161183 2015-03-18 17:04:33Z twu $";
#ifdef HAVE_CONFIG_H
#include <config.h>
#endif
@@ -2297,13 +2297,13 @@ Substring_genomiclength (T this) {
}
Chrpos_T
-Substring_chrstart (T this) {
- return this->genomicstart - this->chroffset;
+Substring_alignstart_chr (T this) {
+ return this->alignstart - this->chroffset;
}
Chrpos_T
-Substring_chrend (T this) {
- return this->genomicend - this->chroffset;
+Substring_alignend_chr (T this) {
+ return this->alignend - this->chroffset;
}
diff --git a/src/substring.h b/src/substring.h
index 3878c72..38651e3 100644
--- a/src/substring.h
+++ b/src/substring.h
@@ -1,4 +1,4 @@
-/* $Id: substring.h 154591 2014-12-04 02:00:32Z twu $ */
+/* $Id: substring.h 161183 2015-03-18 17:04:33Z twu $ */
#ifndef SUBSTRING_INCLUDED
#define SUBSTRING_INCLUDED
@@ -169,9 +169,9 @@ extern Chrpos_T
Substring_genomiclength (T this);
extern Chrpos_T
-Substring_chrstart (T this);
+Substring_alignstart_chr (T this);
extern Chrpos_T
-Substring_chrend (T this);
+Substring_alignend_chr (T this);
extern double
Substring_chimera_prob (T this);
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-med/gmap.git
More information about the debian-med-commit
mailing list