查看: 271|回复: 3

osgearth如何实现heatmap热力图

[复制链接]
  • TA的每日心情
    开心
    2020-3-20 17:50
  • 签到天数: 1 天

    [LV.1]初来乍到

    发表于 2019-11-22 19:57:04 | 显示全部楼层 |阅读模式
    数据:
    var points = [
            {"lng":116.418261, "lat" : 39.921984, "count" : 50},
            { "lng":116.423332,"lat" : 39.916532,"count" : 51 },
            { "lng":116.419787,"lat" : 39.930658,"count" : 15 },
            { "lng":116.418455,"lat" : 39.920921,"count" : 40 },
            { "lng":116.418843,"lat" : 39.915516,"count" : 100 },
            { "lng":116.42546,"lat" : 39.918503,"count" : 6 },
            { "lng":116.423289,"lat" : 39.919989,"count" : 18 },
            { "lng":116.418162,"lat" : 39.915051,"count" : 80 },
            { "lng":116.422039,"lat" : 39.91782,"count" : 11 },
            { "lng":116.41387,"lat" : 39.917253,"count" : 7 },
            { "lng":116.41773,"lat" : 39.919426,"count" : 42 },
            { "lng":116.421107,"lat" : 39.916445,"count" : 4 },
            { "lng":116.417521,"lat" : 39.917943,"count" : 27 },
            { "lng":116.419812,"lat" : 39.920836,"count" : 23 },
            { "lng":116.420682,"lat" : 39.91463,"count" : 60 },
            { "lng":116.415424,"lat" : 39.924675,"count" : 8 },
            { "lng":116.419242,"lat" : 39.914509,"count" : 15 },
            { "lng":116.422766,"lat" : 39.921408,"count" : 25 },
            { "lng":116.421674,"lat" : 39.924396,"count" : 21 },
            { "lng":116.427268,"lat" : 39.92267,"count" : 1 },
            { "lng":116.417721,"lat" : 39.920034,"count" : 51 },
            { "lng":116.412456,"lat" : 39.92667,"count" : 7 },
            { "lng":116.420432,"lat" : 39.919114,"count" : 11 },
            { "lng":116.425013,"lat" : 39.921611,"count" : 35 },
            { "lng":116.418733,"lat" : 39.931037,"count" : 22 },
            { "lng":116.419336,"lat" : 39.931134,"count" : 4 },
            { "lng":116.413557,"lat" : 39.923254,"count" : 5 },
            { "lng":116.418367,"lat" : 39.92943,"count" : 3 },
            { "lng":116.424312,"lat" : 39.919621,"count" : 100 },
            { "lng":116.423874,"lat" : 39.919447,"count" : 87 },
            { "lng":116.424225,"lat" : 39.923091,"count" : 32 },
            { "lng":116.417801,"lat" : 39.921854,"count" : 44 },
            { "lng":116.417129,"lat" : 39.928227,"count" : 21 },
            { "lng":116.426426,"lat" : 39.922286,"count" : 80 },
            { "lng":116.421597,"lat" : 39.91948,"count" : 32 },
            { "lng":116.423895,"lat" : 39.920787,"count" : 26 },
            { "lng":116.423563,"lat" : 39.921197,"count" : 17 },
            { "lng":116.417982,"lat" : 39.922547,"count" : 17 },
            { "lng":116.426126,"lat" : 39.921938,"count" : 25 },
            { "lng":116.42326,"lat" : 39.915782,"count" : 100 },
            { "lng":116.419239,"lat" : 39.916759,"count" : 39 },
            { "lng":116.417185,"lat" : 39.929123,"count" : 11 },
            { "lng":116.417237,"lat" : 39.927518,"count" : 9 },
            { "lng":116.417784,"lat" : 39.915754,"count" : 47 },
            { "lng":116.420193,"lat" : 39.917061,"count" : 52 },
            { "lng":116.422735,"lat" : 39.915619,"count" : 100 },
            { "lng":116.418495,"lat" : 39.915958,"count" : 46 },
            { "lng":116.416292,"lat" : 39.931166,"count" : 9 },
            { "lng":116.419916,"lat" : 39.924055,"count" : 8 },
            { "lng":116.42189,"lat" : 39.921308,"count" : 11 },
            { "lng":116.413765,"lat" : 39.929376,"count" : 3 },
            { "lng":116.418232,"lat" : 39.920348,"count" : 50 },
            { "lng":116.417554,"lat" : 39.930511,"count" : 15 },
            { "lng":116.418568,"lat" : 39.918161,"count" : 23 },
            { "lng":116.413461,"lat" : 39.926306,"count" : 3 },
            { "lng":116.42232,"lat" : 39.92161,"count" : 13 },
            { "lng":116.4174,"lat" : 39.928616,"count" : 6 },
            { "lng":116.424679,"lat" : 39.915499,"count" : 21 },
            { "lng":116.42171,"lat" : 39.915738,"count" : 29 },
            { "lng":116.417836,"lat" : 39.916998,"count" : 99 },
            { "lng":116.420755,"lat" : 39.928001,"count" : 10 },
            { "lng":116.414077,"lat" : 39.930655,"count" : 14 },
            { "lng":116.426092,"lat" : 39.922995,"count" : 16 },
            { "lng":116.41535,"lat" : 39.931054,"count" : 15 },
            { "lng":116.413022,"lat" : 39.921895,"count" : 13 },
            { "lng":116.415551,"lat" : 39.913373,"count" : 17 },
            { "lng":116.421191,"lat" : 39.926572,"count" : 1 },
            { "lng":116.419612,"lat" : 39.917119,"count" : 9 },
            { "lng":116.418237,"lat" : 39.921337,"count" : 54 },
            { "lng":116.423776,"lat" : 39.921919,"count" : 26 },
            { "lng":116.417694,"lat" : 39.92536,"count" : 17 },
            { "lng":116.415377,"lat" : 39.914137,"count" : 19 },
            { "lng":116.417434,"lat" : 39.914394,"count" : 43 },
            { "lng":116.42588,"lat" : 39.922622,"count" : 27 },
            { "lng":116.418345,"lat" : 39.919467,"count" : 8 },
            { "lng":116.426883,"lat" : 39.917171,"count" : 3 },
            { "lng":116.423877,"lat" : 39.916659,"count" : 34 },
            { "lng":116.415712,"lat" : 39.915613,"count" : 14 },
            { "lng":116.419869,"lat" : 39.931416,"count" : 12 },
            { "lng":116.416956,"lat" : 39.925377,"count" : 11 },
            { "lng":116.42066,"lat" : 39.925017,"count" : 38 },
            { "lng":116.416244,"lat" : 39.920215,"count" : 91 },
            { "lng":116.41929,"lat" : 39.915908,"count" : 54 },
            { "lng":116.422116,"lat" : 39.919658,"count" : 21 },
            { "lng":116.4183,"lat" : 39.925015,"count" : 15 },
            { "lng":116.421969,"lat" : 39.913527,"count" : 3 },
            { "lng":116.422936,"lat" : 39.921854,"count" : 24 },
            { "lng":116.41905,"lat" : 39.929217,"count" : 12 },
            { "lng":116.424579,"lat" : 39.914987,"count" : 57 },
            { "lng":116.42076,"lat" : 39.915251,"count" : 70 },
            { "lng":116.425867,"lat" : 39.918989,"count" : 8 }];
  • TA的每日心情
    开心
    2020-3-20 17:50
  • 签到天数: 1 天

    [LV.1]初来乍到

     楼主| 发表于 2020-3-26 16:33:34 | 显示全部楼层
    方法复杂main.cpp   osg::ref_ptr<osg::Node> dt = osgDB::readRefNodeFile("dumptruck.osg");

        // display a solid version of the dump truck
        osg::ref_ptr<osg:ositionAttitudeTransform> solidModel = new osg::PositionAttitudeTransform;
        solidModel->setPosition(osg::Vec3f(7.0f, -2.0f, 7.0f));
        solidModel->addChild(dt);

        // generate the 3D heatmap surface to display
        osg::ref_ptr<Heatmap> hm = new Heatmap(30, 30, 10, 30, 30, 1.0, 0.25);
        float data[30][30];
        for (int x=0; x < 30; ++x)
            for (int y=0; y < 30; ++y)
                data[y][x] = (double)rand() / RAND_MAX;
        hm->setData((float*)data, 10.0, 1.0, 0.25);

        // add a transparent version of the truck to the scene also
        osg::ref_ptr<osg::PositionAttitudeTransform> transparentTruck = new osg::PositionAttitudeTransform;
        transparentTruck->setPosition(osg::Vec3f(7.0f, -25.0f, 7.0f));

        // set the states of the truck so that it actually appears transparently and nicely lit.
        osg::StateSet *state = transparentTruck->getOrCreateStateSet();
        state->setMode(GL_BLEND, osg::StateAttribute::ON | osg::StateAttribute::OVERRIDE);
        state->setAttribute(new osg::BlendFunc(GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA), osg::StateAttribute::ON | osg::StateAttribute::OVERRIDE);
        osg::Material* material = new osg::Material;
        material->setAmbient(osg::Material::FRONT_AND_BACK,osg::Vec4(0.2f,0.2f,0.2f,0.3f));
        material->setDiffuse(osg::Material::FRONT_AND_BACK,osg::Vec4(0.8f,0.8f,0.8f,0.3f));
        state->setAttribute(material,osg::StateAttribute::ON | osg::StateAttribute::OVERRIDE);
        osg:ightModel *lm = new osg::LightModel();
        lm->setTwoSided(true);
        state->setAttribute(lm, osg::StateAttribute::ON | osg::StateAttribute::OVERRIDE);
        (transparentTruck.get())->addChild(dt.get());

        // place the heatmap and a transparent dump truck in the transparent geometry group
        osg::ref_ptr<osg::Group> transparentModel = new osg::Group;
        (transparentModel.get())->addChild(hm.get());
        (transparentModel.get())->addChild(transparentTruck.get());
  • TA的每日心情
    开心
    2020-3-20 17:50
  • 签到天数: 1 天

    [LV.1]初来乍到

     楼主| 发表于 2020-3-26 16:34:06 | 显示全部楼层
    #ifndef HEATMAP_H
    #define HEATMAP_H

    #include <osg/Geode>
    #include <osg/Uniform>
    #include <osg/Texture2D>
    #include <osg/Texture1D>

    class Heatmap : public osg::Geode
    {
    public:
        Heatmap(float width, float depth, float maxheight, unsigned int K, unsigned int N, float maximum, float transparency);
        ~Heatmap();

        void setData(float *buffer, float maxheight, float maximum, float transparency);

    protected:
        unsigned int m_K;
        unsigned int m_N;
        float *m_data;
        osg::ref_ptr<osg::Image> m_img2;
        osg::ref_ptr<osg::Texture2D> m_tex2;

        osg::ref_ptr<osg::Image> colorimg;
        osg::ref_ptr<osg::Texture1D> colortex;

        osg::Uniform *maximumUniform;
        osg::Uniform *maxheightUniform;
        osg::Uniform *transparencyUniform;
    };

    #endif // #ifndef HEATMAP_H



    /*
    * 3D Heat map using vertex displacement mapping
    * Rendered using depth peeling in fragment shader.
    */

    #include <osg/Geode>
    #include <osg/Geometry>
    #include <osg/Vec3>
    #include <osg/MatrixTransform>
    #include <osg/PositionAttitudeTransform>
    #include <osg/LightModel>
    #include <osg/io_utils>
    #include <osg/Material>
    #include <osgDB/ReadFile>
    #include <osgDB/WriteFile>
    #include <osgGA/TrackballManipulator>
    #include <osgViewer/Viewer>
    #include <osg/Math>
    #include <iostream>

    #define _USE_MATH_DEFINES
    #include <math.h>

    #include <osg/TexEnv>
    #include <osg/TexMat>
    #include <osg/Depth>
    #include <osg/ShapeDrawable>
    #include <osg/Texture1D>
    #include <osg/Texture2D>
    #include <osg/PolygonMode>
    #include <osg/PolygonOffset>

    #include "HeatMap.h"
    #include "DepthPeeling.h"

    ///////////////////////////////////////////////////////////////////////////
    // in-line GLSL source code

    static const char *VertexShader = {
        "#version 120\n"
        "uniform float maximum;\n"
        "uniform float maxheight;\n"
        "uniform float transparency;\n"
        "uniform sampler1D colortex;\n"
        "uniform sampler2D datatex;\n"
        "in vec2  xypos;\n"
        "void main(void)\n"
        "{\n"
        "    float foo;\n"
        "    float tmp     = min(texture2D(datatex, xypos).x / maximum, 1.0);\n"
        "    gl_Position   = gl_ModelViewProjectionMatrix * (gl_Vertex + vec4(0.0, 0.0, maxheight * tmp, 0.0));\n"
        "    vec4 color    = texture1D(colortex, tmp);\n"
        "    color.w       = color.w * transparency;\n"
        "    gl_FrontColor = color;\n"
        "    gl_BackColor  = color;\n"
        "}\n"
    };

    static const char *FragmentShader =
    {
        "#version 120\n"
        "bool depthpeeling();\n"
        "void main(void)\n"
        "{\n"
        "  if( depthpeeling() ) discard;\n"
        "  gl_FragColor = gl_Color;\n"
        "}\n"
    };

    /**
    * Overloaded Geometry class to return predefined bounds
    */
    class MyGeometry : public osg::Geometry
    {
    public:
        MyGeometry(osg::BoundingBox bounds)
        {
            m_bounds = bounds;
            m_bsphere = osg::BoundingSphere(bounds);
        }

        // an attempt to return a reasonable bounding box. Still does not prevent clipping of the heat map.
        virtual const osg::BoundingBox& getBoundingBox() const {return m_bounds;}
        virtual osg::BoundingBox computeBoundingBox() const {return m_bounds;}
        virtual const osg::BoundingSphere& getBound() const {return m_bsphere;}

    protected:
        osg::BoundingBox m_bounds;
        osg::BoundingSphere m_bsphere;
    };

    Heatmap::Heatmap(float width, float depth, float maxheight, unsigned int K, unsigned int N, float maximum, float transparency)
    {
        m_K = K;
        m_N = N;
        const int O = 4;

        // create Geometry object to store all the vertices primitives.
        osg::Geometry *meshGeom  = new MyGeometry(osg::BoundingBox(osg::Vec3(-width/2, -depth/2, 0), osg::Vec3(width/2, depth/2, maxheight)));

        // we use a float attribute array storing texcoords
        osg::Vec2Array* xypositions = new osg::Vec2Array();
        xypositions->setName("xypos");

        // create vertex coordinates
        osg::Vec3Array* vertices = new osg::Vec3Array();
        osg::Vec3 off(-width/2, -depth/2, 0);
        for (unsigned int y=0; y < O*N; y++) {
            if (y % 2 == 0)
            {
                for (unsigned int x=0; x < O*K; x++) {
                    vertices->push_back(osg::Vec3(width*x/(O*K-1), depth*y/(O*N-1), 0.0)+off);
                    xypositions->push_back(osg::Vec2(((float)x+0.5f)/(O*K),((float)y+0.5f)/(O*N)));
                }
            }
            else
            {
                vertices->push_back(osg::Vec3(0, depth*y/(O*N-1), 0.0)+off);
                xypositions->push_back(osg::Vec2(0.5f/(O*K),((float)y+0.5f)/(O*N)));
                for (unsigned int x=0; x < O*K-1; x++) {
                    vertices->push_back(osg::Vec3(width*(0.5+x)/(O*K-1), depth*y/(O*N-1), 0.0)+off);
                    xypositions->push_back(osg::Vec2(((float)(x+0.5f)+0.5f)/(O*K),((float)y+0.5f)/(O*N)));
                }
                vertices->push_back(osg::Vec3(width, depth*y/(O*N-1), 0.0)+off);
                xypositions->push_back(osg::Vec2(1.0f-0.5f/(O*K),((float)y+0.5f)/(O*N)));
            }
        }
        xypositions->setBinding(osg::Array::BIND_PER_VERTEX);
        xypositions->setNormalize(false);

        meshGeom->setVertexAttribArray(6, xypositions);
        meshGeom->setVertexArray(vertices);

        // generate several tri strips to form a mesh
        GLuint *indices = new GLuint[4*O*K];
        for (unsigned int y=0; y < O*N-1; y++) {
            if (y % 2 == 0)
            {
                int base  = (y/2) * (O*K+O*K+1);
                int base2 = (y/2) * (O*K+O*K+1) + O*K;
                int i=0; for (unsigned int x=0; x < O*K; x++) { indices[i++] = base2+x; indices[i++] = base+x;}
                indices[i++] = base2+O*K;
                meshGeom->addPrimitiveSet(new osg:rawElementsUInt(osg:rimitiveSet::TRIANGLE_STRIP, i, indices));
            }
            else
            {
                int base = (y/2) * (O*K+O*K+1) + O*K;
                int base2 = (y/2) * (O*K+O*K+1) + O*K + O*K+1;
                int i=0; for (unsigned int x=0; x < O*K; x++) { indices[i++] = base+x; indices[i++] = base2+x;}
                indices[i++] = base+O*K;
                meshGeom->addPrimitiveSet(new osg::DrawElementsUInt(osg::PrimitiveSet::TRIANGLE_STRIP, i, indices));
            }
        }
        delete[] indices;

        // create vertex and fragment shader
        osg::Program* program = new osg::Program;
        program->setName( "mesh" );
        program->addBindAttribLocation("xypos", 6);
        program->addShader( new osg::Shader( osg::Shader::VERTEX, VertexShader ) );
        program->addShader( new osg::Shader( osg::Shader::FRAGMENT, DepthPeeling::PeelingShader ) );
        program->addShader( new osg::Shader( osg::Shader::FRAGMENT, FragmentShader ) );

        // create a 1D texture for color lookups
        colorimg = new osg::Image();
        colorimg->allocateImage(5, 1, 1, GL_RGBA, GL_UNSIGNED_BYTE);
        unsigned char *data = colorimg->data();
        *data++ =   0; *data++ =   0; *data++ = 255; *data++ =   0;  // fully transparent blue
        *data++ =   0; *data++ = 255; *data++ = 255; *data++ = 255;  // turquoise
        *data++ =   0; *data++ = 255; *data++ =   0; *data++ = 255;  // green
        *data++ = 255; *data++ = 255; *data++ =   0; *data++ = 255;  // yellow
        *data++ = 255; *data++ =   0; *data++ =   0; *data++ = 255;  // red
        colortex = new osg::Texture1D(colorimg.get());
        colortex->setFilter(osg::Texture::MIN_FILTER, osg::Texture:INEAR);
        colortex->setFilter(osg::Texture::MAG_FILTER, osg::Texture::LINEAR);
        colortex->setWrap(osg::Texture::WRAP_S, osg::Texture::CLAMP_TO_EDGE);
        colortex->setResizeNonPowerOfTwoHint(false);

        // create a 2D texture for data lookups
        m_img2 = new osg::Image();
        m_img2->allocateImage(K, N, 1, GL_LUMINANCE, GL_FLOAT);
        m_img2->setInternalTextureFormat(GL_RGB32F_ARB);
        m_data = (float*)m_img2.get()->data();
        m_tex2 = new osg::Texture2D(m_img2.get());
        m_tex2.get()->setResizeNonPowerOfTwoHint(false);
        m_tex2.get()->setFilter(osg::Texture::MIN_FILTER, osg::Texture::LINEAR);
        m_tex2.get()->setFilter(osg::Texture::MAG_FILTER, osg::Texture::LINEAR);
        m_tex2.get()->setWrap(osg::Texture::WRAP_S, osg::Texture::CLAMP_TO_EDGE);
        m_tex2.get()->setWrap(osg::Texture::WRAP_T, osg::Texture::CLAMP_TO_EDGE);

        // set render states
        osg::StateSet *meshstate = meshGeom->getOrCreateStateSet();
        meshstate->setMode(GL_BLEND,  osg::StateAttribute::ON);
        meshstate->setMode(GL_LIGHTING, osg::StateAttribute::OFF);
        meshstate->setAttributeAndModes(program, osg::StateAttribute::ON);
        meshstate->setTextureAttributeAndModes(0,colortex.get(),osg::StateAttribute::ON);
        meshstate->setTextureAttributeAndModes(1,m_tex2.get(),osg::StateAttribute::ON);

        // uniforms for height and color scaling
        maximumUniform = new osg::Uniform( "maximum", (float)maximum );
        maxheightUniform = new osg::Uniform( "maxheight", (float)maxheight );
        transparencyUniform = new osg::Uniform( "transparency", (float)transparency);

        osg::Uniform* texUniform = new osg::Uniform(osg::Uniform::SAMPLER_1D, "colortex");
        texUniform->set(0);
        osg::Uniform* texUniform2 = new osg::Uniform(osg::Uniform::SAMPLER_2D, "datatex");
        texUniform2->set(1);
        meshstate->addUniform(texUniform);
        meshstate->addUniform(texUniform2);
        meshstate->addUniform(maximumUniform);
        meshstate->addUniform(maxheightUniform);
        meshstate->addUniform(transparencyUniform);

        // add the geometries to the geode.
        meshGeom->setUseDisplayList(false);
        addDrawable(meshGeom);
    }

    void Heatmap::setData(float *buffer, float maxheight, float maximum, float transparency)
    {
        memcpy(m_data, buffer, m_N*m_K*sizeof(float));

        maximumUniform->set( maximum );
        maxheightUniform->set( maxheight );
        transparencyUniform->set ( transparency );

        m_img2.get()->dirty();
    }

    Heatmap::~Heatmap()
    {
    }
  • TA的每日心情
    开心
    2020-3-20 17:50
  • 签到天数: 1 天

    [LV.1]初来乍到

     楼主| 发表于 2020-3-26 16:50:01 | 显示全部楼层
    本帖最后由 liyihongcug 于 2020-5-25 16:01 编辑

    还是有较大问题, osgearth实现应该在球体一个矩形平面绘图 。应该输入数据,还是有差别.
    反馈思考还是搞定,借助插件。
    解决这个的思路:如果从原生来做,画圆调色等比较慢速度----  参考Baidu地图确实Js都不见得快(heatmap.js)
    反复考虑 采用空间数据库来做比较合适i
    空间数据库可以更新地理数据和数值大小样式自己定义。
    Osgearth本身是支持OGC的,这样只需找打heatmap的插件就可以。完美解决。
    唯一动态设置的需要对样式的间隔进行手动设置
    您需要登录后才可以回帖 登录 | 注册

    本版积分规则

    OSG中国官方论坛-有您OSG在中国才更好

    网站简介:osgChina是国内首个三维相关技术开源社区,旨在为国内更多的技术开发人员提供最前沿的技术资讯,为更多的三维从业者提供一个学习、交流的技术平台。

    联系我们

    • 工作时间:09:00--18:00
    • 反馈邮箱:1315785073@qq.com
    快速回复 返回顶部 返回列表