PathMap.h

Go to the documentation of this file.

 
#ifndef IMPBFF_PATHMAP_H
#define IMPBFF_PATHMAP_H
 
#include <IMP/bff/bff_config.h>
 
#include <stdlib.h>     /* malloc, free, rand */
#include <limits>
 
#include <cmath>
#include <algorithm>
#include <unordered_set>
#include <queue>
#include <vector>
#include <utility>  /* std::pair */
#include <Eigen/Dense>
 
#include <IMP/Object.h>
#include <IMP/Particle.h>
#include <IMP/core/XYZR.h>
#include <IMP/em/SampledDensityMap.h>
 
#include <IMP/em/MRCReaderWriter.h>
#include <IMP/em/XplorReaderWriter.h>
#include <IMP/em/EMReaderWriter.h>
#include <IMP/em/SpiderReaderWriter.h>
 
#include <IMP/bff/PathMapHeader.h>
#include <IMP/bff/PathMapTile.h>
#include <IMP/bff/PathMapTileEdge.h>
 
IMPBFF_BEGIN_NAMESPACE
 
class PathMapTile;
 
 
class IMPBFFEXPORT PathMap : public IMP::em::SampledDensityMap {
 
friend class PathMapTile;
friend class AV;
 
private:
 
    // used in path search
    std::vector<bool>  visited;
    std::vector<bool>  edge_computed;
    std::vector<float> cost;
 
protected:
 
    std::vector<PathMapTile> tiles;
    PathMapHeader pathMapHeader_;
    std::vector<int> offsets_;
    std::vector<PathMapTileEdge>& get_edges(int tile_idx);
 
public:
 
 
    void update_tiles(
        float obstacle_threshold=-1.0, 
        bool binarize=true, 
        float obstacle_penalty=TILE_PENALTY_DEFAULT,
        bool reset_tile_edges=true
    );
 
    void resize(unsigned int nvox);
 
    void set_data(double *input, int n_input, 
        float obstacle_threshold=-1, bool binarize=true, 
        float obstacle_penalty=TILE_PENALTY_DEFAULT);
 
    std::vector<int> get_neighbor_idx_offsets(double neighbor_radius = -1){
        if(neighbor_radius < 0){
            auto pmh = get_path_map_header();
            neighbor_radius = pmh->get_neighbor_radius();
        }
        const int nn = ceil(neighbor_radius);
        const double nr2 = neighbor_radius * neighbor_radius;
 
        const IMP::em::DensityHeader* header = get_header();
        int nx = header->get_nx();
        int ny = header->get_ny();
        int nx_ny = nx * ny;
 
        std::vector<int> offsets;
        for(int z = -nn; z < nn; z += 1) {
            double dz2 = z * z;
            int oz = z * nx_ny;
            for(int y = -nn; y < nn; y++) {
                int oy = y * nx;
                double dz2_dy2 = dz2 + y * y;
                for(int x = -nn; x < nn; x++) {
                    int ox = x;
                    int dz2_dy2_dx2 = dz2_dy2 + x * x;
                    if(dz2_dy2_dx2 <= nr2){
                        int d;                  // edge_cost is a float stored in an 32bit int
                        float *p = (float*) &d; // Make a float pointer point at the integer
                        *p = sqrt((float) dz2_dy2_dx2); // Pretend that the integer is a float and store the value
                        int tile_offset = oz + oy + ox;
                        offsets.emplace_back(z);
                        offsets.emplace_back(y);
                        offsets.emplace_back(x);
                        offsets.emplace_back(tile_offset);
                        offsets.emplace_back(d);
                    }
                }
            }
        }
        return offsets;
    }
 
 
    int get_dim_index_by_voxel(long index, int dim);
 
    const PathMapHeader *get_path_map_header() const { return &pathMapHeader_; }
 
 
    PathMapHeader *get_path_map_header_writable() { return &pathMapHeader_; }
 
 
    void set_path_map_header(PathMapHeader &path_map_header, float resolution = -1.0);
 
 
    std::vector<float> get_tile_values(
            int value_type = PM_TILE_COST,
            std::pair<float, float> bounds = std::pair<float, float>(
                    {std::numeric_limits<float>::min(),
                     std::numeric_limits<float>::max()}),
            const std::string &feature_name=""
    );
 
    void get_tile_values(
            float **output, int *nx, int *ny, int *nz,
            int value_type = PM_TILE_COST,
            std::pair<float, float> bounds = std::pair<float, float>(
                    {std::numeric_limits<float>::min(),
                     std::numeric_limits<float>::max()}),
            const std::string &feature_name=""
    );
 
    std::vector<PathMapTile>& get_tiles();
 
 
    void fill_sphere(IMP::algebra::Vector3D r0, double radius, double value, bool inverse=true);
 
    void find_path(long path_begin_idx, long path_end_idx = -1, int heuristic_mode = 0);
 
 
    void find_path_dijkstra(long path_begin_idx, long path_end_idx = -1);
    
    
    void find_path_astar(long path_begin_idx, long path_end_idx = -1);
 
    std::vector<IMP::algebra::Vector4D> get_xyz_density();
    
    void get_xyz_density(double** output, int* n_output1, int* n_output2);
 
    void sample_obstacles(double extra_radius=0.0);
 
    explicit PathMap(
            PathMapHeader &header,
            std::string name = "PathMap%1%",
            IMP::em::KernelType kt = IMP::em::BINARIZED_SPHERE,
            float resolution = -1.0
    );
 
};
 
 
IMPEMEXPORT
void write_path_map(
        PathMap *m,
        std::string filename,
        int value_type,
        const std::pair<float, float> bounds = std::pair<float, float>(
                std::numeric_limits<float>::min(),
                std::numeric_limits<float>::max()
        ),
        const std::string &feature_name = ""
);
 
IMPBFF_END_NAMESPACE
 
#endif //IMPBFF_PATHMAP_H