On the fly Darstellung von Grenzlängen und Richtungswinkeln

--Rainer Gareis 12:19, 26. Aug 2010 (CEST) Darstellung der Grenzlängen und Winkel in der Liegenschaftskarte Um in einem Tooltip-fenster die Grenzlänge(Gl) und den Richtungswinkel(Rw) zu sehen, muß man eine neue Tabelle mit oid und Geometrie erzeugen. In dieser Tabelle steht dann die Flurstücksgeometrie aufgeschlüsselt als Line und nicht mehr als Multipoygon.

Tabelle erzeugen:

 CREATE TABLE grenzabschnitte2( 
 the_geom geometry,
 CONSTRAINT enforce_dims_the_geom CHECK (ndims(the_geom) = 2),
 CONSTRAINT enforce_geotype_the_geom CHECK (geometrytype(the_geom) = 'LINESTRING'::text OR the_geom IS NULL),
 CONSTRAINT enforce_srid_the_geom CHECK (srid(the_geom) = 2398)
 )
 WITH (OIDS=TRUE);

um die Tabelle füllen zu können benötigt man eine neue Funktion "linefrompoly"

CREATE OR REPLACE FUNCTION linefrompoly(geometry)
 RETURNS geometry AS
$BODY$SELECT 
     geomfromtext(
		 replace(
			 replace(
				 replace(
					 replace(
						 replace(
 							 replace(
								 replace(
									 asText($1),'MULTIPOLYGON','MULTILINESTRING'
								 ),'POLYGON','LINESTRING'
							 ), '((', '('
						 ), '))', ')'
					 ), '(((', '(('
				 ), ')))', '))'
			 ), '))),(((', ')),(('
		 ), srid($1)
	 )$BODY$
    LANGUAGE 'sql' IMMUTABLE STRICT;

mit folgendem Befehl wird dann die Tabelle gefüllt.

 insert into grenzabschnitte2 SELECT makeline(pointn(foo.linestring, foo.count1), pointn(foo.linestring, foo.count1 + 1)) AS the_geom
  FROM ( SELECT generate_series(1, npoints(foofoo.linestring)) AS count1, foofoo.linestring
          FROM ( SELECT linefrompoly(alkobj_e_fla.the_geom) AS linestring
                  FROM alkobj_e_fla where folie = '001') foofoo) foo
 WHERE (foo.count1 + 1) <= npoints(foo.linestring) 

Gist-Index erzeugen, damit die räumliche Abfrage auch schnell durchgeführt wird:

CREATE INDEX grenzabschnitte_gist
 ON grenzabschnitte2
 USING gist
 (the_geom);

(Dauer ca. 15 Minuten für einen Kompletten Landkreis)

Jetzt kann man sich noch eine Sicht erzeugen, in der die Länge und der Richtingswinkel berechnet wird.

CREATE OR REPLACE VIEW grenzabschnitte AS 
SELECT grenzabschnitte2.oid, grenzabschnitte2.the_geom, round(length(grenzabschnitte2.the_geom)::numeric, 2) AS laenge, 
       CASE
           WHEN ((y(endpoint(grenzabschnitte2.the_geom)) - y(startpoint(grenzabschnitte2.the_geom))) = 0 AND (x(endpoint(grenzabschnitte2.the_geom)) > x(startpoint(grenzabschnitte2.the_geom))))THEN '100.0000 gon'::text
           WHEN ((y(endpoint(grenzabschnitte2.the_geom)) - y(startpoint(grenzabschnitte2.the_geom))) = 0 AND (x(endpoint(grenzabschnitte2.the_geom)) < x(startpoint(grenzabschnitte2.the_geom))))THEN '300.0000 gon'::text
           WHEN round((degrees(atan((x(endpoint(grenzabschnitte2.the_geom)) - x(startpoint(grenzabschnitte2.the_geom))) / (y(endpoint(grenzabschnitte2.the_geom)) - y(startpoint(grenzabschnitte2.the_geom))))) * 1.1111111::double precision)::numeric, 4)::double precision < 0::double precision THEN (round((degrees(atan((x(endpoint(grenzabschnitte2.the_geom)) - x(startpoint(grenzabschnitte2.the_geom))) / (y(endpoint(grenzabschnitte2.the_geom)) - y(startpoint(grenzabschnitte2.the_geom))))) * 1.1111111::double precision)::numeric, 4)::double precision + 400::double precision) || ' gon'::text
           ELSE round((degrees(atan((x(endpoint(grenzabschnitte2.the_geom)) - x(startpoint(grenzabschnitte2.the_geom))) / (y(endpoint(grenzabschnitte2.the_geom)) - y(startpoint(grenzabschnitte2.the_geom))))) * 1.1111111::double precision)::numeric, 4)::double precision::text || ' gon'::text
       END AS winkel
  FROM grenzabschnitte2;

