Commit 62face74 by Maarten L. Hekkelman

catch when another thread is still running

parent 3928be49
...@@ -1211,98 +1211,110 @@ DSSPImpl::DSSPImpl(const Structure& s, int min_poly_proline_stretch_length) ...@@ -1211,98 +1211,110 @@ DSSPImpl::DSSPImpl(const Structure& s, int min_poly_proline_stretch_length)
std::thread ta(std::bind(&CalculateAccessibilities, std::ref(mResidues), std::ref(mStats))); std::thread ta(std::bind(&CalculateAccessibilities, std::ref(mResidues), std::ref(mStats)));
auto& db = s.getFile().data(); try
for (auto r: db["struct_conn"].find(cif::Key("conn_type_id") == "disulf"))
{ {
std::string asym1, asym2; auto& db = s.getFile().data();
int seq1, seq2; for (auto r: db["struct_conn"].find(cif::Key("conn_type_id") == "disulf"))
cif::tie(asym1, seq1, asym2, seq2) = r.get("ptnr1_label_asym_id", "ptnr1_label_seq_id", "ptnr2_label_asym_id", "ptnr2_label_seq_id"); {
std::string asym1, asym2;
int seq1, seq2;
cif::tie(asym1, seq1, asym2, seq2) = r.get("ptnr1_label_asym_id", "ptnr1_label_seq_id", "ptnr2_label_asym_id", "ptnr2_label_seq_id");
auto r1 = findRes(asym1, seq1); auto r1 = findRes(asym1, seq1);
if (r1 == mResidues.end()) if (r1 == mResidues.end())
throw std::runtime_error("Invalid file, missing residue for SS bond"); throw std::runtime_error("Invalid file, missing residue for SS bond");
auto r2 = findRes(asym2, seq2); auto r2 = findRes(asym2, seq2);
if (r2 == mResidues.end()) if (r2 == mResidues.end())
throw std::runtime_error("Invalid file, missing residue for SS bond"); throw std::runtime_error("Invalid file, missing residue for SS bond");
mSSBonds.emplace_back(&*r1, &*r2); mSSBonds.emplace_back(&*r1, &*r2);
} }
if (cif::VERBOSE) std::cerr << "."; if (cif::VERBOSE) std::cerr << ".";
CalculateHBondEnergies(mResidues); CalculateHBondEnergies(mResidues);
if (cif::VERBOSE) std::cerr << "."; if (cif::VERBOSE) std::cerr << ".";
CalculateBetaSheets(mResidues, mStats); CalculateBetaSheets(mResidues, mStats);
if (cif::VERBOSE) std::cerr << "."; if (cif::VERBOSE) std::cerr << ".";
CalculateAlphaHelices(mResidues, mStats); CalculateAlphaHelices(mResidues, mStats);
if (cif::VERBOSE) std::cerr << "."; if (cif::VERBOSE) std::cerr << ".";
CalculatePPHelices(mResidues, mStats, m_min_poly_proline_stretch_length); CalculatePPHelices(mResidues, mStats, m_min_poly_proline_stretch_length);
if (cif::VERBOSE) std::cerr << std::endl; if (cif::VERBOSE) std::cerr << std::endl;
if (cif::VERBOSE > 1) if (cif::VERBOSE > 1)
{
for (auto& r: mResidues)
{ {
auto& m = r.mM; for (auto& r: mResidues)
char helix[5] = { };
for (HelixType helixType: { HelixType::rh_3_10, HelixType::rh_alpha, HelixType::rh_pi, HelixType::rh_pp })
{ {
switch (r.GetHelixFlag(helixType)) auto& m = r.mM;
char helix[5] = { };
for (HelixType helixType: { HelixType::rh_3_10, HelixType::rh_alpha, HelixType::rh_pi, HelixType::rh_pp })
{ {
case Helix::Start: helix[static_cast<int>(helixType)] = '>'; break; switch (r.GetHelixFlag(helixType))
case Helix::Middle: helix[static_cast<int>(helixType)] = helixType == HelixType::rh_pp ? 'P' : '3' + static_cast<int>(helixType); break; {
case Helix::StartAndEnd: helix[static_cast<int>(helixType)] = 'X'; break; case Helix::Start: helix[static_cast<int>(helixType)] = '>'; break;
case Helix::End: helix[static_cast<int>(helixType)] = '<'; break; case Helix::Middle: helix[static_cast<int>(helixType)] = helixType == HelixType::rh_pp ? 'P' : '3' + static_cast<int>(helixType); break;
case Helix::None: helix[static_cast<int>(helixType)] = ' '; break; case Helix::StartAndEnd: helix[static_cast<int>(helixType)] = 'X'; break;
case Helix::End: helix[static_cast<int>(helixType)] = '<'; break;
case Helix::None: helix[static_cast<int>(helixType)] = ' '; break;
}
} }
auto id = m.asymID() + ':' + std::to_string(m.seqID()) + '/' + m.compoundID();
std::cerr << id << std::string(12 - id.length(), ' ')
<< char(r.mSecondaryStructure) << ' '
<< helix
<< std::endl;
} }
auto id = m.asymID() + ':' + std::to_string(m.seqID()) + '/' + m.compoundID();
std::cerr << id << std::string(12 - id.length(), ' ')
<< char(r.mSecondaryStructure) << ' '
<< helix
<< std::endl;
} }
}
// finish statistics // finish statistics
mStats.nrOfSSBridges = mSSBonds.size(); mStats.nrOfSSBridges = mSSBonds.size();
mStats.nrOfIntraChainSSBridges = 0; mStats.nrOfIntraChainSSBridges = 0;
int ssBondNr = 0; int ssBondNr = 0;
for (const auto& [a, b]: mSSBonds) for (const auto& [a, b]: mSSBonds)
{ {
if (a == b) if (a == b)
throw std::runtime_error("first and second residue are the same"); {
std::cerr << "In the SS bonds list, the residue " << a->mM << " is bonded to itself" << std::endl;
if (a->mM.asymID() == b->mM.asymID() and NoChainBreak(a, b)) continue;
++mStats.nrOfIntraChainSSBridges; }
// throw std::runtime_error("first and second residue are the same");
a->mSSBridgeNr = b->mSSBridgeNr = ++ssBondNr; if (a->mM.asymID() == b->mM.asymID() and NoChainBreak(a, b))
} ++mStats.nrOfIntraChainSSBridges;
mStats.nrOfHBonds = 0; a->mSSBridgeNr = b->mSSBridgeNr = ++ssBondNr;
for (auto& r: mResidues) }
{
auto donor = r.mHBondDonor;
for (int i = 0; i < 2; ++i) mStats.nrOfHBonds = 0;
for (auto& r: mResidues)
{ {
if (donor[i].residue != nullptr and donor[i].energy < kMaxHBondEnergy) auto donor = r.mHBondDonor;
for (int i = 0; i < 2; ++i)
{ {
++mStats.nrOfHBonds; if (donor[i].residue != nullptr and donor[i].energy < kMaxHBondEnergy)
auto k = donor[i].residue->mNumber - r.mNumber; {
if (k >= -5 and k <= 5) ++mStats.nrOfHBonds;
mStats.nrOfHBondsPerDistance[k + 5] += 1; auto k = donor[i].residue->mNumber - r.mNumber;
if (k >= -5 and k <= 5)
mStats.nrOfHBondsPerDistance[k + 5] += 1;
}
} }
} }
} }
catch (const std::exception& ex)
{
ta.join();
throw;
}
ta.join(); ta.join();
} }
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment