Skip to content

D
dssp
Commit 5fb6af0d by Maarten L. Hekkelman
Options

writing sheet info

parent 0cedf7b1
Showing with 307 additions and 307 deletions
src/dssp-io.cpp
...@@ -350,8 +350,8 @@ void writeSheets(cif::datablock &db, const dssp &dssp) ...@@ -350,8 +350,8 @@ void writeSheets(cif::datablock &db, const dssp &dssp)
// write out the struct_sheet, since all info is available now // write out the struct_sheet, since all info is available now
auto &struct_sheet = db["struct_sheet"]; auto &struct_sheet = db["struct_sheet"];
auto &struct_sheet_range = db["struct_sheet_range"]; auto &struct_sheet_range = db["struct_sheet_range"];
auto &struct_sheet_order = db["struct_sheet_order"]; // auto &struct_sheet_order = db["struct_sheet_order"];
auto &struct_sheet_hbond = db["struct_sheet_hbond"]; // auto &struct_sheet_hbond = db["struct_sheet_hbond"];
// auto &pdbx_struct_sheet_hbond = db["pdbx_struct_sheet_hbond"]; // auto &pdbx_struct_sheet_hbond = db["pdbx_struct_sheet_hbond"];
// auto &pdbx_poly_seq_scheme = db["pdbx_poly_seq_scheme"]; // auto &pdbx_poly_seq_scheme = db["pdbx_poly_seq_scheme"];
...@@ -359,7 +359,7 @@ void writeSheets(cif::datablock &db, const dssp &dssp) ...@@ -359,7 +359,7 @@ void writeSheets(cif::datablock &db, const dssp &dssp)
{ {
auto sheetID = cif::cif_id_for_number(sheetNr - 1); auto sheetID = cif::cif_id_for_number(sheetNr - 1);
std::map<std::tuple<int,int>,int> sheetOrder; // std::map<std::tuple<int,int>,int> sheetOrder;
struct_sheet.emplace({ struct_sheet.emplace({
{ "id", sheetID }, { "id", sheetID },
...@@ -403,344 +403,344 @@ void writeSheets(cif::datablock &db, const dssp &dssp) ...@@ -403,344 +403,344 @@ void writeSheets(cif::datablock &db, const dssp &dssp)
{ "end_auth_asym_id", end.auth_asym_id() }, { "end_auth_asym_id", end.auth_asym_id() },
{ "end_auth_seq_id", end.auth_seq_id() } }); { "end_auth_seq_id", end.auth_seq_id() } });
// loop over the ladders for this strand and write out the struct_sheet_order // // loop over the ladders for this strand and write out the struct_sheet_order
for (auto &res : strand) // for (auto &res : strand)
{ // {
for (int i : { 0, 1 }) // for (int i : { 0, 1 })
{ // {
const auto &[bp, ladder, parallel] = res.bridge_partner(i); // const auto &[bp, ladder, parallel] = res.bridge_partner(i);
// if (not bp or bp.type() != dssp::structure_type::Strand) // // if (not bp or bp.type() != dssp::structure_type::Strand)
if (not bp) // if (not bp)
continue; // continue;
auto s1 = res.strand(); // auto s1 = res.strand();
auto s2 = bp.strand(); // auto s2 = bp.strand();
if (s1 > s2) // if (s1 > s2)
std::swap(s1, s2); // std::swap(s1, s2);
sheetOrder[{ s1, s2 }] += parallel ? 1 : -1; // sheetOrder[{ s1, s2 }] += parallel ? 1 : -1;
} // }
} // }
} }
for (const auto &[k, parallel_test] : sheetOrder) // for (const auto &[k, parallel_test] : sheetOrder)
{ // {
const auto &[s1, s2] = k; // const auto &[s1, s2] = k;
if (parallel_test == 0) // if (parallel_test == 0)
continue; // continue;
bool parallel = parallel_test > 0; // bool parallel = parallel_test > 0;
std::string strandID_1 = cif::cif_id_for_number(s1 - 1);
std::string strandID_2 = cif::cif_id_for_number(s2 - 1);
struct_sheet_order.emplace({
{ "sheet_id", sheetID },
{ "range_id_1", strandID_1 },
{ "range_id_2", strandID_2 },
{ "sense", parallel > 0 ? "parallel" : "anti-parallel" }
});
res_list strand1, strand2;
for (auto &&[strandTuple, strand] : strands)
{
const auto &[sheet_nr, strand_nr] = strandTuple;
if (sheet_nr != sheetNr)
continue;
if (strand_nr == s1)
strand1 = strand;
else if (strand_nr == s2)
{
strand2 = strand;
break;
}
}
assert(not(strand1.empty() or strand2.empty()));
int beg1SeqID = 0, beg2SeqID = 0, end1SeqID = 0, end2SeqID = 0;
std::string beg1AtomID, beg2AtomID, end1AtomID, end2AtomID;
if (parallel)
{
// I. a d II. a d parallel
// \ /
// b e b e <= the residues forming the bridge
// / \ ..
// c f c f
for (auto b : strand1)
{
for (int i : { 0, 1 })
{
const auto &[e, ladder1, parallel1] = b.bridge_partner(i);
auto esi = std::find(strand2.begin(), strand2.end(), e);
if (esi == strand2.end())
continue;
auto bi = std::find(dssp.begin(), dssp.end(), b);
assert(bi != dssp.end() and bi != dssp.begin());
auto a = *std::prev(bi);
auto c = *std::next(bi);
auto ei = std::find(dssp.begin(), dssp.end(), e);
assert(ei != dssp.end() and ei != dssp.begin());
auto d = *std::prev(ei);
auto f = *std::next(ei);
if (test_bond(e, a) and test_bond(c, e)) // case I.
{
beg1SeqID = a.seq_id();
beg2SeqID = e.seq_id();
beg1AtomID = "O";
beg2AtomID = "N";
}
else if (test_bond(b, d) and test_bond(f, b)) // case II.
{
beg1SeqID = b.seq_id();
beg2SeqID = d.seq_id();
beg1AtomID = "N";
beg2AtomID = "O";
}
break;
}
if (beg1SeqID) // std::string strandID_1 = cif::cif_id_for_number(s1 - 1);
break; // std::string strandID_2 = cif::cif_id_for_number(s2 - 1);
}
std::reverse(strand1.begin(), strand1.end()); // struct_sheet_order.emplace({
std::reverse(strand2.begin(), strand2.end()); // { "sheet_id", sheetID },
// { "range_id_1", strandID_1 },
// { "range_id_2", strandID_2 },
// { "sense", parallel > 0 ? "parallel" : "anti-parallel" }
// });
for (auto b : strand1) // res_list strand1, strand2;
{
for (int i : { 0, 1 })
{
const auto &[e, ladder1, parallel1] = b.bridge_partner(i);
auto esi = std::find(strand2.begin(), strand2.end(), e);
if (esi == strand2.end())
continue;
auto bi = std::find(dssp.begin(), dssp.end(), b);
assert(bi != dssp.end() and bi != dssp.begin());
auto a = *std::next(bi);
auto c = *std::prev(bi);
auto ei = std::find(dssp.begin(), dssp.end(), e);
assert(ei != dssp.end() and ei != dssp.begin());
auto d = *std::next(ei);
auto f = *std::prev(ei);
if (test_bond(a, e) and test_bond(e, c)) // case I.
{
end1SeqID = a.seq_id();
end2SeqID = e.seq_id();
end1AtomID = "N";
end2AtomID = "O";
}
else if (test_bond(d, b) and test_bond(b, f)) // case II.
{
end1SeqID = b.seq_id();
end2SeqID = d.seq_id();
end1AtomID = "O";
end2AtomID = "N";
}
break; // for (auto &&[strandTuple, strand] : strands)
} // {
// const auto &[sheet_nr, strand_nr] = strandTuple;
// if (sheet_nr != sheetNr)
// continue;
if (end1SeqID) // if (strand_nr == s1)
break; // strand1 = strand;
} // else if (strand_nr == s2)
} // {
else // strand2 = strand;
{ // break;
// III. a <- f IV. a f antiparallel // }
// // }
// b e b <-> e <= the residues forming the bridge
//
// c -> d c d
std::reverse(strand2.begin(), strand2.end()); // assert(not(strand1.empty() or strand2.empty()));
for (auto b : strand1) // int beg1SeqID = 0, beg2SeqID = 0, end1SeqID = 0, end2SeqID = 0;
{ // std::string beg1AtomID, beg2AtomID, end1AtomID, end2AtomID;
for (int i : { 0, 1 })
{
const auto &[e, ladder1, parallel1] = b.bridge_partner(i);
auto esi = std::find(strand2.begin(), strand2.end(), e);
if (esi == strand2.end())
continue;
auto bi = std::find(dssp.begin(), dssp.end(), b);
assert(bi != dssp.end() and bi != dssp.begin());
auto a = *std::prev(bi);
auto c = *std::next(bi);
auto ei = std::find(dssp.begin(), dssp.end(), e);
assert(ei != dssp.end() and ei != dssp.begin());
auto d = *std::prev(ei);
auto f = *std::next(ei);
if (test_bond(f, a) and test_bond(c, d)) // case III.
{
beg1SeqID = a.seq_id();
beg2SeqID = f.seq_id();
beg1AtomID = "O";
beg2AtomID = "N";
}
else if (test_bond(b, e) and test_bond(e, b)) // case IV.
{
beg1SeqID = b.seq_id();
beg2SeqID = e.seq_id();
beg1AtomID = "N";
beg2AtomID = "O";
}
break; // if (parallel)
} // {
// // I. a d II. a d parallel
// // \ /
// // b e b e <= the residues forming the bridge
// // / \ ..
// // c f c f
if (beg1SeqID) // for (auto b : strand1)
break; // {
} // for (int i : { 0, 1 })
// {
// const auto &[e, ladder1, parallel1] = b.bridge_partner(i);
// auto esi = std::find(strand2.begin(), strand2.end(), e);
// if (esi == strand2.end())
// continue;
// auto bi = std::find(dssp.begin(), dssp.end(), b);
// assert(bi != dssp.end() and bi != dssp.begin());
// auto a = *std::prev(bi);
// auto c = *std::next(bi);
// auto ei = std::find(dssp.begin(), dssp.end(), e);
// assert(ei != dssp.end() and ei != dssp.begin());
// auto d = *std::prev(ei);
// auto f = *std::next(ei);
// if (test_bond(e, a) and test_bond(c, e)) // case I.
// {
// beg1SeqID = a.seq_id();
// beg2SeqID = e.seq_id();
// beg1AtomID = "O";
// beg2AtomID = "N";
// }
// else if (test_bond(b, d) and test_bond(f, b)) // case II.
// {
// beg1SeqID = b.seq_id();
// beg2SeqID = d.seq_id();
// beg1AtomID = "N";
// beg2AtomID = "O";
// }
// break;
// }
// if (beg1SeqID)
// break;
// }
std::reverse(strand1.begin(), strand1.end()); // std::reverse(strand1.begin(), strand1.end());
std::reverse(strand2.begin(), strand2.end()); // std::reverse(strand2.begin(), strand2.end());
for (auto b : strand1) // for (auto b : strand1)
{ // {
for (int i : { 0, 1 }) // for (int i : { 0, 1 })
{ // {
const auto &[e, ladder1, parallel1] = b.bridge_partner(i); // const auto &[e, ladder1, parallel1] = b.bridge_partner(i);
auto esi = std::find(strand2.begin(), strand2.end(), e); // auto esi = std::find(strand2.begin(), strand2.end(), e);
if (esi == strand2.end()) // if (esi == strand2.end())
continue; // continue;
auto bi = std::find(dssp.begin(), dssp.end(), b); // auto bi = std::find(dssp.begin(), dssp.end(), b);
assert(bi != dssp.end() and bi != dssp.begin()); // assert(bi != dssp.end() and bi != dssp.begin());
auto a = *std::next(bi); // auto a = *std::next(bi);
auto c = *std::prev(bi); // auto c = *std::prev(bi);
auto ei = std::find(dssp.begin(), dssp.end(), e); // auto ei = std::find(dssp.begin(), dssp.end(), e);
assert(ei != dssp.end() and ei != dssp.begin()); // assert(ei != dssp.end() and ei != dssp.begin());
auto d = *std::next(ei); // auto d = *std::next(ei);
auto f = *std::prev(ei); // auto f = *std::prev(ei);
if (test_bond(a, f) and test_bond(d, c)) // case III. // if (test_bond(a, e) and test_bond(e, c)) // case I.
{ // {
end1SeqID = a.seq_id(); // end1SeqID = a.seq_id();
end2SeqID = f.seq_id(); // end2SeqID = e.seq_id();
end1AtomID = "N"; // end1AtomID = "N";
end2AtomID = "O"; // end2AtomID = "O";
} // }
else if (test_bond(b, e) and test_bond(e, b)) // case IV. // else if (test_bond(d, b) and test_bond(b, f)) // case II.
{ // {
end1SeqID = b.seq_id(); // end1SeqID = b.seq_id();
end2SeqID = e.seq_id(); // end2SeqID = d.seq_id();
end1AtomID = "N"; // end1AtomID = "O";
end2AtomID = "O"; // end2AtomID = "N";
} // }
// break;
// }
// if (end1SeqID)
// break;
// }
// }
// else
// {
// // III. a <- f IV. a f antiparallel
// //
// // b e b <-> e <= the residues forming the bridge
// //
// // c -> d c d
// std::reverse(strand2.begin(), strand2.end());
// for (auto b : strand1)
// {
// for (int i : { 0, 1 })
// {
// const auto &[e, ladder1, parallel1] = b.bridge_partner(i);
// auto esi = std::find(strand2.begin(), strand2.end(), e);
// if (esi == strand2.end())
// continue;
// auto bi = std::find(dssp.begin(), dssp.end(), b);
// assert(bi != dssp.end() and bi != dssp.begin());
// auto a = *std::prev(bi);
// auto c = *std::next(bi);
// auto ei = std::find(dssp.begin(), dssp.end(), e);
// assert(ei != dssp.end() and ei != dssp.begin());
// auto d = *std::prev(ei);
// auto f = *std::next(ei);
// if (test_bond(f, a) and test_bond(c, d)) // case III.
// {
// beg1SeqID = a.seq_id();
// beg2SeqID = f.seq_id();
// beg1AtomID = "O";
// beg2AtomID = "N";
// }
// else if (test_bond(b, e) and test_bond(e, b)) // case IV.
// {
// beg1SeqID = b.seq_id();
// beg2SeqID = e.seq_id();
// beg1AtomID = "N";
// beg2AtomID = "O";
// }
// break;
// }
// if (beg1SeqID)
// break;
// }
break; // std::reverse(strand1.begin(), strand1.end());
} // std::reverse(strand2.begin(), strand2.end());
if (end1SeqID) // for (auto b : strand1)
break; // {
} // for (int i : { 0, 1 })
} // {
// const auto &[e, ladder1, parallel1] = b.bridge_partner(i);
// auto esi = std::find(strand2.begin(), strand2.end(), e);
// if (esi == strand2.end())
// continue;
// auto bi = std::find(dssp.begin(), dssp.end(), b);
// assert(bi != dssp.end() and bi != dssp.begin());
// auto a = *std::next(bi);
// auto c = *std::prev(bi);
// auto ei = std::find(dssp.begin(), dssp.end(), e);
// assert(ei != dssp.end() and ei != dssp.begin());
// auto d = *std::next(ei);
// auto f = *std::prev(ei);
// if (test_bond(a, f) and test_bond(d, c)) // case III.
// {
// end1SeqID = a.seq_id();
// end2SeqID = f.seq_id();
// end1AtomID = "N";
// end2AtomID = "O";
// }
// else if (test_bond(b, e) and test_bond(e, b)) // case IV.
// {
// end1SeqID = b.seq_id();
// end2SeqID = e.seq_id();
// end1AtomID = "N";
// end2AtomID = "O";
// }
// break;
// }
// if (end1SeqID)
// break;
// }
// }
// struct_sheet_hbond.emplace({
// { "sheet_id", sheetID },
// { "range_id_1", strandID_1 },
// { "range_id_2", strandID_2 },
// { "range_1_beg_label_seq_id", beg1SeqID },
// { "range_1_beg_label_atom_id", beg1AtomID },
// { "range_2_beg_label_seq_id", beg2SeqID },
// { "range_2_beg_label_atom_id", beg2AtomID },
// { "range_1_end_label_seq_id", end1SeqID },
// { "range_1_end_label_atom_id", end1AtomID },
// { "range_2_end_label_seq_id", end2SeqID },
// { "range_2_end_label_atom_id", end2AtomID } });
// }
struct_sheet_hbond.emplace({ // // Data items in the PDBX_STRUCT_SHEET_HBOND category record details
{ "sheet_id", sheetID }, // // about the hydrogen bonding between residue ranges in a beta sheet.
{ "range_id_1", strandID_1 }, // // This category is provided for cases where only a single hydrogen
{ "range_id_2", strandID_2 }, // // bond is used to register the two residue ranges. Category
{ "range_1_beg_label_seq_id", beg1SeqID }, // // STRUCT_SHEET_HBOND should be used when the initial and terminal
{ "range_1_beg_label_atom_id", beg1AtomID }, // // hydrogen bonds for strand pair are known.
{ "range_2_beg_label_seq_id", beg2SeqID },
{ "range_2_beg_label_atom_id", beg2AtomID },
{ "range_1_end_label_seq_id", end1SeqID },
{ "range_1_end_label_atom_id", end1AtomID },
{ "range_2_end_label_seq_id", end2SeqID },
{ "range_2_end_label_atom_id", end2AtomID } });
}
// Data items in the PDBX_STRUCT_SHEET_HBOND category record details // // // Okay, so that means we should write entries in that category only
// about the hydrogen bonding between residue ranges in a beta sheet. // // // when the ladder is size 1. Right? But those are called BetaBridges
// This category is provided for cases where only a single hydrogen // // // in
// bond is used to register the two residue ranges. Category
// STRUCT_SHEET_HBOND should be used when the initial and terminal
// hydrogen bonds for strand pair are known.
// // Okay, so that means we should write entries in that category only
// // when the ladder is size 1. Right? But those are called BetaBridges
// // in
// // {
// {
// // auto b1 = pdbx_poly_seq_scheme.find_first("asym_id"_key == strand1.front().asym_id() and "seq_id"_key == beg1SeqID);
// // auto e1 = pdbx_poly_seq_scheme.find_first("asym_id"_key == strand1.front().asym_id() and "seq_id"_key == end1SeqID);
// // auto b2 = pdbx_poly_seq_scheme.find_first("asym_id"_key == strand2.front().asym_id() and "seq_id"_key == beg2SeqID);
// // auto e2 = pdbx_poly_seq_scheme.find_first("asym_id"_key == strand2.front().asym_id() and "seq_id"_key == end2SeqID);
// auto b1 = pdbx_poly_seq_scheme.find_first("asym_id"_key == strand1.front().asym_id() and "seq_id"_key == beg1SeqID); // // if (not(b1 and e1 and b2 and e2))
// auto e1 = pdbx_poly_seq_scheme.find_first("asym_id"_key == strand1.front().asym_id() and "seq_id"_key == end1SeqID); // // {
// auto b2 = pdbx_poly_seq_scheme.find_first("asym_id"_key == strand2.front().asym_id() and "seq_id"_key == beg2SeqID); // // if (cif::VERBOSE > 0)
// auto e2 = pdbx_poly_seq_scheme.find_first("asym_id"_key == strand2.front().asym_id() and "seq_id"_key == end2SeqID); // // std::cerr << "error looking up strand ends" << std::endl;
// // continue;
// // }
// if (not(b1 and e1 and b2 and e2)) // // pdbx_struct_sheet_hbond.emplace({
// { // // { "sheet_id", sheetID },
// if (cif::VERBOSE > 0) // // { "range_id_1", strandID_1 },
// std::cerr << "error looking up strand ends" << std::endl; // // { "range_id_2", strandID_2 },
// continue;
// }
// pdbx_struct_sheet_hbond.emplace({ // // { "range_1_label_atom_id", beg1AtomID },
// { "sheet_id", sheetID }, // // { "range_1_label_comp_id", b1["mon_id"].as<std::string>() },
// { "range_id_1", strandID_1 }, // // { "range_1_label_asym_id", b1["asym_id"].as<std::string>() },
// { "range_id_2", strandID_2 }, // // { "range_1_label_seq_id", b1["seq_id"].as<std::string>() },
// // { "range_1_PDB_ins_code", b1["pdb_ins_code"].as<std::string>() },
// { "range_1_label_atom_id", beg1AtomID }, // // { "range_1_auth_atom_id", beg1AtomID },
// { "range_1_label_comp_id", b1["mon_id"].as<std::string>() }, // // { "range_1_auth_comp_id", b1["auth_mon_id"].as<std::string>() },
// { "range_1_label_asym_id", b1["asym_id"].as<std::string>() }, // // { "range_1_auth_asym_id", b1["pdb_strand_id"].as<std::string>() },
// { "range_1_label_seq_id", b1["seq_id"].as<std::string>() }, // // { "range_1_auth_seq_id", b1["auth_seq_num"].as<std::string>() },
// { "range_1_PDB_ins_code", b1["pdb_ins_code"].as<std::string>() }, // // { "range_2_label_atom_id", beg2AtomID },
// { "range_1_auth_atom_id", beg1AtomID }, // // { "range_2_label_comp_id", b2["mon_id"].as<std::string>() },
// { "range_1_auth_comp_id", b1["auth_mon_id"].as<std::string>() }, // // { "range_2_label_asym_id", b2["asym_id"].as<std::string>() },
// { "range_1_auth_asym_id", b1["pdb_strand_id"].as<std::string>() }, // // { "range_2_label_seq_id", b2["seq_id"].as<std::string>() },
// { "range_1_auth_seq_id", b1["auth_seq_num"].as<std::string>() }, // // { "range_2_PDB_ins_code", b2["pdb_ins_code"].as<std::string>() },
// { "range_2_label_atom_id", beg2AtomID }, // // { "range_2_auth_atom_id", beg2AtomID },
// { "range_2_label_comp_id", b2["mon_id"].as<std::string>() }, // // { "range_2_auth_comp_id", b2["auth_mon_id"].as<std::string>() },
// { "range_2_label_asym_id", b2["asym_id"].as<std::string>() }, // // { "range_2_auth_asym_id", b2["pdb_strand_id"].as<std::string>() },
// { "range_2_label_seq_id", b2["seq_id"].as<std::string>() }, // // { "range_2_auth_seq_id", b2["auth_seq_num"].as<std::string>() },
// { "range_2_PDB_ins_code", b2["pdb_ins_code"].as<std::string>() }, // // });
// { "range_2_auth_atom_id", beg2AtomID }, // // }
// { "range_2_auth_comp_id", b2["auth_mon_id"].as<std::string>() },
// { "range_2_auth_asym_id", b2["pdb_strand_id"].as<std::string>() },
// { "range_2_auth_seq_id", b2["auth_seq_num"].as<std::string>() },
// });
// }
} }
} }
......
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