Added computation of PDBDistanceGrid for all atoms of the backbone (not just alpha-carbons)
--- a/examples/grid/pdbdistance-vineyard.cpp Sat Jan 20 13:39:51 2007 -0500
+++ b/examples/grid/pdbdistance-vineyard.cpp Wed Jan 24 17:50:45 2007 -0500
@@ -30,11 +30,12 @@
if (argc < 5)
{
- std::cout << "Usage: pdbdistance FILENAME LASTFRAME LASTSUBFRAME OUTFILENAME" << std::endl;
+ std::cout << "Usage: pdbdistance FILENAME LASTFRAME LASTSUBFRAME OUTFILENAME [CAs_ONLY]" << std::endl;
std::cout << " FILENAME - prefix of the filenames of the PDB frames" << std::endl;
std::cout << " LASTFRAME - the last frame number" << std::endl;
std::cout << " LASTSUBFRAME - the last subframe number" << std::endl;
std::cout << " OUTFILENAME - filename prefix for the resulting vineyards" << std::endl;
+ std::cout << " CAs_ONLY - only use alpha carbons [1 = true, 0 = false, default: 1]" << std::endl;
std::cout << std::endl;
std::cout << "Computes a vineyard of the pairwise distance function for a sequence of PDB frames." << std::endl;
std::cout << "Frames are in files FILENAME#1_#2.pdb, where #1 goes from 0 to LASTFRAME, " << std::endl;
@@ -45,11 +46,14 @@
int lastframe; std::istringstream(argv[2]) >> lastframe;
int lastsubframe; std::istringstream(argv[3]) >> lastsubframe;
std::string outfilename = argv[4];
+ bool cas_only = true;
+ if (argc > 5)
+ std::istringstream(argv[5]) >> cas_only;
// Compute initial filtration
int f = 0; int sf = 0;
std::ifstream in(frame_filename(infilename, f, sf++).c_str());
- Grid2DVineyard v(new PDBDistanceGrid(in));
+ Grid2DVineyard v(new PDBDistanceGrid(in, cas_only));
in.close();
std::cout << "Filtration generated, size: " << v.filtration()->size() << std::endl;
v.compute_pairing();
@@ -61,7 +65,7 @@
std::string fn = frame_filename(infilename, f, sf++);
std::cout << "Processing " << fn << std::endl;
in.open(fn.c_str());
- v.compute_vineyard(new PDBDistanceGrid(in));
+ v.compute_vineyard(new PDBDistanceGrid(in, cas_only));
in.close();
if (sf == lastsubframe) { sf = 0; ++f; }
}
--- a/examples/grid/pdbdistance.h Sat Jan 20 13:39:51 2007 -0500
+++ b/examples/grid/pdbdistance.h Wed Jan 24 17:50:45 2007 -0500
@@ -28,23 +28,35 @@
PDBDistanceGrid()
{}
- PDBDistanceGrid(std::istream& in)
+ PDBDistanceGrid(std::istream& in, bool ca_only = true)
{
- load_stream(in);
+ load_stream(in, ca_only);
}
- void load_stream(std::istream& in)
+ void load_stream(std::istream& in, bool ca_only = true)
{
dsrpdb::Protein p(in);
- std::vector<dsrpdb::Point> CAs(ca_coordinates_begin(p), ca_coordinates_end(p));
- std::cout << "CAs created, size: " << CAs.size() << std::endl;
+ typedef std::vector<dsrpdb::Point> PointVector;
+ PointVector coordinates;
+ if (ca_only)
+ {
+ PointVector v(ca_coordinates_begin(p), ca_coordinates_end(p));
+ coordinates.swap(v);
+ }
+ else
+ {
+ PointVector v(backbone_coordinates_begin(p), backbone_coordinates_end(p));
+ coordinates.swap(v);
+ }
- Grid2D::change_dimensions(CAs.size(), CAs.size());
- for (Grid2D::CoordinateIndex i = 0; i < CAs.size(); ++i)
- for (Grid2D::CoordinateIndex j = 0; j < CAs.size(); ++j)
+ std::cout << "Coordinatess created, size: " << coordinates.size() << std::endl;
+
+ Grid2D::change_dimensions(coordinates.size(), coordinates.size());
+ for (Grid2D::CoordinateIndex i = 0; i < coordinates.size(); ++i)
+ for (Grid2D::CoordinateIndex j = 0; j < coordinates.size(); ++j)
{
if (i < j)
- Grid2D::operator()(i,j) = distance(CAs[i], CAs[j]);
+ Grid2D::operator()(i,j) = distance(coordinates[i], coordinates[j]);
else
Grid2D::operator()(i,j) = 0;
}