[Shapelib] Calculating areas

jakob at lanstorp.com jakob at lanstorp.com
Thu Aug 26 09:31:32 EDT 2004


I am trying to calculate an area with coordinates extracted from a shape file with shapelib. It does work fine except for polygons with holes in.

I can not id the inner ring, and always get the ring (just ring not outer or inner) returned for all segments of my shape. If I can id the inner ring I can substract that area from the total. I enclose a method in Object Pascal for calculating areas of a shape:

function TReadShape.CalcArea(ft: integer): double;
// This function calculate the areal of a single shape object
// ft : current shape
// Return value is area in m^2
var
  psShape               : PSHPObject;
  nSegment              : integer;
  nVertices             : integer;
  area                  : double;
  areaSegment           : double;
  i,j                   : integer;
  nVertexInSegment      : integer;
  nVertexInLastSegment  : integer;
  index                 : integer;
  VertexInSegment       : array of TVertexInSegment;
begin
  psShape   := SHPReadObject(hSHP, ft);    // Reference to current shape feature
  nSegment  := psShape^.nParts;           // Number of segment
  nVertices := psShape^.nVertices;

  area := 0;

  // If polygon has only one segment
  if nSegment = 1 then begin
    for j := 0 to nVertices -2 do begin             //*** Loop every vertex in current segment
      area := area + ((psShape^.padfX[j+1] * psShape^.padfY[j]) - (psShape^.padfX[j] * psShape^.padfY[j+1])) / 2;
    end;
  end;

  // If polygon has two or more segments
  if nSegment > 1 then begin
    //*** First get the number of vertexes in each segment
    SetLength(VertexInSegment,nSegment);
    for i := 0 to nSegment-1 do begin
      if i = nSegment-1 then VertexInSegment[i].count := nVertices - psShape^.panPartStart[i]   //*** Count vertexes in last segment
      else VertexInSegment[i].count := psShape^.panPartStart[i+1] - psShape^.panPartStart[i];   //*** count number of vertexes in every segment
    end;

    area                  := 0;
    areaSegment           := 0;
    nVertexInLastSegment  := 0;

    // loop segments to calc every segment to add up to the total area
    for i := 0 to nSegment-1 do begin                               //*** Loop every segment
      nVertexInSegment := VertexInSegment[i].count;
      //*** Try to get type of segment, i.e. inner ring or outer ring
      ShowMessage('SHPPartTypeName(psShape.panPartType(i)): ' + SHPPartTypeName(psShape^.panPartType[i]));

      for j := 0 to nVertexInSegment -2 do begin                    //*** Loop every vertex in current segment
        index := nVertexInLastSegment + j;                          //*** Get index of next vertex
        areaSegment := areaSegment + ((psShape^.padfX[index+1] * psShape^.padfY[index]) - (psShape^.padfX[index] * psShape^.padfY[index+1])) / 2;
      end;
      area := area + areaSegment;                    // sum up all segments
      areaSegment := 0;
      nVertexInLastSegment := nVertexInSegment;      // This is needed to jump from on segment to the next when finding vertex number      areaSegment := 0;
    end;
  end;

  Result := Abs(area);      // Get absolute value. The area can be negativ if nodes arranged in clockwise direction
end;

Thx,
Jake



More information about the Shapelib mailing list