Added computation of PDBDistanceGrid for all atoms of the backbone (not just alpha-carbons)
authorDmitriy Morozov <morozov@cs.duke.edu>
Wed, 24 Jan 2007 17:50:45 -0500
changeset 14 d68c28a18f37
parent 13 b1535347e37d
child 15 ef3bf83302ca
Added computation of PDBDistanceGrid for all atoms of the backbone (not just alpha-carbons)
examples/grid/pdbdistance-vineyard.cpp
examples/grid/pdbdistance.h
--- 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;
 				}