Commit edc492ee authored by Andrius Merkys's avatar Andrius Merkys

New upstream version 0.8.9.2

parent eeaed758
2018-04-01 Paul Emsley <pemsley at mrc dash lmb dot ac dot uk>
* Release 0.8.9.2
2018-03-18 Paul Emsley <pemsley at mrc dash lmb dot ac dot uk>
* Release 0.8.9.1
......
-----
Release 0.8.9.2
o FEATURE: LINK info added to header browser
o FEATURE: External angle restraints are now available (corresponding to REFMAC type
0 restraints) - mini-rsr included [Robbie Joosten]
o FEATURE: raster3d output now includes labels [Vito Calderone]
o FEATURE: API for user-defined unimodal torsion restraints (for pyranoses) [Juan Hermoso]
o FEATURE: "Variable" bond widths i.e. as you zoom out, the bond lines get thinner
[Erec Stebbins]
(set-use-variable-bond-thickness 1)
o CHANGE: Iron-Sulfur clusters are now bonded using the dictionary (and so appear "cubic")
[Hannah Bridges]
o CHANGE: Edit -> Copy Fragment now has a "Move it Here" check-button
o CHANGE: Old interface to PRODRG removed [Elinor Breiner]
o CHANGE: Jiggle-fit speed up by using backbone mode if possible
o CHANGE: Forwards-fitting a terminal residue now correctly positions the O atom
of the current residue that corresponds to the new residue
o CHANGE: Planar peptide restraints now apply to PTRANS restraints as well as TRANS
restraints
o CHANGE: "Linear" refine no longer exists, the refinement button now uses the same
interface as sphere refine. Residues now have environments but "finding peptide
links by residue number, no matter how far the residues are away from each
other" has been removed
o CHANGE: Squares to cirlces in the Ramachandran Plot
o CHANGE: Default is now not to generate a monomer on reading a dictionary cif file
o CHANGE: More functions added to the coot_extended python module
o CHANGE: RCrane restored
o CHANGE: Add link-by-torsion-to-NAG-core-NAG-SER.tab so that add_linked_residue()
can link a NAG to a SER
o CHANGE: CA Zone -> Mainchain builds in both directions
o CHANGE: Alignment gap penalty changed to -3.0
o CHANGE: Merging a ligand molecule now adds the ligand to the closest chain - unless
a merging molecule ligand has been specified
o BUG-FIX: the H atom in a peptide is now moved along with the other peptide atoms
on a cis->trans conversion (and vice versa)
o BUG-FIX: align-and-mutate occasional reversing of residues fixed [Rob Kirchdoerfer]
O BUG-FIX: Using PDB-REDO now loads an anomalous map also (if possible) [Robbie Joosten]
o BUG-FIX: Fix typo in template-keybindings-to-preferences that causes initial
installation of key-bindings to fail sometimes [Colin Palmer].
o BUG-FIX: The slow update on deleting deviant extra distance restraints has been
been replaced with a faster version [JED]
o BUG-FIX: db_mainchain() fragment/chain indexing bug fixed [Oliver Clarke]
o BUG-FIX: db_mainchain() no longer crashes when using negative residue numbers
[Oliver Clarke]
o BUG-FIX: Delete All Carbohydrate now works with carbohydrate with insertion codes
-----
Release 0.8.9.1
......
......@@ -13,7 +13,8 @@ AM_CPPFLAGS = \
$(MMDB_CXXFLAGS)
libcoot_analysis_la_SOURCES = bfkurt.cc mogul.cc kolmogorov.hh kolmogorov.cc stats.cc \
cablam.hh cablam.cc typed-distances.cc typed-distances.hh
cablam.hh cablam.cc typed-distances.cc typed-distances.hh \
b-factor-histogram.cc b-factor-histogram.hh
libcoot_analysis_la_LIBADD = \
$(top_builddir)/coot-utils/libcoot-coord-utils.la \
......@@ -27,7 +28,7 @@ libcoot_analysis_la_LDFLAGS = $(SHARED_LDFLAGS)
bin_PROGRAMS = coot-bfactan
check_PROGRAMS = test-mogul test-cablam peptide-projection test-align
check_PROGRAMS = test-mogul test-cablam peptide-projection test-align improper-dihedrals
coot_bfactan_SOURCES = bfactan.cc
......@@ -46,6 +47,11 @@ peptide_projection_SOURCES = peptide-projection.cc
peptide_projection_LDADD = libcoot-analysis.la
improper_dihedrals_SOURCES = improper-dihedrals.cc
improper_dihedrals_LDADD = $(top_builddir)/coords/libcoot-coords.la \
libcoot-analysis.la
test_align_LDADD = libcoot-analysis.la
coot_bfactan_LDADD = \
......
......@@ -37,7 +37,8 @@ build_triplet = @build@
host_triplet = @host@
bin_PROGRAMS = coot-bfactan$(EXEEXT)
check_PROGRAMS = test-mogul$(EXEEXT) test-cablam$(EXEEXT) \
peptide-projection$(EXEEXT) test-align$(EXEEXT)
peptide-projection$(EXEEXT) test-align$(EXEEXT) \
improper-dihedrals$(EXEEXT)
subdir = analysis
DIST_COMMON = $(pkginclude_HEADERS) $(srcdir)/Makefile.am \
$(srcdir)/Makefile.in
......@@ -113,7 +114,7 @@ libcoot_analysis_la_DEPENDENCIES = \
$(top_builddir)/utils/libcoot-utils.la $(am__DEPENDENCIES_1) \
$(am__DEPENDENCIES_1)
am_libcoot_analysis_la_OBJECTS = bfkurt.lo mogul.lo kolmogorov.lo \
stats.lo cablam.lo typed-distances.lo
stats.lo cablam.lo typed-distances.lo b-factor-histogram.lo
libcoot_analysis_la_OBJECTS = $(am_libcoot_analysis_la_OBJECTS)
libcoot_analysis_la_LINK = $(LIBTOOL) --tag=CXX $(AM_LIBTOOLFLAGS) \
$(LIBTOOLFLAGS) --mode=link $(CXXLD) $(AM_CXXFLAGS) \
......@@ -127,6 +128,10 @@ coot_bfactan_DEPENDENCIES = libcoot-analysis.la \
$(top_builddir)/mini-mol/libcoot-mini-mol.la \
$(top_builddir)/utils/libcoot-utils.la $(am__DEPENDENCIES_1) \
$(am__DEPENDENCIES_1) $(am__DEPENDENCIES_1)
am_improper_dihedrals_OBJECTS = improper-dihedrals.$(OBJEXT)
improper_dihedrals_OBJECTS = $(am_improper_dihedrals_OBJECTS)
improper_dihedrals_DEPENDENCIES = \
$(top_builddir)/coords/libcoot-coords.la libcoot-analysis.la
am_peptide_projection_OBJECTS = peptide-projection.$(OBJEXT)
peptide_projection_OBJECTS = $(am_peptide_projection_OBJECTS)
peptide_projection_DEPENDENCIES = libcoot-analysis.la
......@@ -163,11 +168,13 @@ LINK = $(LIBTOOL) --tag=CC $(AM_LIBTOOLFLAGS) $(LIBTOOLFLAGS) \
--mode=link $(CCLD) $(AM_CFLAGS) $(CFLAGS) $(AM_LDFLAGS) \
$(LDFLAGS) -o $@
SOURCES = $(libcoot_analysis_la_SOURCES) $(coot_bfactan_SOURCES) \
$(peptide_projection_SOURCES) $(test_align_SOURCES) \
$(test_cablam_SOURCES) $(test_mogul_SOURCES)
$(improper_dihedrals_SOURCES) $(peptide_projection_SOURCES) \
$(test_align_SOURCES) $(test_cablam_SOURCES) \
$(test_mogul_SOURCES)
DIST_SOURCES = $(libcoot_analysis_la_SOURCES) $(coot_bfactan_SOURCES) \
$(peptide_projection_SOURCES) $(test_align_SOURCES) \
$(test_cablam_SOURCES) $(test_mogul_SOURCES)
$(improper_dihedrals_SOURCES) $(peptide_projection_SOURCES) \
$(test_align_SOURCES) $(test_cablam_SOURCES) \
$(test_mogul_SOURCES)
HEADERS = $(pkginclude_HEADERS)
ETAGS = etags
CTAGS = ctags
......@@ -389,7 +396,8 @@ AM_CPPFLAGS = \
$(MMDB_CXXFLAGS)
libcoot_analysis_la_SOURCES = bfkurt.cc mogul.cc kolmogorov.hh kolmogorov.cc stats.cc \
cablam.hh cablam.cc typed-distances.cc typed-distances.hh
cablam.hh cablam.cc typed-distances.cc typed-distances.hh \
b-factor-histogram.cc b-factor-histogram.hh
libcoot_analysis_la_LIBADD = \
$(top_builddir)/coot-utils/libcoot-coord-utils.la \
......@@ -410,6 +418,10 @@ test_cablam_LDADD = $(top_builddir)/coords/libcoot-coords.la \
peptide_projection_SOURCES = peptide-projection.cc
peptide_projection_LDADD = libcoot-analysis.la
improper_dihedrals_SOURCES = improper-dihedrals.cc
improper_dihedrals_LDADD = $(top_builddir)/coords/libcoot-coords.la \
libcoot-analysis.la
test_align_LDADD = libcoot-analysis.la
coot_bfactan_LDADD = \
libcoot-analysis.la \
......@@ -544,6 +556,9 @@ clean-checkPROGRAMS:
coot-bfactan$(EXEEXT): $(coot_bfactan_OBJECTS) $(coot_bfactan_DEPENDENCIES)
@rm -f coot-bfactan$(EXEEXT)
$(CXXLINK) $(coot_bfactan_OBJECTS) $(coot_bfactan_LDADD) $(LIBS)
improper-dihedrals$(EXEEXT): $(improper_dihedrals_OBJECTS) $(improper_dihedrals_DEPENDENCIES)
@rm -f improper-dihedrals$(EXEEXT)
$(CXXLINK) $(improper_dihedrals_OBJECTS) $(improper_dihedrals_LDADD) $(LIBS)
peptide-projection$(EXEEXT): $(peptide_projection_OBJECTS) $(peptide_projection_DEPENDENCIES)
@rm -f peptide-projection$(EXEEXT)
$(CXXLINK) $(peptide_projection_OBJECTS) $(peptide_projection_LDADD) $(LIBS)
......@@ -563,9 +578,11 @@ mostlyclean-compile:
distclean-compile:
-rm -f *.tab.c
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/b-factor-histogram.Plo@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/bfactan.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/bfkurt.Plo@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/cablam.Plo@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/improper-dihedrals.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/kolmogorov.Plo@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/mogul.Plo@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/peptide-projection.Po@am__quote@
......
#include <iostream>
#include <cmath>
#include "b-factor-histogram.hh"
// set n_atoms, n_bins, b_max
coot::b_factor_histogram::b_factor_histogram(mmdb::Manager *mol) {
init();
b_max = -1.0;
n_atoms = 0;
for(int imod = 1; imod<=mol->GetNumberOfModels(); imod++) {
mmdb::Model *model_p = mol->GetModel(imod);
if (model_p) {
int n_chains = model_p->GetNumberOfChains();
for (int ichain=0; ichain<n_chains; ichain++) {
mmdb::Chain *chain_p = model_p->GetChain(ichain);
int nres = chain_p->GetNumberOfResidues();
for (int ires=0; ires<nres; ires++) {
mmdb::Residue *residue_p = chain_p->GetResidue(ires);
int n_atoms_in_residue = residue_p->GetNumberOfAtoms();
for (int iat=0; iat<n_atoms_in_residue; iat++) {
mmdb::Atom *at = residue_p->GetAtom(iat);
const float &b = at->tempFactor;
if (b >= 0.0) {
n_atoms++;
if (b > b_max) {
b_max = b;
}
}
}
}
}
}
}
if (n_atoms > 0) {
n_bins = get_n_bins(); // use b_max and n_atoms
}
b_vector.resize(n_bins);
for(int imod = 1; imod<=mol->GetNumberOfModels(); imod++) {
mmdb::Model *model_p = mol->GetModel(imod);
if (model_p) {
int n_chains = model_p->GetNumberOfChains();
for (int ichain=0; ichain<n_chains; ichain++) {
mmdb::Chain *chain_p = model_p->GetChain(ichain);
int nres = chain_p->GetNumberOfResidues();
for (int ires=0; ires<nres; ires++) {
mmdb::Residue *residue_p = chain_p->GetResidue(ires);
int n_atoms_in_residue = residue_p->GetNumberOfAtoms();
for (int iat=0; iat<n_atoms_in_residue; iat++) {
mmdb::Atom *at = residue_p->GetAtom(iat);
const float &b = at->tempFactor;
if (b >= 0.0) {
int bin_idx = b_to_bin(b);
b_vector[bin_idx].push_back(b);
}
}
}
}
}
}
}
coot::b_factor_histogram::b_factor_histogram(mmdb::Manager *mol, int atom_selection_handle) {
init();
b_max = -1.0;
n_atoms = 0;
mmdb::Atom **atom_selection = 0;
int n_selection_atoms;
mol->GetSelIndex(atom_selection_handle, atom_selection, n_selection_atoms);
for (int i=0; i<n_selection_atoms; i++) {
mmdb::Atom *at = atom_selection[i];
const float &b = at->tempFactor;
if (b >= 0.0) {
n_atoms++;
if (b > b_max) {
b_max = b;
}
}
}
if (n_atoms > 0) {
n_bins = get_n_bins(); // use b_max and n_atoms
}
b_vector.resize(n_bins);
for (int i=0; i<n_selection_atoms; i++) {
mmdb::Atom *at = atom_selection[i];
const float &b = at->tempFactor;
if (b >= 0.0) {
int bin_idx = b_to_bin(b);
b_vector[bin_idx].push_back(b);
}
}
}
int
coot::b_factor_histogram::b_to_bin(const float &b) const {
int ibin = 0;
double f = b/static_cast<double>(b_max);
ibin = static_cast<int>(f*n_bins);
// std::cout << "b_to_bin() ibin " << ibin << " for f " << f << " b_max " << b_max
// << std::endl;
if (ibin >= n_bins)
ibin = n_bins - 1;
return ibin;
}
int
coot::b_factor_histogram::get_n_bins() const {
float f_n_bins = 40.0f;
int nb = static_cast<int>(f_n_bins);
return nb;
}
std::vector<std::pair<double, double> >
coot::b_factor_histogram::get_data() const {
std::vector<std::pair<double, double> > v(n_bins);
for (std::size_t i=0; i<b_vector.size(); i++) {
int freq = b_vector[i].size();
double frac = static_cast<double>(i)/static_cast<double>(n_bins);
double b_val = b_max * frac;
std::pair<double, double> p(b_val, freq);
v.push_back(p);
}
return v;
}
std::vector<std::pair<double, double> >
coot::b_factor_histogram::get_model() const {
int n_model_points = 100; // smooth enough?
std::vector<std::pair<double, double> > v(n_model_points+1); // so that we capture the biggest b.
if (false) {
std::cout << "in get_model() with alpha_estimate " << alpha_estimate << std::endl;
std::cout << "in get_model() with beta_estimate " << beta_estimate << std::endl;
}
double sf = 1.0;
double sum = 0.0;
for (int i=0; i<=n_model_points; i++) {
double frac = static_cast<double>(i)/static_cast<double>(n_model_points);
double b = b_max * frac;
double b_model = ig(b + shift_estimate);
sum += b_model;
}
sf = 1.0/sum;
double sf_2 = 0.5 * static_cast<double>(b_max) / static_cast<double>(get_n_bins());
for (int i=0; i<=n_model_points; i++) {
double frac = static_cast<double>(i)/static_cast<double>(n_model_points);
double b = b_max * frac;
double b_model = ig(b + shift_estimate);
std::pair<double, double> p(b, b_model * sf * sf_2 * n_atoms);
// std::cout << " " << b << " " << b_model << std::endl;
v.push_back(p);
}
return v;
}
double
coot::b_factor_histogram::ig(const double &x) const {
if (x <= 0.0) {
return 0.0;
} else {
double g_a = Gamma(alpha_estimate);
double part_1 = (std::pow(beta_estimate, alpha_estimate) * std::pow(x, -alpha_estimate-1.0))/g_a;
double part_2 = std::exp(-beta_estimate/x);
if (false)
std::cout << " ig: for x " << x << " g_a: " << g_a << " " << part_1 << " " << part_2 << " returning "
<< part_1 * part_2 << std::endl;
return part_1 * part_2;
}
}
double
coot::b_factor_histogram::Gamma(const double &b) const {
#ifdef HAVE_CXX11
return std::tgamma(b);
#else
std::cout << "ERROR:: Missing CXX11" << std::endl;
return 1.0;
#endif
}
void
coot::b_factor_histogram::init() {
alpha_estimate = 0.0;
beta_estimate = 0.0;
shift_estimate = 0.0;
}
void
coot::b_factor_histogram::model() {
// This is pure invese gamma modelling - not shifted inverse gamma.
// model this bf distribution with parameters alpha, beta
// bf is shifted from b-factor:
// bf = b-factor + shift
// f(bf) = beta^alpha/Gamma(alpha) * (1/)^(alpha+1) exp(-beta/bf)
// where Gamma(alpha) is
// Gamma(alpha) = (n-1)!
// \mathcal{IG}(x|\alpha,\beta) = \frac{\beta^\alpha x^{-\alpha-1}}{|\Gamma{\alpha}}exp(\frac{-\beta}{x})
// Gaussian approximation to the moments, mean (mu) and variance (nu)
// \mu = \frac{\beta}{\alpha-1}
// \nu = \frac{\beta^2}{(\alpha-1)^2(\alpha-2))
// \alpha(hat) = \frac{\mu^2}{\nu} + 2
// \beta(hat) = \mu(\frac{\mu^2}{\nu}+1)
double sum = 0.0;
double sum_sq = 0.0;
int n = 0;
for (std::size_t ibin=0; ibin<b_vector.size(); ibin++) {
for (std::size_t i=0; i<b_vector[ibin].size(); i++) {
const float &b = b_vector[ibin][i];
sum += b;
sum_sq += b*b;
n++;
}
}
double n_d = static_cast<double>(n);
double mu = sum/n_d;
double nu = sum_sq/n_d - mu * mu;
if (nu < 0.0) nu = 0;
alpha_estimate = mu * mu / nu + 2.0;
beta_estimate = mu * (mu*mu/nu+1);
// optimize_estimates(); too complicated to make work at the moment.
}
#include <fstream>
#include "kolmogorov.hh"
#include "utils/coot-utils.hh"
std::vector<double>
coot::b_factor_histogram::select_from_model() const {
unsigned int n_points = 200;
std::vector<double> v;
for (;;) {
double b_phi = b_max * static_cast<double>(util::random())/static_cast<double>(RAND_MAX);
double pr = ig(b_phi + shift_estimate);
double r = static_cast<double>(util::random())/static_cast<double>(RAND_MAX);
bool state = false;
if (pr > r) state = true;
// std::cout << " select_from_model " << b_phi << " " << pr << " " << state << std::endl;
if (pr > r) {
v.push_back(b_phi);
if (v.size() == n_points)
break;
}
}
std::ofstream f("bfm.tab");
for (std::size_t i=0; i<v.size(); i++) {
f << i << " " << v[i] << "\n";
}
f.close();
return v;
}
void
coot::b_factor_histogram::optimize_estimates() {
// this doesn't work
double alpha_orig = alpha_estimate;
std::vector<double> data_1;
std::vector<double> data_2;
for (std::size_t i=0; i<b_vector.size(); i++) {
const std::vector<float> &bv = b_vector[i];
for (std::size_t j=0; j<bv.size(); j++) {
data_1.push_back(bv[j]);
}
}
std::cout << "alpha_orig " << alpha_orig << std::endl;
double l = 0.99;
for (std::size_t i=0; i<20; i++) {
double f = static_cast<double>(i) * 0.05;
double m_factor = (1.0-l) + 2.0 * l * f;
alpha_estimate = m_factor * alpha_orig;
std::vector<std::pair<double, double> > m = get_model();
// how do I turn a probability model into a set of data?
// data_2 = some function of m
data_2 = select_from_model();
std::pair<double, double> kl_divergences = nicholls::get_KL(data_1, data_2);
std::cout << "f " << f << " l " << l << " alpha " << alpha_estimate << " k-l div: "
<< kl_divergences.first << " " << kl_divergences.second
<< std::endl;
}
// restore alpha_estimate
alpha_estimate = alpha_orig;
}
#include <vector>
#include <mmdb2/mmdb_manager.h>
namespace coot {
class b_factor_histogram {
int n_bins;
int n_atoms;
float b_max;
void init();
int get_n_atoms() const;
int get_n_bins() const;
int b_to_bin(const float &b) const;
std::vector<std::vector<float> > b_vector;
float alpha_estimate;
float beta_estimate;
float shift_estimate;
double ig(const double &x) const;
double Gamma(const double &b) const;
void optimize_estimates();
std::vector<double> select_from_model() const;
public:
b_factor_histogram(mmdb::Manager *mol);
b_factor_histogram(mmdb::Manager *mol, int atom_selection_handle);
void model();
std::vector<std::pair<double, double> > get_data() const;
std::vector<std::pair<double, double> > get_model() const; // call model() first
};
}
#include <string>
#include "coords/mmdb-extras.h"
#include "coords/mmdb.h"
struct atom_vector_map_t {
coot::atom_quad q1;
coot::atom_quad q2;
int res_no;
std::string chain_id;
std::string fn;
};
void print_tors(std::map<std::string, std::vector<atom_vector_map_t> > &atom_vector_map) {
std::map<std::string, std::vector<atom_vector_map_t> >::const_iterator it;
for(it=atom_vector_map.begin(); it!=atom_vector_map.end(); it++) {
const std::string &key = it->first;
const std::vector<atom_vector_map_t> &vv = it->second;
for (std::size_t i=0; i<vv.size(); i++) {
const atom_vector_map_t &qp = vv[i];
double tors_1 = qp.q1.torsion();
double tors_2 = qp.q2.torsion();
std::cout << "pseudo-torsion " << " " << qp.fn << " " << qp.chain_id << " " << qp.res_no << " "
<< key << " N " << tors_1 << " C " << tors_2 << "\n";