diff --git a/src/utility/ssw/ssw.c b/src/utility/ssw/ssw.c index d73bdfea..0300cbe2 100755 --- a/src/utility/ssw/ssw.c +++ b/src/utility/ssw/ssw.c @@ -56,7 +56,7 @@ * * Created by Mengyao Zhao on 6/22/10. * Copyright 2010 Boston College. All rights reserved. - * Version 1.2.4 + * Version 1.2.5 * Last revision by Mengyao Zhao on 2022-Apr-17. * * The lazy-F loop implementation was derived from SWPS3, which is @@ -66,8 +66,6 @@ * BSD licensed under Micharl Farrar. */ -// #include - #include #include #include @@ -310,7 +308,9 @@ static alignment_end* sw_sse2_byte (const int8_t* ref, _mm_store_si128(pvHStore + j, vH); vH = _mm_subs_epu8(vH, vGapO); vF = _mm_subs_epu8(vF, vGapE); - if (UNLIKELY(! _mm_movemask_epi8(_mm_cmpgt_epi8(vF, vH)))) goto end; + vTemp = _mm_subs_epu8(vF, vH); + vTemp = _mm_cmpeq_epi8 (vTemp, vZero); + if (UNLIKELY(_mm_movemask_epi8(vTemp) == 0xffff)) goto end; } } @@ -711,7 +711,7 @@ static cigar* banded_sw (const int8_t* ref, op = 'D'; break; default: - // Rcpp::warning("Trace back error: %d.\n", direction_line[temp1 - 1]); + // fprintf(stderr, "Trace back error: %d.\n", direction_line[temp1 - 1]); free(direction); free(h_c); free(e_b); @@ -837,7 +837,7 @@ s_align* ssw_align (const s_profile* prof, r->cigarLen = 0; r->flag = 0; if (maskLen < 15) { - // Rcpp::warning("When maskLen < 15, the function ssw_align doesn't return 2nd best alignment information.\n"); + // fprintf(stderr, "When maskLen < 15, the function ssw_align doesn't return 2nd best alignment information.\n"); } // Find the alignment scores and ending positions @@ -848,7 +848,7 @@ s_align* ssw_align (const s_profile* prof, bests = sw_sse2_word(ref, 0, refLen, readLen, weight_gapO, weight_gapE, prof->profile_word, -1, maskLen); word = 1; } else if (bests[0].score == 255) { - // Rcpp::warning("Please set 2 to the score_size parameter of the function ssw_init, otherwise the alignment results will be incorrect.\n"); + // fprintf(stderr, "Please set 2 to the score_size parameter of the function ssw_init, otherwise the alignment results will be incorrect.\n"); free(r); return NULL; } @@ -856,7 +856,7 @@ s_align* ssw_align (const s_profile* prof, bests = sw_sse2_word(ref, 0, refLen, readLen, weight_gapO, weight_gapE, prof->profile_word, -1, maskLen); word = 1; }else { - // Rcpp::warning("Please call the function ssw_init before ssw_align.\n"); + // fprintf(stderr, "Please call the function ssw_init before ssw_align.\n"); free(r); return NULL; } @@ -888,7 +888,7 @@ s_align* ssw_align (const s_profile* prof, r->read_begin1 = r->read_end1 - bests_reverse[0].read; if (UNLIKELY(r->score1 > bests_reverse[0].score)) { // banded_sw result will miss a small part - // Rcpp::warning("Warning: The alignment path of one pair of sequences may miss a small part. [ssw.c ssw_align]\n"); + // fprintf(stderr, "Warning: The alignment path of one pair of sequences may miss a small part. [ssw.c ssw_align]\n"); r->flag = 2; } free(bests_reverse); diff --git a/src/utility/ssw/ssw_cpp.cpp b/src/utility/ssw/ssw_cpp.cpp index d5e3a7ec..fe7dc998 100755 --- a/src/utility/ssw/ssw_cpp.cpp +++ b/src/utility/ssw/ssw_cpp.cpp @@ -1,6 +1,6 @@ // ssw_cpp.cpp // Created by Wan-Ping Lee -// Last revision by Mengyao Zhao on 2017-05-30 +// Last revision by Mengyao Zhao on 2023-Apr-21 #include "ssw_cpp.h" #include "ssw.h" @@ -290,25 +290,18 @@ int Aligner::SetReferenceSequence(const char* seq, const int& length) { int len = 0; if (translation_matrix_) { - // calculate the valid length - //int calculated_ref_length = static_cast(strlen(seq)); - //int valid_length = (calculated_ref_length > length) - // ? length : calculated_ref_length; - int valid_length = length; // delete the current buffer CleanReferenceSequence(); // allocate a new buffer - translated_reference_ = new int8_t[valid_length]; + translated_reference_ = new int8_t[length]; - len = TranslateBase(seq, valid_length, translated_reference_); + len = TranslateBase(seq, length, translated_reference_); } else { // nothing } reference_length_ = len; return len; - - } int Aligner::TranslateBase(const char* bases, const int& length, @@ -326,7 +319,7 @@ int Aligner::TranslateBase(const char* bases, const int& length, } -bool Aligner::Align(const char* query, const Filter& filter, +uint16_t Aligner::Align(const char* query, const Filter& filter, Alignment* alignment, const int32_t maskLen) const { if (!translation_matrix_) return false; @@ -351,18 +344,18 @@ bool Aligner::Align(const char* query, const Filter& filter, alignment->Clear(); ConvertAlignment(*s_al, query_len, alignment); alignment->mismatches = CalculateNumberMismatch(&*alignment, translated_reference_, translated_query, query_len); - + uint16_t align_flag = s_al->flag; // Free memory delete [] translated_query; align_destroy(s_al); init_destroy(profile); - return true; + return align_flag; } -bool Aligner::Align(const char* query, const char* ref, const int& ref_len, +uint16_t Aligner::Align(const char* query, const char* ref, const int& ref_len, const Filter& filter, Alignment* alignment, const int32_t maskLen) const { if (!translation_matrix_) return false; @@ -392,6 +385,7 @@ bool Aligner::Align(const char* query, const char* ref, const int& ref_len, alignment->Clear(); ConvertAlignment(*s_al, query_len, alignment); alignment->mismatches = CalculateNumberMismatch(&*alignment, translated_ref, translated_query, query_len); + uint16_t align_flag = s_al->flag; // Free memory delete [] translated_query; @@ -399,7 +393,7 @@ bool Aligner::Align(const char* query, const char* ref, const int& ref_len, align_destroy(s_al); init_destroy(profile); - return true; + return align_flag; } void Aligner::Clear(void) { diff --git a/src/utility/ssw/ssw_cpp.h b/src/utility/ssw/ssw_cpp.h index 4901352e..c312fc1a 100755 --- a/src/utility/ssw/ssw_cpp.h +++ b/src/utility/ssw/ssw_cpp.h @@ -1,6 +1,6 @@ // ssw_cpp.h // Created by Wan-Ping Lee -// Last revision by Mengyao Zhao on 2017-05-30 +// Last revision by Mengyao Zhao on 2023-Apr-21 #ifndef COMPLETE_STRIPED_SMITH_WATERMAN_CPP_H_ #define COMPLETE_STRIPED_SMITH_WATERMAN_CPP_H_ @@ -109,7 +109,7 @@ class Aligner { // and replaced. // @param seq The reference bases; // [NOTICE] It is not necessary null terminated. - // @param length The length of bases will be be built. + // @param length The length of bases will be built. // @return The length of the built bases. // ========= int SetReferenceSequence(const char* seq, const int& length); @@ -134,9 +134,9 @@ class Aligner { // @param maskLen The distance between the optimal and suboptimal alignment ending position will >= maskLen. We suggest to // use readLen/2, if you don't have special concerns. Note: maskLen has to be >= 15, otherwise this function // will NOT return the suboptimal alignment information. - // @return True: succeed; false: fail. + // @return If the alignment path is accurate (or has missing part). 0: accurate; 1: banded_sw is totally failed; 2: banded_sw returned path has missing part // ========= - bool Align(const char* query, const Filter& filter, Alignment* alignment, const int32_t maskLen) const; + uint16_t Align(const char* query, const Filter& filter, Alignment* alignment, const int32_t maskLen) const; // ========= // @function Align the query againt the reference. @@ -151,9 +151,9 @@ class Aligner { // @param maskLen The distance between the optimal and suboptimal alignment ending position will >= maskLen. We suggest to // use readLen/2, if you don't have special concerns. Note: maskLen has to be >= 15, otherwise this function // will NOT return the suboptimal alignment information. - // @return True: succeed; false: fail. + // @return If the alignment path is accurate (or has missing part). 0: accurate; 1: banded_sw is totally failed; 2: banded_sw returned path has missing part // ========= - bool Align(const char* query, const char* ref, const int& ref_len, + uint16_t Align(const char* query, const char* ref, const int& ref_len, const Filter& filter, Alignment* alignment, const int32_t maskLen) const; // @function Clear up all containers and thus the aligner is disabled.