-- Vytvoreni vrstvy bunkry

CREATE TABLE f10.bunkry(
				gid SERIAL PRIMARY KEY,
				osm_id INTEGER,
				name TEXT,
				ref TEXT,
				way GEOMETRY);
				
INSERT INTO f10.bunkry (osm_id,ref,way) 
	SELECT  osm_id, name, way 
	FROM czech_point 
	WHERE military = 'bunker' AND name LIKE '%/%';
	
INSERT INTO f10.bunkry (osm_id,name,way)  
	SELECT  osm_id, name, way 
	FROM czech_point 
	WHERE military = 'bunker' AND (name NOT LIKE '%/%' OR name IS NULL);

INSERT INTO f10.bunkry (osm_id,name,ref,way)  
	SELECT  osm_id, name,ref, ST_Centroid(way) 
	FROM czech_polygon 
	WHERE military = 'bunker';
	
SELECT Populate_geometry_columns('f10.bunkry'::regclass);
CREATE INDEX bunkry_gist ON f10.bunkry USING GIST (way);

-- Vytvoreni vrstvy ujezdy

CREATE TABLE f10.ujezdy AS
	SELECT name,military,way 
	  FROM czech_polygon 
	 WHERE military = 'danger_area' AND name LIKE '%újezd%';
 
ALTER TABLE f10.ujezdy ADD COLUMN gid serial;
ALTER TABLE f10.ujezdy ADD PRIMARY KEY (gid);
 
SELECT Populate_geometry_columns('f10.ujezdy'::regclass);
CREATE INDEX f10_ujezdy_gist ON f10.ujezdy USING GIST (way);
 
-- Vytvoreni vrstvy lesu

CREATE TABLE f10.lesni_porosty AS
 	SELECT osm_id,name,way 
 	  FROM czech_polygon 
 	 WHERE landuse = 'forest';
 
ALTER TABLE f10.lesni_porosty ADD COLUMN gid SERIAL; 
ALTER TABLE f10.lesni_porosty ADD PRIMARY KEY (gid);
 
SELECT Populate_geometry_columns('f10.lesni_porosty'::regclass);
CREATE INDEX lesni_porosty_gist ON f10.lesni_porosty USING GIST (way);

--------------------------------------------------------------------------------------
-- Validace nevalidnich polygonu lesu

CREATE TABLE f10.valid AS 
SELECT gid, ST_BuildArea(ST_Collect(way)) AS way
  FROM
	(SELECT gid, ST_BuildArea(ST_Union(ring.way,ST_StartPoint(ring.way))) AS way
		FROM
			(SELECT DISTINCT gid, ST_ExteriorRing((ST_DumpRings(way)).geom) AS way 
				FROM f10.lesni_porosty WHERE NOT ST_IsValid(way)) AS ring) AS after_union
GROUP BY gid;

ALTER TABLE f10.valid ADD PRIMARY KEY (gid);

SELECT Populate_Geometry_Columns('f10.valid'::regclass);
CREATE INDEX valid_gist ON f10.valid USING GIST (way);

-- Kontrola zvalidovanych polygonu
SELECT n.gid,
	ST_NRings(n.way) AS r_neval,
	ST_NRings(v.way) AS r_val,
	ST_NRings(n.way) - ST_NRings(v.way) AS roz,
	ROUND( ( (ST_Area(n.way) - ST_Area(v.way)) / ST_Area(n.way) )*100) AS pl_proc,
	ST_GeometryType(v.way) AS type,
	ST_IsValid(n.way) AS b,
	ST_IsValid(v.way) AS a,
	ST_IsValidReason(n.way) AS reason
		
  FROM f10.lesni_porosty AS n
  JOIN f10.valid AS v
    ON v.gid = n.gid
 ORDER BY n.gid;
 
 
-- Update lesu
ALTER TABLE f10.lesni_porosty DROP CONSTRAINT enforce_geotype_way;

UPDATE f10.lesni_porosty
   SET way = f10.valid.way
  FROM f10.valid
 WHERE f10.lesni_porosty.gid = f10.valid.gid;
 
UPDATE f10.lesni_porosty
   SET way = ST_Multi(way);
   
SELECT Populate_geometry_columns('f10.lesni_porosty'::regclass);

------------------------------------------------------------------------------
-- Prekryt lesu - zjisteni intersectu


--CREATE AGGREGATE f10.array_accum (
--    sfunc = array_append,
--    basetype = anyelement,
--    stype = anyarray,
--   initcond = '{}'
--);

