AV.h

Go to the documentation of this file.

 
#ifndef IMPBFF_AV_H
#define IMPBFF_AV_H
 
#include <IMP/bff/bff_config.h>
 
#include <IMP/decorator_macros.h>
#include <IMP/Decorator.h>
#include <IMP/algebra/Vector3D.h>
#include <IMP/algebra/Transformation3D.h>
#include <IMP/core/XYZ.h>
#include <IMP/log.h>
 
 
#include <IMP/bff/PathMap.h>
 
#include <string>
#include <cmath>
#include <vector>
#include <limits>
#include <iostream>            // std::cout, std::cout, std::flush
 
#include <IMP/bff/internal/json.h>
#include <IMP/bff/internal/InverseSampler.h>
// requires C++14
// #include <boost/histogram.hpp> // make_histogram, regular, weight, indexed
#include <IMP/bff/internal/Histogram.h>
 
IMPBFF_BEGIN_NAMESPACE
 
 
typedef enum{
    DYE_PAIR_DISTANCE_E,            
    DYE_PAIR_DISTANCE_MEAN,         
    DYE_PAIR_DISTANCE_MP,           
    DYE_PAIR_EFFICIENCY,            
    DYE_PAIR_DISTANCE_DISTRIBUTION, 
    DYE_PAIR_XYZ_DISTANCE           
} DyePairMeasures;
 
 
 
class IMPBFFEXPORT AVPairDistanceMeasurement{
 
public:
 
    double distance = -1;
    double error_neg = -1;
    double error_pos = -1;
    double forster_radius = 52.0;
    int distance_type = IMP::bff::DYE_PAIR_DISTANCE_MEAN;
    std::string position_1;
    std::string position_2;
 
    std::string get_json();
 
    double score_model(double model);
 
};
 
 
 
 
class IMPBFFEXPORT AV : public IMP::core::Gaussian {
 
private:
 
    IMP::bff::PathMap* av_map_ = nullptr;
 
protected:
 
    IMP::algebra::VectorD<9> get_parameter() const {
        IMP::algebra::VectorD<9> v;
        v[0] = get_linker_length();
        v[1] = get_radius1();
        v[2] = get_radius2();
        v[3] = get_radius3();
        v[4] = get_linker_width();
        v[5] = get_allowed_sphere_radius();
        v[6] = get_contact_volume_thickness();
        v[7] = get_contact_volume_trapped_fraction();
        v[8] = get_simulation_grid_resolution();
        return v;
    }
 
public:
 
    IMP::bff::PathMapHeader create_path_map_header();
 
    void init_path_map();
 
    static FloatKey get_av_key(unsigned int i){
        IMP_USAGE_CHECK(i < 9, "Out of range av feature");
        static const FloatKey k[] = {
                FloatKey("linker_length"),
                FloatKey("radius1"),
                FloatKey("radius2"),
                FloatKey("radius3"),
                FloatKey("linker_width"),
                FloatKey("allowed_sphere_radius"),
                FloatKey("contact_volume_thickness"),
                FloatKey("contact_volume_trapped_fraction"),
                FloatKey("simulation_grid_resolution")
        };
        return k[i];
    }
 
    static void do_setup_particle(Model *m, ParticleIndex pi,
                                ParticleIndex pi_source,
                                double linker_length = 20.0,
                                const algebra::Vector3D radii = algebra::Vector3D(3.5, 0, 0),
                                double linker_width = 0.5,
                                double allowed_sphere_radius = 1.5,
                                double contact_volume_thickness = 0.0,
                                double contact_volume_trapped_fraction = -1,
                                double simulation_grid_resolution = 1.5) {
        if (!IMP::core::Gaussian::get_is_setup(m, pi)) {
            IMP::core::Gaussian::setup_particle(m, pi);
        }
        m->add_attribute(get_av_key(0), pi, linker_length);
        m->add_attribute(get_av_key(1), pi, radii[0]);
        m->add_attribute(get_av_key(2), pi, radii[1]);
        m->add_attribute(get_av_key(3), pi, radii[2]);
        m->add_attribute(get_av_key(4), pi, linker_width);
        m->add_attribute(get_av_key(5), pi, allowed_sphere_radius);
        m->add_attribute(get_av_key(6), pi, contact_volume_thickness);
        m->add_attribute(get_av_key(7), pi, contact_volume_trapped_fraction);
        m->add_attribute(get_av_key(8), pi, simulation_grid_resolution);
        m->add_attribute(get_particle_key(0), pi, pi_source);
    }
 
    static ParticleIndexKey get_particle_key(unsigned int i) {
        static const ParticleIndexKey k[1] = {
            ParticleIndexKey("Source particle")
        };
        return k[i];
    }
 
    ParticleIndex get_particle_index() const {
        return Decorator::get_particle_index();
    }
 
