# Interpolation of non-rectangular data failed

I want to interpolate a type of data which is on a triangular lattice in order to make DensityPlot faster(ListDensityPlot is so slow). However, the interoplation failed with error, no matter setting InterpolationOrder->1 or InterpolationOrder->All

Here is a less dense data

data = Import["https://pastebin.com/raw/6XmFDzmf"];
plotData =
Partition[
StringCases[data, x : NumberString :> InternalStringToDouble@x],
3];
Interpolation[plotData]


it will raise several errors

Interpolation::udeg: Interpolation on unstructured grids is currently only supported for InterpolationOrder->1 or InterpolationOrder->All. Order will be reduced to 1.

Interpolation::femimq: The element mesh has insufficient quality of 0.. A quality estimate below 0. may be caused by a wrong ordering of element incidents or self-intersecting elements.

Interpolation::fememtlq: The quality 0.of the underlying mesh is too low. The quality needs to be larger than 0..

What does it mean? How to interpolate such data or general non-rectangular data?

PS: the density plot is

The problem is that for regular (but not rectangular) meshes the Delaunay mesh is unstable. It's a bug, also mentioned in this question, which was about rectangular grids. The same workaround works here -- just jiggle the points a tiny bit around their perfect lattice positions:

epsilon = 10^-7;
jiggledPlotData = {#1 + RandomReal[epsilon {-1, 1}], #2 +
RandomReal[epsilon {-1, 1}], #3} & @@@ plotData;
reg = ConvexHullMesh[jiggledPlotData[[All, 1 ;; 2]]];
f = Interpolation[jiggledPlotData, InterpolationOrder -> 1];
DensityPlot[f[x, y], {x, y} \[Element] reg]

• Thank you so much! You saved my day! Dec 16, 2017 at 6:00

A better method than jiggle original grid which I learned recently is dealing with mesh explicitly.

First, we need a package

Needs["NDSolveFEM"];


We can generate Delaunay mesh by

originalMesh = ToElementMesh[plotData[[;; , 1 ;; 2]]]


but this will throw a warning message

ToElementMesh::femimq: The element mesh has insufficient quality of 0.. A quality estimate below 0. may be caused by a wrong ordering of element incidents or self-intersecting elements.

this is the same femimq problem, though the mesh looks fine at first sight as(enlarged)

The ElementMesh provides powerful tools to examing what is going on.

originalMesh["Quality"][[1]] // Histogram


The histogram shows most of the triangle has good quality close to 1. But there are plenty of triangle has quality near 0 which are bad enough.

using highlightBadMeshElement(code at the end), we can show the location of these bad triangles.

highlightBadMeshElement[originalMesh, 0.2]


They are all located at the edge shown by red marker.

We can delete these redundent bad triangles and form a refined mesh using refineMesh(code at the end)

mesh = refineMesh[originalMesh, 0.2];
f = ElementMeshInterpolation[{mesh}, plotData[[;; , -1]]];
DensityPlot[f[x, y], {x, y} \[Element] mesh]


gives

functions:

highlightBadMeshElement[mesh_, qualityThreshold_] := Module[{},
pos = Position[mesh["Quality"], _?(# <= qualityThreshold &)];
Show[
mesh["Wireframe"[
"ElementMeshDirective" ->
Directive[EdgeForm[GrayLevel[.6]], FaceForm[]]]],
mesh["Wireframe"[pos, "MeshElement" -> "MeshElements",
"ElementMeshDirective" -> Directive[EdgeForm[Red], FaceForm[]]]]
, Boxed -> False]]

refineMesh[mesh_, qualityThreshold_] := Module[{},
qualityList = First@mesh["Quality"];
pos = Flatten@Position[qualityList, _?(# > qualityThreshold &)];
If[Length@mesh["MeshElements"] != 1,
Message[refineMesh::MoreThanOneTypeOfElement,
"More than one type of element. 1";
`