CREATE TABLE f10.lesy_intersect AS 

	SELECT g1.gid AS g1_gid,
	--f10.array_accum(g2.gid) AS g2_gid,
	-- pouze pro informaci
	ST_Union(ST_Buffer(ST_Intersection(g1.way, g2.way),0.0)) AS inter
  	  FROM f10.lesni_porosty g1, f10.lesni_porosty g2
	 WHERE g1.gid < g2.gid 
	 	AND st_intersects(g1.way, g2.way) 
  		AND NOT ST_Touches(g1.way, g2.way)
		AND g1.gid != 40686
 GROUP BY g1.gid
 ORDER BY g1.gid
 
-- samostatne insertovan intersect pro polygon gid 40686 kvuli hlaseni vyjimky unionem
INSERT INTO f10.lesy_intersect
SELECT g1.gid AS g1_gid,
	--f10.array_accum(g2.gid) AS g2_gid,
	-- pouze pro informaci
       ST_Union(ST_Buffer(ST_Intersection(
                            ST_SnapToGrid(g1.way, 0.001), g2.way
                            ), 0.0)) AS inter
  FROM f10.lesni_porosty g1, f10.lesni_porosty g2
 WHERE g1.gid < g2.gid
       AND ST_Intersects(g1.way, g2.way)
       AND NOT ST_Touches(g1.way, g2.way)
       AND g1.gid = 40686
 GROUP BY g1.gid;

-- Prekryt lesu - odecteni zvetseneho intersectu
CREATE TABLE f10.lesy_no_intersect_buffer AS
	SELECT A.gid, ST_Multi(ST_Difference(A.way,ST_Buffer(B.inter,0.01)) ) AS way
	  FROM f10.lesni_porosty AS A
	  JOIN f10.lesy_intersect AS B
		 ON A.gid = B.g1_gid;


-- kontrola

SELECT	gid, 
			ST_GeometryType(way),
			ST_Area(way),
			ST_NPoints(way),
			ST_IsValidReason(way) 
  FROM f10.lesy_no_intersect_buffer 
 WHERE ST_GeometryType(way) <> 'ST_MultiPolygon' 
 	 OR ST_Area(way) < 10 
 	 OR NOT ST_IsValid(way);

-- vytvoreni zalohy zvalidovanych lesu
CREATE TABLE f10.lesni_porosty_winter AS -- With INTERsect
SELECT * FROM f10.lesni_porosty;

-- update lesu
UPDATE f10.lesni_porosty 
   SET way = A.way
  FROM f10.lesy_no_intersect_buffer AS A
 WHERE A.gid = f10.lesni_porosty.gid
   AND ST_GeometryType(A.way) = 'ST_MultiPolygon' -- pro jistotu
   AND ST_Area(A.way) > 1.0
   AND ST_IsValid(A.way); -- pro jistotu

DELETE FROM f10.lesni_porosty
 USING f10.lesy_no_intersect_buffer AS A
 WHERE A.gid = f10.lesni_porosty.gid
   AND ST_Area(A.way) <= 1.0;
   
-- kontrola porovnanim plochy 
SELECT
(
 SELECT SUM(ST_Area(way))
   FROM f10.lesni_porosty
  WHERE ST_Within(way, 
                  ST_SetSRID(ST_MakeBox2D(
                        ST_Point(1500000, 6300000),
                        ST_Point(1560000, 6360000)
                  ) , 900913) )
) AS new,
(
 SELECT ST_Area(ST_Union(way))
   FROM f10.lesni_porosty_winter
  WHERE ST_Within(way,
                  ST_SetSRID(ST_MakeBox2D(
                        ST_Point(1500000, 6300000),
                        ST_Point(1560000, 6360000)
                  ) , 900913) )
) AS old;

-- Zjisteni dalsich zbyvajicich intersectu, jejich odecteni a nasledny update se provede opakovanim postupu.

----------------------------------------------------------------------------------------------------
-- Dotazy


-- Kolik bunkrů je v lese?}


SELECT COUNT(*)
  FROM f10.bunkry AS bunkry 
  JOIN f10.lesni_porosty AS lesy 
    ON ST_Within(bunkry.way, lesy.way);

-- count
-- -----
--    127
-- (1 row)


-- Který újezd má největší rozlohu a kolik?}


SELECT name, ROUND(ST_Area(way)/1e6) AS area
  FROM f10.ujezdy AS ujezdy
 ORDER BY area DESC
 LIMIT 1;



--           name           | area 
-- -------------------------+------
--  Vojenský újezd Hradiště |  809 
-- (1 row)

