Skip to content

Commit

Permalink
Merge pull request #362 from rachelse/master
Browse files Browse the repository at this point in the history 
Change in scoremultimer reflected in multimercluster
  • Loading branch information
@sooyoung-cha
sooyoung-cha authored Oct 21, 2024
2 parents 6b00c5d + 3d85d5c commit e1f38a1
Show file tree
Hide file tree
Showing 2 changed files with 45 additions and 22 deletions.
2 changes: 1 addition & 1 deletion README.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -325,7 +325,7 @@ The `easy-multimercluster` module is designed for multimer-level structural clus

Make sure chain names in PDB/mmcIF files does not contain underscores(_).

foldseek easy-multimercluster example/ clu tmp --multimer-tm-threshold 0.65 --chain-tm-threshold 0.5 --interface-tm-threshold 0.65
foldseek easy-multimercluster example/ clu tmp --multimer-tm-threshold 0.65 --chain-tm-threshold 0.5 --interface-lddt-threshold 0.65

#### Output MultimerCluster
##### Tab-separated multimercluster
Expand Down
65 changes: 44 additions & 21 deletions src/strucclustutils/filtermultimer.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -551,13 +551,13 @@ localThreads = std::max(std::min((size_t)par.threads, alnDbr.getSize()), (size_t
getlookupInfo(tLookupFile, tChainKeyToComplexIdMap, tComplexes, tComplexIdToIdx);
getComplexResidueLength(tDbr, tComplexes);
}
std::vector<unsigned int> qComplexOrder(qComplexes.size());
for (size_t qComplexIdx = 0; qComplexIdx < qComplexes.size(); qComplexIdx++) {
qComplexOrder[qComplexIdx] = qComplexIdx;
}
std::sort(qComplexOrder.begin(), qComplexOrder.end(), [&qComplexes](unsigned int lhs, unsigned int rhs) {
return qComplexes[lhs].chainKeys.size() > qComplexes[rhs].chainKeys.size();
});
// std::vector<unsigned int> qComplexOrder(qComplexes.size());
// for (size_t qComplexIdx = 0; qComplexIdx < qComplexes.size(); qComplexIdx++) {
// qComplexOrder[qComplexIdx] = qComplexIdx;
// }
// std::sort(qComplexOrder.begin(), qComplexOrder.end(), [&qComplexes](unsigned int lhs, unsigned int rhs) {
// return qComplexes[lhs].chainKeys.size() > qComplexes[rhs].chainKeys.size();
// });

#pragma omp parallel num_threads(localThreads)
{
Expand All @@ -579,8 +579,8 @@ localThreads = std::max(std::min((size_t)par.threads, alnDbr.getSize()), (size_t

Matcher::result_t res;
#pragma omp for schedule(dynamic, 1)
// for (size_t qComplexIdx = 0; qComplexIdx < qComplexes.size(); qComplexIdx++) {
for (size_t qComplexIdx : qComplexOrder) {
for (size_t qComplexIdx = 0; qComplexIdx < qComplexes.size(); qComplexIdx++) {
// for (size_t qComplexIdx : qComplexOrder) {
progress.updateProgress();
Complex qComplex = qComplexes[qComplexIdx];
unsigned int qComplexId = qComplex.complexId;
Expand All @@ -590,13 +590,9 @@ localThreads = std::max(std::min((size_t)par.threads, alnDbr.getSize()), (size_t
unsigned int qChainAlnId = alnDbr.getId(qChainKey);
unsigned int qChainDbId = qDbr->sequenceReader->getId(qChainKey);
// Handling monomer as singleton
// if (qChainAlnId == NOT_AVAILABLE_CHAIN_KEY) {
// char *outpos = Itoa::u32toa_sse2(qComplexId, buffer);
// result.append(buffer, (outpos - buffer - 1));
// result.push_back('\n');
// resultWriter.writeData(result.c_str(), result.length(), qComplexId, thread_idx);
// break;
// }
if (qChainAlnId == NOT_AVAILABLE_CHAIN_KEY) {
break;
}

char *qcadata = qStructDbr.getData(qChainDbId, thread_idx);
size_t qCaLength = qStructDbr.getEntryLen(qChainDbId);
Expand All @@ -619,11 +615,11 @@ localThreads = std::max(std::min((size_t)par.threads, alnDbr.getSize()), (size_t
unsigned int tComplexId = tChainKeyToComplexIdMap.at(tChainKey);
unsigned int tComplexIdx = tComplexIdToIdx.at(tComplexId);
std::vector<unsigned int> tChainKeys = tComplexes[tComplexIdx].chainKeys;
//if target is monomer, break to be singleton
// unsigned int tChainAlnId = alnDbr.getId(tChainKey);
// if (tChainAlnId == NOT_AVAILABLE_CHAIN_KEY) {
// continue;
// }
//if target is monomer, but user doesn't want, continue
unsigned int tChainAlnId = alnDbr.getId(tChainKey);
if (tChainAlnId == NOT_AVAILABLE_CHAIN_KEY) {
continue;
}
float u[3][3];
float t[3];
fillUArr(retComplex.uString, u);
Expand Down Expand Up @@ -735,6 +731,33 @@ localThreads = std::max(std::min((size_t)par.threads, alnDbr.getSize()), (size_t
tmpBuff = filterToBuffer(cmplfiltcrit, tmpBuff);
resultWrite5.writeAdd(buffer2, tmpBuff - buffer2, thread_idx);
}
if (selectedAssIDs.size() == 0) {
float t[3];
float u[3][3];
for (int i=0; i < 3; i++) {
t[i] = 0.0;
}
for (int i=0; i < 3; i++) {
for (int j=0; j < 3; j++) {
u[i][j] = 0.0;
}
}
ComplexFilterCriteria cmplfiltcrit(qComplexId, 1.0, 1.0, t, u);
cmplfiltcrit.qCov = 1.0;
cmplfiltcrit.tCov = 1.0;
cmplfiltcrit.interfaceLddt = 1.0;
resultWrite5.writeStart(thread_idx);
char * tmpBuff = Itoa::u32toa_sse2(qComplexId, buffer2);
tmpBuff = filterToBuffer(cmplfiltcrit, tmpBuff);
resultWrite5.writeAdd(buffer2, tmpBuff - buffer2, thread_idx);

char *outpos = Itoa::u32toa_sse2(qComplexId, buffer);
result.append(buffer, (outpos - buffer - 1));
result.push_back('\n');
}
// if (qComplexId == 1) {
// Debug(Debug::WARNING)<< "hi\n";
// }
resultWriter.writeData(result.c_str(), result.length(), qComplexId, thread_idx);
resultWrite5.writeEnd(qComplexId, thread_idx);
result.clear();
Expand Down

0 comments on commit e1f38a1

Please sign in to comment.