AnalysisData data_;
AnalysisDataAverageModulePointer avem_;
-
-
+
+
Selection selec;
std::vector< std::vector < RVec > > trajectory;
std::vector< double > noise;
std::vector< int > index;
- double epsi = 0.1; // should be selectable
- int etalon_frame = 0; // should be selectable
- int frames = 0;
+ double epsi = 0.1; // should be selectable
+ int etalon_frame = 0; // should be selectable
+ int frames = 0;
};
.store(&fnDist_).defaultBasename("avedist")
.description("Average distances from reference group"));
- options->addOption(SelectionOption("reference")
+ options->addOption(SelectionOption("reference")
.store(&refsel_).required()
.description("Reference group to calculate distances from"));*/
options->addOption(SelectionOption("select")
.store(&selec).required()
.description("Atoms that are considered as part of the excluded volume"));
- /* options->addOption(DoubleOption("cutoff").store(&cutoff_)
- .description("Cutoff for distance calculation (0 = no cutoff)"));
+ /* options->addOption(DoubleOption("cutoff").store(&cutoff_)
+ .description("Cutoff for distance calculation (0 = no cutoff)"));
- settings->setFlag(TrajectoryAnalysisSettings::efRequireTop);*/
+ settings->setFlag(TrajectoryAnalysisSettings::efRequireTop);*/
}
void
{
index.resize(0);
ConstArrayRef< int > atomind = selec.atomIndices();
- for (ConstArrayRef< int >::iterator ai = atomind.begin(); (ai < atomind.end()); ai++) {
+ for (ConstArrayRef< int >::iterator ai = atomind.begin(); (ai < atomind.end()); ai++)
+ {
index.push_back(*ai);
}
}
{
trajectory.resize(frames + 1);
trajectory[frames].resize(index.size());
- for (int i = 0; i < index.size(); i++) {
+ for (int i = 0; i < index.size(); i++)
+ {
trajectory[frames][i] = fr.x[index[i]];
}
frames++;
void
AnalysisTemplate::finishAnalysis(int /*nframes*/)
{
- real *w_rls;
- rvec *x1, *x2, temp;
+ real *w_rls;
+ rvec *x1, *x2, temp;
std::vector< bool > flag;
- double noise_mid = 9999;
- int count = 0;
+ double noise_mid = 9999;
+ int count = 0;
snew(x1, index.size());
snew(x2, index.size());
flag.resize(index.size(), true);
noise.resize(index.size(), 0);
- while (noise_mid > epsi) {
+ while (noise_mid > epsi)
+ {
noise_mid = 0;
- count = 0;
+ count = 0;
noise.assign(index.size(), 0);
- for (int i = 0; i < index.size(); i++) {
- if (flag[i]) {
+ for (int i = 0; i < index.size(); i++)
+ {
+ if (flag[i])
+ {
w_rls[i] = 1;
- } else {
+ }
+ else
+ {
w_rls[i] = 0;
}
}
- for (int i = 0; i < frames; i++) {
- for (int j = 0; j < trajectory[i].size(); j++) {
+ for (int i = 0; i < frames; i++)
+ {
+ for (int j = 0; j < trajectory[i].size(); j++)
+ {
copy_rvec(trajectory[i][j].as_vec(), x2[j]);
copy_rvec(trajectory[etalon_frame][j].as_vec(), x1[j]); // i think its nessesary to do because of lots of reset_x
}
reset_x(index.size(), NULL, index.size(), NULL, x1, w_rls);
reset_x(index.size(), NULL, index.size(), NULL, x2, w_rls);
do_fit(index.size(), w_rls, x1, x2);
- if (i % 100 == 0) {
+ if (i % 100 == 0)
+ {
std::cout << i << " / " << frames << "\n";
}
- for (int j = 0; j < index.size(); j++) {
+ for (int j = 0; j < index.size(); j++)
+ {
rvec_sub(x1[j], x2[j], temp);
noise[j] += norm(temp);
}
}
- for (int i = 0; i < index.size(); i++) {
+ for (int i = 0; i < index.size(); i++)
+ {
noise[i] /= frames;
}
- for (int i = 0; i < index.size(); i++) {
- if (flag[i]) {
+ for (int i = 0; i < index.size(); i++)
+ {
+ if (flag[i])
+ {
count++;
noise_mid += noise[i];
}
}
noise_mid /= count;
- for (int i = 0; i < index.size(); i++) {
- if (noise[i] > noise_mid) {
+ for (int i = 0; i < index.size(); i++)
+ {
+ if (noise[i] > noise_mid)
+ {
flag[i] = false;
}
}
}
- for (int i = 0; i < noise.size(); i++) {
- if (noise[i] <= epsi) {
+ for (int i = 0; i < noise.size(); i++)
+ {
+ if (noise[i] <= epsi)
+ {
//fprintf(pFile, "%d ", i + 1);
// numbers from index, +1 for starting from 1 not 0
}