Nun benötigt man nur noch einen abfragbaren Layer ohne Klassen.

Mit freundlicher Unterstützung durch STEFAN RAHN ist dieses Kunststück gelungen.

--Markus Hentschel 08:29, 5. Okt 2010 (CEST) Da zwei benachbarte Flurstücke auch zwei gemeinsame Grenzabschnitte mitbringen, muss man die doppelten noch löschen. Dafür bietet sich folgender Befehl an:

delete from grenzabschnitte2 where oid in (
select a.oid from grenzabschnitte2 a left join grenzabschnitte2 b 
on length(a.the_geom)=length(b.the_geom) 
and a.the_geom=b.the_geom 
and x(pointn(a.the_geom,1)) > x(pointn(b.the_geom,1))
)

Problem: Es gibt zu viele Linien, die Abfrage dauert unendlich. Lösung: ein Script. In diesem Script kann dann auch das Füllen der Tabelle passieren. Das Script geht kilometerquadratweise vor und löscht die doppelten Linien immer nur in diesem Bereich, wodurch die Abfragezeit auf der Datenbank deutlich verringert wird.

Dazu muss eine weitere Tabelle angelegt werden, die ebenfalls nur ein Geometriefeld enthält. Angenommen, sie heißt grenzabschnitte1. Dann sieht das Script so aus:

#!/bin/ksh

# Füllt die Tabelle grenzabschnitte2 mit den Grenzlängen aus den Flurstücksgeometrien

cd /usr/local/pgsql/bin

LOGFILE=/pfad/zu/load_grenzabschnitte.log; export LOGFILE
DBNAME=kvwmapsp; export DBNAME
DBUSER=kvwmap; export DBUSER

# Definition des zu kachelnden Gebiets im Ziel-Koordinatenreferenzsystem S42/83
min_rechts=4517000   # Rechtswert der Ecke links unten
min_hoch=5982000     # Hochwert der Ecke links unten
max_rechts=70000     # maximale Ausdehnung des Bereichs in Metern (min_rechts + max_rechts = Rechtswert rechts oben)
max_hoch=59000	      # maximale Ausdehnung des Bereichs in Metern (min_hoch + max_hoch = Hochwert rechts oben)
kachel=1000          # Größe der Kachel in Metern

echo "############## Beginn `date +%c` ################" >> $LOGFILE 2>&1

./psql -d $DBNAME -U $DBUSER -c "

   -- Leeren der temporären Tabelle
   TRUNCATE grenzabschnitte1;

   -- Vorhandene Indices löschen
   DROP INDEX g1ixgeom;

   -- temporäre Tabelle füllen
   INSERT INTO grenzabschnitte1 
     SELECT makeline(pointn(foo.linestring, foo.count1), pointn(foo.linestring, foo.count1 + 1)) AS the_geom
     FROM ( SELECT generate_series(1, npoints(foofoo.linestring)) AS count1, foofoo.linestring
            FROM ( SELECT linefrompoly(alkobj_e_fla.the_geom) AS linestring
                   FROM alkobj_e_fla where folie = '001') foofoo) foo
     WHERE (foo.count1 + 1) <= npoints(foo.linestring);

   -- Tabelle analysieren für einen schnellen Zugriff
   VACUUM ANALYZE grenzabschnitte1;

   -- Indices wieder aufbauen für einen noch schnelleren Zugriff
   CREATE INDEX g1ixgeom ON grenzabschnitte1 USING gist (the_geom);

"  # Ende SQL

