@property
    def areal(self):
        """Compute areal distortion of the flatmap.

        Areal distortion is calculated at each triangle as the log2 ratio of
        the triangle area in the flatmap to the area in the reference surface.
        Distortion values are then resampled onto the vertices.

        Thus a value of 0 indicates the areas are equal (no distortion), a 
        value of +1 indicates that the area in the flatmap is 2x the area
        in the reference surface (expansion), and a value of -1 indicates
        that the area in the flatmap is 1/2x the area in the reference
        surface (compression).

        See: http://brainvis.wustl.edu/wiki/index.php/Caret:Operations/Morphing

        Returns
        -------
        vertratios : 1D ndarray, shape (total_verts,)
            Areal distortion at each vertex.
        """
        def area(pts, polys):
            ppts = pts[polys]
            cross = np.cross(ppts[:,1] - ppts[:,0], ppts[:,2] - ppts[:,0])
            return np.sqrt((cross**2).sum(-1))

        refarea = area(self.ref, self.polys)
        flatarea = area(self.flat, self.polys)
        tridists = np.log2(flatarea/refarea)
        
        vertratios = np.zeros((len(self.ref),))
        vertratios[self.polys[:,0]] += tridists
        vertratios[self.polys[:,1]] += tridists
        vertratios[self.polys[:,2]] += tridists
        vertratios /= np.bincount(self.polys.ravel())
        vertratios = np.nan_to_num(vertratios)
        vertratios[vertratios==0] = 1
        return vertratios
dataset/ETHPy150Open gallantlab/pycortex/cortex/polyutils.py/Distortion.areal/original