    ParticleIndex get_particle_index(unsigned int i) const {
        return get_model()->get_attribute(get_particle_key(i), get_particle_index());
    }
 
    Particle* get_particle() const {
        return Decorator::get_particle();
    }
 
    Particle* get_particle(unsigned int i) const {
        return get_model()->get_particle(get_particle_index(i));
    }
    IMP_DECORATOR_METHODS(AV, IMP::core::Gaussian);
 
    IMP_DECORATOR_SETUP_1(AV, IMP::ParticleIndex, pi_source);
    IMP_DECORATOR_GET_SET(linker_length, get_av_key(0), Float, Float);
    IMP_DECORATOR_GET_SET(radius1, get_av_key(1), Float, Float);
    IMP_DECORATOR_GET_SET(radius2, get_av_key(2), Float, Float);
    IMP_DECORATOR_GET_SET(radius3, get_av_key(3), Float, Float);
    IMP_DECORATOR_GET_SET(linker_width, get_av_key(4), Float, Float);
    IMP_DECORATOR_GET_SET(allowed_sphere_radius, get_av_key(5), Float, Float);
    IMP_DECORATOR_GET_SET(contact_volume_thickness, get_av_key(6), Float, Float);
    IMP_DECORATOR_GET_SET(contact_volume_trapped_fraction, get_av_key(7), Float, Float);
    IMP_DECORATOR_GET_SET(simulation_grid_resolution, get_av_key(8), Float, Float);
 
 
    IMP::algebra::Vector3D get_radii(){
        return IMP::algebra::Vector3D({get_radius1(), get_radius3(), get_radius3()});
    }
 
 
    bool get_parameters_are_optimized() const {
        return
            get_particle()->get_is_optimized(get_av_key(0)) &&
            get_particle()->get_is_optimized(get_av_key(1)) &&
            get_particle()->get_is_optimized(get_av_key(2)) &&
            get_particle()->get_is_optimized(get_av_key(3)) &&
            get_particle()->get_is_optimized(get_av_key(6)) &&
            get_particle()->get_is_optimized(get_av_key(7));
    }
 
    void set_av_parameters_are_optimized(bool tf) const {
        get_particle()->set_is_optimized(get_av_key(0), tf);
        get_particle()->set_is_optimized(get_av_key(1), tf);
        get_particle()->set_is_optimized(get_av_key(2), tf);
        get_particle()->set_is_optimized(get_av_key(3), tf);
        get_particle()->set_is_optimized(get_av_key(6), tf);
        get_particle()->set_is_optimized(get_av_key(7), tf);
    }
 
    void set_av_parameter(const nlohmann::json &j);
 
 
 
    algebra::Vector3D get_derivatives() const {
        return get_model()->get_coordinate_derivatives(get_particle_index());
    }
 
    static bool get_is_setup(Model *m, ParticleIndex pi) {
        return m->get_has_attribute(get_av_key(2), pi);
    }
 
    IMP::bff::PathMap* get_map() const;
 
    void resample(bool shift_xyz=true);
 
    IMP::algebra::Vector3D get_mean_position(bool include_source=true) const;
 
    IMP::algebra::Vector3D get_source_coordinates() const;
 
    Particle* get_source() const;
 
};
 
IMP_DECORATORS(AV, AVs, ParticlesTemp);
 
template<typename T>
T inline fret_efficiency(T distance, double forster_radius){
    double rda_r0_6 = std::pow(distance / forster_radius, 6.0);
    return 1. / (1. + rda_r0_6);
}
 
template<typename T>
T inline distance_fret(double fret_efficiency, double forster_radius){
    return forster_radius * std::pow(1. / fret_efficiency - 1.0, 1. / 6.);
}
 
IMPBFFEXPORT double av_distance(
        const AV& a,
        const AV& b,
        double forster_radius = 52.0,
        int distance_type = DYE_PAIR_DISTANCE_MEAN,
        int n_samples = 10000
);
 
// Draw random points in AV. Returns (x,y,z,d) vector
IMPBFFEXPORT std::vector<double> av_random_points(
        const AV& av1,
        int n_samples=10000
);
 
IMPBFFEXPORT std::vector<double> av_random_distances(
        const AV& av1,
        const AV& av2,
        int n_samples=10000
);
 
 
IMPBFFEXPORT std::vector<double> av_distance_distribution(
        const AV& av1,
        const AV& av2,
        std::vector<double> axis,
        //double start, double stop, int n_bins, // for boost histogram
        int n_samples=10000
);
 
 
IMPBFFEXPORT IMP::ParticleIndex search_labeling_site(
        const IMP::core::Hierarchy& hier,
        std::string json_str = "",
        const nlohmann::json &json_data = nlohmann::json()
);
 
 
IMPBFF_END_NAMESPACE
 
#endif /* IMPBFF_AV_H */