revert changes

src/dssp-io.cpp

@@ -379,15 +379,10 @@ void writeSheets(cif::datablock &db, const dssp &dssp)
 	int lastSheet = -1;
 	int lastSheet = -1;
 	for (const auto &[sheetNr, strand] : strands)
 	for (const auto &[sheetNr, strand] : strands)
 	{
 	{
-		if (sheetNr == lastSheet)
-			continue;
-
-		auto id = cif::cif_id_for_number(sheetNr);
-
-		if (not struct_sheet.exists(cif::key("id") == id))
+		if (sheetNr != lastSheet)
 		{
 		{
 			struct_sheet.emplace({
 			struct_sheet.emplace({
-				{ "id", id },
+				{ "id", cif::cif_id_for_number(sheetNr) },
 				{ "number_strands",
 				{ "number_strands",
 					std::count_if(strands.begin(), strands.end(), [nr = sheetNr](std::tuple<int, res_list> const &s)
 					std::count_if(strands.begin(), strands.end(), [nr = sheetNr](std::tuple<int, res_list> const &s)
 						{ return std::get<0>(s) == nr; })
 						{ return std::get<0>(s) == nr; })

src/dssp.cpp

@@ -768,15 +768,14 @@ double CalculateHBondEnergy(residue &inDonor, residue &inAcceptor)
 
 
 // --------------------------------------------------------------------
 // --------------------------------------------------------------------
 
 
-void CalculateHBondEnergies(std::vector<residue> &inResidues, queue_type &q1, cif::progress_bar *progress)
+void CalculateHBondEnergies(std::vector<residue> &inResidues, std::vector<std::tuple<uint32_t, uint32_t>> &q)
 {
 {
-	for (;;)
-	{
-		const auto &[i, j] = q1.pop();
-
-		if (i == 0 and j == 0)
-			break;
+	std::unique_ptr<cif::progress_bar> progress;
+	if (cif::VERBOSE == 0 or cif::VERBOSE == 1)
+		progress.reset(new cif::progress_bar(q.size(), "calculate hbond energies"));
 
 
+	for (const auto &[i, j] : q)
+	{
 		auto &ri = inResidues[i];
 		auto &ri = inResidues[i];
 		auto &rj = inResidues[j];
 		auto &rj = inResidues[j];
 
 
@@ -1554,46 +1553,48 @@ void DSSP_impl::calculateSecondaryStructure()
 		mSSBonds.emplace_back(&*r1, &*r2);
 		mSSBonds.emplace_back(&*r1, &*r2);
 	}
 	}
 
 
-	// Calculate the HBond energies
+	// Prefetch the c-alpha positions. No, really, that might be the trick
 
 
-	queue_type q1, q2;
-	std::vector<std::tuple<uint32_t,uint32_t>> near;
+	std::vector<point> cAlphas;
+	cAlphas.reserve(mResidues.size());
+	for (auto &r : mResidues)
+		cAlphas.emplace_back(r.mCAlpha);
 
 
 	std::unique_ptr<cif::progress_bar> progress;
 	std::unique_ptr<cif::progress_bar> progress;
 	if (cif::VERBOSE == 0 or cif::VERBOSE == 1)
 	if (cif::VERBOSE == 0 or cif::VERBOSE == 1)
-		progress.reset(new cif::progress_bar((mResidues.size() * (mResidues.size() - 1) / 2), "calculate hbond energies"));
+		progress.reset(new cif::progress_bar(mResidues.size() - 1, "calculate distances"));
 
 
-	std::thread hbond_thread(std::bind(&CalculateHBondEnergies, std::ref(mResidues), std::ref(q1), progress.get()));
+	// Calculate the HBond energies
+	std::vector<std::tuple<uint32_t,uint32_t>> near;
 
 
 	for (uint32_t i = 0; i + 1 < mResidues.size(); ++i)
 	for (uint32_t i = 0; i + 1 < mResidues.size(); ++i)
 	{
 	{
-		auto &ri = mResidues[i];
+		auto cai = cAlphas[i];
 
 
 		for (uint32_t j = i + 1; j < mResidues.size(); ++j)
 		for (uint32_t j = i + 1; j < mResidues.size(); ++j)
 		{
 		{
-			auto &rj = mResidues[j];
+			auto caj = cAlphas[j];
+
+			if (distance_sq(cai, caj) > (kMinimalCADistance * kMinimalCADistance))
+				continue;
 
 
-			if (distance_sq(ri.mCAlpha, rj.mCAlpha) < kMinimalCADistance * kMinimalCADistance)
-			{
-				q1.push({ i, j });
 			near.emplace_back(i, j);
 			near.emplace_back(i, j);
 		}
 		}
-		}
-	}
 
 
-	q1.push({0, 0});
+		if (progress)
+			progress->consumed(1);
+	}
 
 
-	hbond_thread.join();
+	if (cif::VERBOSE > 0)
+		std::cerr << "Considering " << near.size() << " pairs of residues" << std::endl;
 
 
-	// std::thread bsheet_thread(std::bind(&CalculateBetaSheets, std::ref(mResidues), std::ref(mStats), std::ref(near)));
+	progress.reset(nullptr);
 
 
-	// CalculateHBondEnergies(mResidues);
+	CalculateHBondEnergies(mResidues, near);
 	CalculateBetaSheets(mResidues, mStats, near);
 	CalculateBetaSheets(mResidues, mStats, near);
 	CalculateAlphaHelices(mResidues, mStats);
 	CalculateAlphaHelices(mResidues, mStats);
 	CalculatePPHelices(mResidues, mStats, m_min_poly_proline_stretch_length);
 	CalculatePPHelices(mResidues, mStats, m_min_poly_proline_stretch_length);
 
 
-	// bsheet_thread.join();
-
 	if (cif::VERBOSE > 1)
 	if (cif::VERBOSE > 1)
 	{
 	{
 		for (auto &r : mResidues)
 		for (auto &r : mResidues)