ap_PairwiseCrmsd

ap_PairwiseCrmsd calculates pairwise crmsd values for a set of protein structures (at least two).

This example evaluates crmsd for each pair of proteins twice: on C-alpha atoms and on all backbone atoms

USAGE:
ap_PairwiseCrmsd structureA.pdb structureB.pdb [structureC.pdb ... ]

Categories:

  • core::protocols::PairwiseCrmsd

Output files:

Program source:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
#include <iostream>

#include <core/protocols/PairwiseCrmsd.hh>
#include <core/data/io/Pdb.hh>
#include <core/data/structural/structure_selectors.hh>

#include <utils/exit.hh>

std::string program_info = R"(

ap_PairwiseCrmsd calculates pairwise crmsd values for a set of protein structures (at least two).

This example evaluates crmsd for each pair of proteins twice: on C-alpha atoms and on all backbone atoms

USAGE:
    ap_PairwiseCrmsd structureA.pdb structureB.pdb [structureC.pdb ... ]

)";

/** @brief Calculates crmsd value for a set of protein structures (at least two)
 *
 * CATEGORIES: core::protocols::PairwiseCrmsd
 * KEYWORDS:   PDB input; crmsd; structure selectors
 */
int main(const int argc, const char *argv[]) {

  if(argc < 3) utils::exit_OK_with_message(program_info); // --- complain about missing program parameter

  using core::data::basic::Vec3;
  using namespace core::data::structural; // --- for all AtomSelector types
  using namespace core::data::io;
  using namespace core::protocols;

  std::vector<core::data::structural::Structure_SP> structures;
  std::vector<std::string> tags;
  for (int i = 1; i < argc; ++i) {
    core::data::io::Pdb reader(argv[i],all_true(is_not_alternative,is_not_water), false); // --- note we read all atoms but skip alternate locators and waters
    structures.push_back(reader.create_structure(0));
    tags.push_back(structures.back()->code());
  }

  // ---------- crmsd on C-alpha : this is the
  std::cout <<"# crmsd on alpha carbons:\n";
  std::shared_ptr<AtomSelector> is_CA = std::make_shared<IsCA>();
  PairwiseCrmsd rmsd_ca(structures,is_CA,tags);
  rmsd_ca.crmsd_cutoff(20.0); // crmsd cutoff large anough to get some output
  rmsd_ca.calculate();

  // ---------- crmsd on backbone
  std::cout <<"# crmsd on heavy backbone atoms:\n";
  std::shared_ptr<AtomSelector> is_bb = std::make_shared<IsBB>();
  std::shared_ptr<AtomSelector> not_h = std::make_shared<NotHydrogen>();
  std::shared_ptr<LogicalANDSelector> heavy_bb = std::make_shared<LogicalANDSelector>();
  heavy_bb->add_selector(is_bb);
  heavy_bb->add_selector(not_h);
  PairwiseCrmsd rmsd_bb(structures,heavy_bb,tags);
  rmsd_bb.crmsd_cutoff(20.0); // crmsd cutoff large anough to get some output
  rmsd_bb.calculate();
}
../_images/file_icon.png