All the code here is under the CC-0 Universal Licence which means that you are free to do whatever you want with it (i.e. it is Public Domain code)
This page contains the following examples, written in all the languages supported by chemfiles:
indexes
program read the first frame of a trajectory and fill a
vector called indexes
with the indexes of all the atom obeying the
\(x < 5\) relation.rmsd
program read all the frames from a trajectory, and compute the
root mean square displacement of the 10th atom in the frame.convert
program read all the frames from a trajectory, set the unit
cell and the topology of the frame, and then write it back to another file.
The input file contain here the trajectory of a single water molecule, but
without any topological or cell information. So the program add these
informations to a new trajectory file.#include <iostream>
#include <vector>
#include "chemfiles.hpp"
int main() {
chemfiles::Trajectory traj("filename.xyz");
chemfiles::Frame frame;
traj >> frame;
auto positions = frame.positions();
std::vector<size_t> indexes;
for (size_t i=0; i<frame.natoms(); i++) {
if (positions[i][0] < 5) {
indexes.push_back(i);
}
}
std::cout << "Atoms with x < 5: " << std::endl;
for (auto i : indexes) {
std::cout << " - " << i << std::endl;
}
return 0;
}
#include <iostream>
#include <vector>
#include <numeric>
#include <cmath>
#include "chemfiles.hpp"
int main() {
chemfiles::Trajectory traj("filename.nc");
chemfiles::Frame frame;
std::vector<double> distances;
// Accumulate the distances to the origin of the 10th atom throughtout the
// trajectory
while (!traj.done()) {
traj >> frame;
// Position of the 10th atom
auto pos = frame.positions()[9];
double distance = sqrt(pos[0]*pos[0] + pos[1]*pos[1] + pos[2]*pos[2]);
distances.push_back(distance);
}
double mean = std::accumulate(std::begin(distances), std::end(distances), 0.0) / distances.size();
double rmsd = 0.0;
for (auto dist : distances) {
rmsd += (mean - dist)*(mean - dist);
}
rmsd /= distances.size();
rmsd = sqrt(rmsd);
std::cout << "Root-mean square displacement is: " << rmsd << std::endl;
return 0;
}
#include "chemfiles.hpp"
int main() {
chemfiles::Trajectory input("water.xyz");
chemfiles::Frame frame{};
chemfiles::Topology water_topology{};
// Orthorombic UnitCell with lengths of 20, 15 and 35 A
chemfiles::UnitCell cell(20, 15, 35);
// Create Atoms
chemfiles::Atom O("O");
chemfiles::Atom H("H");
// Fill the topology with one water molecule
water_topology.append(O);
water_topology.append(H);
water_topology.append(H);
water_topology.add_bond(0, 1);
water_topology.add_bond(0, 2);
chemfiles::Trajectory output("water.pdb", "w");
while (!input.done()) {
input >> frame;
// Set the frame cell and topology
frame.cell(cell);
frame.topology(water_topology);
// Write the frame to the output file, using PDB format
output << frame;
}
return 0;
}
#include <stdio.h>
#include <stdlib.h>
#include "chemfiles.h"
int main() {
CHFL_TRAJECTORY* traj = chfl_trajectory_open("tests/files/xyz/helium.xyz", "r");
CHFL_FRAME* frame = chfl_frame(0);
if (traj == NULL)
goto error;
if (chfl_trajectory_read(traj, frame))
goto error;
size_t natoms = 0;
chfl_frame_atoms_count(frame, &natoms);
float (*positions)[3] = (float(*)[3])malloc(natoms*3*sizeof(float));
unsigned* indexes = (unsigned*)malloc(natoms*sizeof(unsigned));
if (positions == NULL || indexes == NULL)
goto error;
for (int i=0; i<natoms; i++) {
indexes[i] = (unsigned)-1;
}
if (chfl_frame_positions(frame, positions, natoms))
goto error;
unsigned last_index = 0;
for (int i=0; i<natoms; i++) {
if (positions[i][0] < 5) {
indexes[last_index] = i;
last_index++;
}
}
printf("Atoms with x < 5:\n");
int i = 0;
while(indexes[i] != (unsigned)-1 && i < natoms) {
printf(" - %d\n", indexes[i]);
i++;
}
printf("Number of atoms: %d\n", i);
chfl_trajectory_close(traj);
chfl_frame_free(frame);
free(positions);
return 0;
error:
printf("Error, cleaning up …\n");
chfl_trajectory_close(traj);
chfl_frame_free(frame);
free(positions);
return 1;
}
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include "chemfiles.h"
int main() {
CHFL_TRAJECTORY* traj = chfl_trajectory_open("filename.nc", "r");
CHFL_FRAME* frame = chfl_frame(0);
float (*positions)[3] = NULL;
double* distances = NULL;
if(traj == NULL || frame == NULL)
goto error;
size_t nsteps = 0;
chfl_trajectory_nsteps(traj, &nsteps);
distances = (double*)malloc(sizeof(double)*nsteps);
if (distances == NULL)
goto error;
// Accumulate the distances to the origin of the 10th atom throughtout the
// trajectory
for (size_t i=0; i<nsteps; i++) {
if(!chfl_trajectory_read(traj, frame))
goto error;
// Only allocate on the first iteration. This assume a constant number
// of particles
size_t natoms = 0;
if (i == 0) {
chfl_frame_atoms_count(frame, &natoms);
positions = (float(*)[3])malloc(natoms*3*sizeof(float));
if (positions == NULL)
goto error;
}
// Position of the 10th atom
chfl_frame_positions(frame, positions, natoms);
double distance = sqrt(positions[9][0]*positions[9][0] +
positions[9][1]*positions[9][1] +
positions[9][2]*positions[9][2]);
distances[i] = distance;
}
double mean = 0;
for (size_t i=0; i<nsteps; i++) {
mean += distances[i];
}
mean /= nsteps;
double rmsd = 0.0;
for (size_t i=0; i<nsteps; i++) {
rmsd += (mean - distances[i])*(mean - distances[i]);
}
rmsd /= nsteps;
rmsd = sqrt(rmsd);
printf("Root-mean square displacement is: %f", rmsd);
// Free the memory
chfl_trajectory_close(traj);
chfl_frame_free(frame);
free(distances);
free(positions);
return 0;
error:
printf("Error, cleaning up …\n");
chfl_trajectory_close(traj);
chfl_frame_free(frame);
free(distances);
free(positions);
return 1;
}
#include <stdio.h>
#include <stdlib.h>
#include "chemfiles.h"
int main() {
CHFL_TRAJECTORY* input = chfl_trajectory_open("water.xyz", "r");
CHFL_TRAJECTORY* output = chfl_trajectory_open("water.pdb", "w");
CHFL_FRAME* frame = chfl_frame(0);
CHFL_TOPOLOGY* water_topology = chfl_topology();
// Orthorombic UnitCell with lengths of 20, 15 and 35 A
CHFL_CELL* cell = chfl_cell(20, 15, 35);
// Create Atoms
CHFL_ATOM* O = chfl_atom("O");
CHFL_ATOM* H = chfl_atom("H");
if (input == NULL || frame == NULL || water_topology == NULL ||
cell == NULL || O == NULL || H == NULL || output == NULL)
goto error;
// Fill the topology with one water molecule
chfl_topology_append(water_topology, O);
chfl_topology_append(water_topology, H);
chfl_topology_append(water_topology, H);
chfl_topology_add_bond(water_topology, 0, 1);
chfl_topology_add_bond(water_topology, 0, 2);
size_t nsteps;
chfl_trajectory_nsteps(input, &nsteps);
for (size_t i=0; i<nsteps; i++) {
chfl_trajectory_read(input, frame);
// Set the frame cell and topology
chfl_frame_set_cell(frame, cell);
chfl_frame_set_topology(frame, water_topology);
// Write the frame to the output file, using PDB format
chfl_trajectory_write(output, frame);
}
chfl_trajectory_close(input);
chfl_trajectory_close(output);
chfl_frame_free(frame);
chfl_cell_free(cell);
chfl_topology_free(water_topology);
chfl_atom_free(O);
chfl_atom_free(H);
return 0;
error:
printf("Error, cleaning up …\n");
chfl_trajectory_close(input);
chfl_trajectory_close(output);
chfl_frame_free(frame);
chfl_cell_free(cell);
chfl_topology_free(water_topology);
chfl_atom_free(O);
chfl_atom_free(H);
return 1;
}