-- Poznámka:
-- Skutečná rozloha je 332 $km^2$. Rozdíl způsoben nezavedením opravy ze zobrazení (platí pro všechny otázky na plochu).


-- Kolik procent vojenského újezdu Hradiště pokrývá les?}


SELECT ROUND(((lesy_area/area)*100)::NUMERIC, 2)
FROM (
 SELECT SUM(ST_Area(ST_Intersection(ujezdy.way, lesy.way))) AS lesy_area,
        ST_Area(ujezdy.way) AS area
   FROM f10.ujezdy AS ujezdy
   JOIN f10.lesni_porosty AS lesy
     ON ujezdy.name = 'Vojensky ujezd Hradiste'
        AND ST_Intersects(ujezdy.way, lesy.way)
  GROUP BY ujezdy.way
) AS lesy_hrad;



-- round
-- -------
--  16.86
-- (1 row)




-- Jaká je vzdálenost mezi nejzápadnějším a nejvýchodnějším bunkrem v ČR?}


SELECT ROUND(ST_Distance(
       (SELECT way 
          FROM f10.bunkry AS bunkry 
         ORDER BY ST_X(bunkry.way) DESC 
         LIMIT 1),
       (SELECT way 
          FROM f10.bunkry 
         ORDER BY ST_X(bunkry.way) ASC 
         LIMIT 1)
        )/1e3) AS vzdalenost;


--  vzdalenost 
-- ------------
--         524
-- (1 row)

-- Poznámka:
-- Jistě je zde také rozdíl oproti skutečnosti způsobený nezavedením opravy ze zobrazení (platí pro všechny otázky na vzdálenost).





-- Který vojenský újezd má největší poměr obvodu ku ploše (poměr uveďte na 5 desetinných míst)?}


SELECT name, ROUND((ST_Perimeter(way)/ST_Area(way))::NUMERIC,5) AS hodnota 
  FROM f10.ujezdy AS ujezdy
 ORDER BY hodnota DESC LIMIT 1;



--           name          | hodnota 
-- ------------------------+----------- 
--  Vojenský újezd Březina | 0.00045
-- (1 row)






-- Vypište souřadnice bunkru, který se nachází nejdál od vojenského újezdu Brdy.}


SELECT ST_X(bunkr.way) AS X, 
       ST_Y(bunkr.way) AS Y 
  FROM (
    SELECT bunkry.way,ST_Distance(bunkry.way, ujezdy.way) AS vzdalenost  
      FROM f10.bunkry AS bunkry 
      JOIN f10.ujezdy AS  ujezdy 
        ON ujezdy.name = 'Vojenský újezd Brdy' 
     ORDER BY vzdalenost DESC LIMIT 1) AS bunkr;



--         x         |        y
-- ------------------+------------------
--  1996460.74801165 | 6437929.75357603
-- (1 row)





-- Kolik lesů má celou svou plochu do 10 km od vojenského újezdu Brdy?}


 SELECT COUNT(*) 
   FROM f10.lesni_porosty AS lesy 
   JOIN ST_Buffer(
          (SELECT way 
             FROM f10.ujezdy 
            WHERE name='Vojenský újezd Brdy'),1e4) AS buffer 
     ON ST_Within(lesy.way, buffer);



--  count
-- -------
--    322
-- (1 row)

-- Poznámky: Příkaz se dlouho se vykonává (josef: 79153,033 ms). Opět problém s jednotkami, tentokrát však už ve znění otázky.


-- Jaká je výměra lesa ve vzdálenosti 10 km kolem nej\-sever\-něj\-šího bunkru?}


SELECT ROUND(SUM(ST_Area(ST_Intersection(lesy.way, buffer))/1e6)) 
  FROM f10.lesni_porosty AS lesy 
  JOIN ST_Buffer(
        (SELECT way 
           FROM f10.bunkry 
           ORDER BY ST_Y(way) DESC LIMIT 1)
        ,1e4) AS buffer
    ON ST_Intersects(lesy.way, buffer);



--  round 
-- -------
--    270
-- (1 row)


-- Jaká je výměra nejmenšího lesa v ČR?}


SELECT ROUND(ST_Area(lesy.way)) AS area 
  FROM f10.lesni_porosty AS lesy 
 ORDER BY area ASC LIMIT 1;



--  area 
-- ------
--   106
-- (1 row)




-- Jaká je celková plocha všech újezdů v ČR?}


SELECT ROUND(SUM(ST_Area(way))/1e6) AS area
  FROM f10.ujezdy;



--  area
-- -------
--   3058
-- (1 row)

 
