VCRA-SZCZERBACKI-PointBasedRendering/splatting.cpp

#include "splatting.h"

#include "analytic.inl"

extern Mesh    *lwo;
extern Color   *buffer;
extern ABuffer *abuffer;

extern std::vector<Splat>       splats;

void loadLwoFile(char *filename) 
{
  glutSetCursor(GLUT_CURSOR_WAIT);

  lwo = new Mesh( filename );
        lwo->setParametrization();
        lwo->subdivide( subdivision );
        splats.resize(lwo->getVertexCount());
        for(int i=0;i<splats.size();++i)
        {
                float maxl = 0.0f,l;
                int valence, left, right;
                bool boundary;
                const std::vector<Mesh::ventry> &vertexlist = lwo->get1Neighbourhood( i, valence, boundary, left, right );
                for(int j=0;j<vertexlist.size();++j)
                {
                        l = (lwo->getVertex(vertexlist[j].m_vertex).m_position - lwo->getVertex(i).m_position).length();
                        maxl = max2(l,maxl);
                }

                splats[i].m_position = lwo->getVertex(i).m_position;
                splats[i].m_normal = lwo->getVertex(i).m_normal;
                splats[i].m_size = maxl;
    splats[i].m_color.set(1.0f,0.85f,0.65f, 1.0f);
        }

  glutSetCursor(GLUT_CURSOR_INHERIT);
}

void poissonDiscDistribution( Vector2 *points, int count )
{
        SS_ASSERT(count > 0);

        float threshold = 1.0f;

        //repeat until we have all subpixels distributed
        int j=0;
        while(j < count)
        {
                //try to place a random sample in a unit circle so that it is not closer
                // to earlier samples than the threshold 
                int tries = 100;
                int k;
                Vector2 p;
                do
                {
                        //make a random sample inside unit circle
                        p = randInCircle();

                        //test whether the sample is closer to points already placed than the threshold
                        for(k=0;k<j;k++)
                        {
                                if( (p - points[k]).length() < threshold )
                                {       //too close, break
                                        tries--;
                                        break;
                                }
                        }
                } while(k < j && tries > 0);

                if( tries == 0 )
                {
                        //if did not succeed, lower the threshold and redistribute all points
                        threshold *= 0.98f;     //the closer the number is to 1, the slower the convergence
                        j = 0;
                }
                else
                {
                        //if succeeded, add the sample
                        points[j] = p;
                        j++;
                }
        }
}

float cumulativeGaussian[gausssamples];

void initCumulativeGaussian()
{
        //init cumulative gaussian
        float x = 0.0f;
        float xstep = prefiltercutoff / (float)gausssamples;
        x += xstep/2.0f;        //calculate x at the center of each slice
        cumulativeGaussian[0] = exp(-0.5f*x*x)*xstep;
        for(int i=1;i<gausssamples;i++)
        {
                x += xstep;
                cumulativeGaussian[i] = cumulativeGaussian[i-1] + exp(-0.5f*x*x)*xstep;
        }
        //normalize
        for(int i=0;i<gausssamples;i++)
        {
                cumulativeGaussian[i] /= cumulativeGaussian[gausssamples-1];
        }
}

float inverseCumulativeGaussian( float r )
{
        int i;
        for(i=0;i<gausssamples-1;i++)
        {
                if( r > cumulativeGaussian[i] && r < cumulativeGaussian[i+1] )
                        break;
        }
        return (float)i/(float)gausssamples;
}

void makeJitter(int subpixels)
{
        SS_ASSERT(subpixels > 0);

        if( prefilterJitter ) delete[] prefilterJitter;
        if( lensJitter) delete[] lensJitter;

        prefilterJitter = new Vector2[subpixels*jitterPatterns];
        lensJitter = new Vector2[subpixels*jitterPatterns];

        svPrint("poisson disk distribution with %d points.",subpixels);
        Timer t;
        t.sample();
        initCumulativeGaussian();
        for(int p=0;p<jitterPatterns;p++)
        {
                poissonDiscDistribution(&prefilterJitter[p*subpixels],subpixels);
                for(int j=0;j<subpixels;j++)
                {
                        Vector2 jit = prefilterJitter[p*subpixels+j];
                        float r = jit.length();
                        if( r > 0.0f ) jit *= 1.0f/r;
                        jit *= inverseCumulativeGaussian(r) * prefiltersize;
                        prefilterJitter[p*subpixels+j] = jit;
                }

                poissonDiscDistribution(&lensJitter[p*subpixels],subpixels);
                for(int j=0;j<subpixels;j++)
                {
                        lensJitter[p*subpixels+j] *= lensRadius;
                }
        }
        t.sample();
        svPrint(" took %f seconds\n",t.getDelta());
}

---