if test $? -eq 0
  then
    echo "`date +%c`   grenzabschnitte1 aktualisiert" >> $LOGFILE 2>&1
    
    
    # Lösche die doppelt vorhandenen wieder raus
    # portionsweiseweise, weil sonst die Prüfung zu lange dauert
    
    # äußere Schleife - erst die Zeilen von unten nach oben ...
    i=0
    while [ $i -lt $max_hoch + $kachel ]
    do
    
          z2=$(($min_hoch + $i))
          
    
    # innere Schleife - ... dann die Kacheln von links nach rechts
      j=0
      while [ $j -lt $max_rechts + $kachel ]
      do
    	z1=$(($min_rechts + $j))
    
    	z3=`expr $z1 + $kachel`
    	z4=`expr $z2 + $kachel`
    
    	echo $z1 " " $z2 ", " $z3 " " $z4
       
            sql="delete from grenzabschnitte1 where oid in (select a.oid from grenzabschnitte1 a left join grenzabschnitte1 b on length(a.the_geom)=length(b.the_geom) and a.the_geom=b.the_geom and x(pointn(a.the_geom,1)) > x(pointn(b.the_geom,1)) where a.the_geom && GeomFromText('POLYGON((${z1} ${z2},${z3} ${z2},${z3} ${z4},${z1} ${z4},${z1} ${z2}))', 2398) and b.the_geom && GeomFromText('POLYGON((${z1} ${z2},${z3} ${z2},${z3} ${z4},${z1} ${z4},${z1} ${z2}))', 2398));"
            #echo $sql
            ./psql -d $DBNAME -U $DBUSER -c "$sql"
      
            j=`expr $j + $kachel`
    
      done # Ende innere Schleife
    
          i=`expr $i + $kachel`
    
    done # Ende äußere Schleife
    
    echo "`date +%c`   grenzabschnitte1 auf doppelte Grenzen geprüft" >> $LOGFILE 2>&1
   
    
    # Übertragung der Grenzlängen in die endgültige Tabelle
    ./psql -d $DBNAME -U $DBUSER -c "
    
       -- Starten der Transaktion
       BEGIN;

       -- Index entfernen
       DROP INDEX g2ixgeom;
       
       -- Leeren der Tabelle
       TRUNCATE grenzabschnitte2;  
       
       -- Füllen der Tabelle
       INSERT INTO grenzabschnitte2 SELECT * FROM grenzabschnitte1;
       
       -- Beenden der Transaktion
       END;

    "  # Ende SQL
    
    if test $? -eq 0
      then
        echo "   grenzabschnitte2 erfolgreich aktualisiert" >> $LOGFILE 2>&1
      else
        echo ">> grenzabschnitte2 konnte nicht aktualisiert werden" >> $LOGFILE 2>&1
    fi
    
    ./psql -d $DBNAME -U $DBUSER -c "
    
       -- Bereinigen der Tabellenstruktur
       VACUUM ANALYZE grenzabschnitte2;

    "  # Ende SQL
    
    if test $? -eq 0
      then
        echo "   VACUUM erfolgreich durchgeführt" >> $LOGFILE 2>&1
      else
        echo ">> VACUUM konnte nicht durchgeführt werden" >> $LOGFILE 2>&1
    fi
    
    ./psql -d $DBNAME -U $DBUSER -c "
    
       CREATE INDEX g2ixgeom
       ON grenzabschnitte2
       USING gist
       (the_geom);
     
    "  # Ende SQL
    
    if test $? -eq 0
      then
        echo "   INDEX erfolgreich angelegt" >> $LOGFILE 2>&1
      else
        echo ">> INDEX konnte nicht angelegt werden" >> $LOGFILE 2>&1
    fi

  else
    echo "`date +%c`>> grenzabschnitte1 konnte nicht aktualisiert werden" >> $LOGFILE 2>&1
    echo "`date +%c`>> grenzabschnitte2 behält den alten Stand" >> $LOGFILE 2>&1
fi

echo " " >> $LOGFILE 2>&1
echo "############## Ende `date +%c` ################" >> $LOGFILE 2>&1
echo " " >> $LOGFILE 2>&1

exit 0