Skip to content

Commit

Permalink
update ssw
Browse files Browse the repository at this point in the history
  • Loading branch information
@ChangqingW
ChangqingW committed Mar 27, 2024
1 parent 1638326 commit 94fa129
Show file tree
Hide file tree
Showing 3 changed files with 24 additions and 30 deletions.
18 changes: 9 additions & 9 deletions src/utility/ssw/ssw.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand All @@ -66,8 +66,6 @@
* BSD licensed under Micharl Farrar.
*/

// #include <Rcpp.h>

#include <stdint.h>
#include <stdlib.h>
#include <stdio.h>
Expand Down Expand Up @@ -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;
}
}

Expand Down Expand Up @@ -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);
Expand Down Expand Up @@ -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
Expand All @@ -848,15 +848,15 @@ 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;
}
}else if (prof->profile_word) {
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;
}
Expand Down Expand Up @@ -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);
Expand Down
24 changes: 9 additions & 15 deletions src/utility/ssw/ssw_cpp.cpp
View file
Original file line number Diff line number Diff line change
@@ -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"
Expand Down Expand Up @@ -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<int>(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,
Expand All @@ -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;
Expand All @@ -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;
Expand Down Expand Up @@ -392,14 +385,15 @@ 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;
delete [] translated_ref;
align_destroy(s_al);
init_destroy(profile);

return true;
return align_flag;
}

void Aligner::Clear(void) {
Expand Down
12 changes: 6 additions & 6 deletions src/utility/ssw/ssw_cpp.h
View file
Original file line number Diff line number Diff line change
@@ -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_
Expand Down Expand Up @@ -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);
Expand All @@ -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.
Expand All @@ -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.
Expand Down

0 comments on commit 94fa129

Please sign in to comment.