Opened 15 years ago

Closed 7 years ago

#52 closed task (fixed)

Calculation of all upstream river segments

Reported by: cedric Owned by: cedric
Priority: major Milestone:
Component: SIG-data Version: EDA2.0
Keywords: Cc: celine

Description

PLSQL function of
all identifiers upstream from one rivernode within a bassin (use pfafstetter)
all segments belonging to a subbassin

Change History (12)

comment:1 Changed 15 years ago by cedric

  • Status changed from new to accepted

which segments are upstream from the current segment ?

-- Vilaine
--- creating a second index to accelerate the requests
CREATE INDEX indexriversegments_wso_id
ON ccm21.riversegments
  USING btree
  (wso_id);

select * from ccm21.riversegments where wso_id=291146;
CREATE INDEX indexcatchments_wso_id
ON ccm21.catchments
  USING btree
  (wso_id);
select * from ccm21.catchments where wso_id=291146 and wso3_id=291828 ;
-- extraction du chiffre immediatement à gauche du segment selectionné
select * from ccm21.riversegments where wso_id=291146 order by cum_len desc;
select round(pfafstette) from ccm21.riversegments where wso_id=291146 order by cum_len desc;
select CAST(round(pfafstette) AS TEXT) from ccm21.riversegments where wso_id=291146 order by cum_len desc;
select substring(CAST(round(pfafstette) AS TEXT),2,2) from ccm21.riversegments where wso_id=291146 order by cum_len desc;

Create an index on all wso1_id
CREATE INDEX indexriversegments_wso1_id
ON ccm21.riversegments
  USING btree
  (wso1_id);
Create an index on all wso1..._id
CREATE INDEX indexcatchments_wso1_id
ON ccm21.catchments
  USING btree
  (wso1_id);
-- segment mer vilaine
select * from ccm21.riversegments where gid=234706;
-- le rang de strahler est 6
-- recupération de l'identifiant du bv correspondant
select wso6_id from ccm21.riversegments r join ccm21.catchments c on r.wso1_id=c.wso1_id where r.gid=234706;
-- 291146
select wso1_id from ccm21.catchments c  where wso6_id = 291146;
-- le m^ en plus compliqué
select wso1_id from ccm21.catchments c  where wso6_id = 
select wso1_id from ccm21.catchments c  where wso6_id =
	(select wso6_id from ccm21.riversegments r 
	join ccm21.catchments c on r.wso1_id=c.wso1_id 
	where r.gid=234706);
-- liste de wso1_id
-- extraction des segments contenus dans ces bassins
select * from ccm21.riversegments where wso1_id in 
(select wso1_id from ccm21.catchments c  where wso6_id = 291146) order by cum_len desc;
-- le m^ en plus compliqué
select * from ccm21.riversegments where wso1_id in 
(select wso1_id from ccm21.catchments c  where wso6_id =
	(select wso6_id from ccm21.riversegments r 
	join ccm21.catchments c on r.wso1_id=c.wso1_id 
	where r.gid=234706)
) ;
-- avec un cas adapté à chaque rang de strahler (marche pas)

select CASE 	WHEN r.strahler=5 THEN (select wso1_id from ccm21.catchments c  where wso5_id =
					(select wso5_id from ccm21.riversegments r 
					join ccm21.catchments c on r.wso1_id=c.wso1_id 
					where r.gid=234706)) 
		WHEN r.strahler=6 THEN (select wso1_id from ccm21.catchments c  where wso6_id =
					(select wso6_id from ccm21.riversegments r 
					join ccm21.catchments c on r.wso1_id=c.wso1_id 
					where r.gid=234706))
		ELSE 0
	END AS wso1_id
	FROM ccm21.riversegments r
	join ccm21.catchments c on r.wso1_id=c.wso1_id 
	where r.gid=234706;		
	
-- marche pas 



-- longueur de la cle pfstette
select character_length(CAST(round(pfafstette) AS TEXT)) from ccm21.riversegments where gid=234706; --1
-- selection de tous les noeuds de la Vilaine en amont du noeud choisi
select gid from ccm21.riversegments
	where wso_id=291146 
	and     CAST(substring(CAST(round(pfafstette) AS TEXT),
				1,
				1 -- voir ci dessus pour le calcul
				) as numeric	)	
	
	>1		;
-- le même en plus compliqué
select gid from ccm21.riversegments
	where wso_id=291146 
	and     CAST(substring(CAST(round(pfafstette) AS TEXT),
				1,
				1 -- voir ci dessus pour le calcul
				) as numeric	)	
	
	> (select pfafstette FROM ccm21.riversegments where gid=234706);

-- le même en plus compliqué
select gid from ccm21.riversegments
	where wso_id=291146 
	and     CAST(substring(CAST(round(pfafstette) AS TEXT),
				1,
				(select character_length(CAST(round(pfafstette) AS TEXT)) from ccm21.riversegments where gid=234706)
				) as numeric	)	
	
	> (select pfafstette FROM ccm21.riversegments where gid=234706);
select max (strahler) from ccm21.riversegments -- 10
-- le même en plus compliqué
select gid from ccm21.riversegments
	where wso1_id in 
	(select wso1_id from ccm21.catchments c  where wso6_id =
		(select wso6_id from ccm21.riversegments r 
		join ccm21.catchments c on r.wso1_id=c.wso1_id 
		where r.gid=234706)
	) 
	and     CAST(substring(CAST(round(pfafstette) AS TEXT),
				1,
				(select character_length(CAST(round(pfafstette) AS TEXT)) from ccm21.riversegments where gid=234706)
				) as numeric	)	
	
	> (select pfafstette FROM ccm21.riversegments where gid=234706);

-- le même en plus compliqué

comment:2 Changed 15 years ago by cedric

C'est vraiment dur
Voilà tout ce que j'ai réussi à programmer pour l'instant, si on lui passe gid_integer ça marche pas, RETURNS TABLE ne renvoit rien,
RETURN integer ne marche pas

SELECT gid, strahler FROM ccm21.riversegments where gid=234706;
DROP  FUNCTION IF EXISTS ccm21.pricatch_from_catch(gid_ numeric);
CREATE OR REPLACE FUNCTION ccm21.pricatch_from_catch(gid_ numeric) RETURNS SETOF ccm21.riversegments AS $$
DECLARE
BEGIN
    RETURN QUERY SELECT * FROM ccm21.riversegments where gid=gid_;
END;
$$
LANGUAGE 'plpgsql' ;
SELECT ccm21.pricatch_from_catch(234706);

Franchement c'est pas gagné, je pense qu'il va falloir écrire les résulats intermédiaires dans des tables crées à cet effet, cependant l'exemple ci dessous montre que ça marche

DROP  FUNCTION IF EXISTS ccm21.pricatch_from_catch(gid_ numeric);
CREATE OR REPLACE FUNCTION ccm21.pricatch_from_catch(gid_ numeric) RETURNS void AS $$
DECLARE
result RECORD;
BEGIN
    SELECT INTO result strahler FROM ccm21.riversegments where gid=gid_;
    RAISE NOTICE 'result is %' ,result;
END;
$$
LANGUAGE 'plpgsql' ;
SELECT ccm21.pricatch_from_catch(234706);
-- message NOTICE:  result is (6)

comment:3 Changed 15 years ago by cedric

One good function

-- it is necessary to create a type to hold the output type
DROP TYPE IF EXISTS wso1_id;
CREATE TYPE wso1_id as (wso1_id int);
-- In the function I tried to create a wso type to simplify the text. Pb you need to specify 
-- a tablename.columname%TYPE which will have to change along with the column, hence the big ugly script below...
DROP  FUNCTION IF EXISTS ccm21.pricatch_from_catch(gid_ numeric);
CREATE OR REPLACE FUNCTION ccm21.pricatch_from_catch(gid_ numeric) RETURNS setof int AS $$
	DECLARE
	result RECORD;
	pri wso1_id%rowtype;
	resultcount RECORD;
	BEGIN
	    SELECT INTO result strahler FROM ccm21.riversegments where gid=gid_;
	    RAISE NOTICE 'result is %' ,result;
	    IF result.strahler=2 THEN 	
		RAISE NOTICE 'WSO_2 used';
		SELECT INTO resultcount count(*) from  ccm21.catchments c  where wso2_id =
			(select wso2_id from ccm21.riversegments r 
			join ccm21.catchments c on r.wso1_id=c.wso1_id 
			where r.gid=gid_);
		RAISE NOTICE 'number of primary catchments in the catchment =%',resultcount; 
		for pri in select wso1_id from ccm21.catchments c  where wso2_id =
			(select wso2_id from ccm21.riversegments r 
			join ccm21.catchments c on r.wso1_id=c.wso1_id 
			where r.gid=gid_) loop
		pri.wso1_id=CAST(pri.wso1_id AS int8);
		return next pri.wso1_id;
		end loop;
		return;	
	    ELSIF result.strahler=3 THEN 
		RAISE NOTICE 'WSO_3 used';
		SELECT INTO resultcount count(*) from  ccm21.catchments c  where wso3_id =
			(select wso3_id from ccm21.riversegments r 
			join ccm21.catchments c on r.wso1_id=c.wso1_id 
			where r.gid=gid_);
		RAISE NOTICE 'number of primary catchments in the catchment =%',resultcount; 
		for pri in select wso1_id from ccm21.catchments c  where wso3_id =
			(select wso3_id from ccm21.riversegments r 
			join ccm21.catchments c on r.wso1_id=c.wso1_id 
			where r.gid=gid_) loop
		pri.wso1_id=CAST(pri.wso1_id AS int8);
		return next pri.wso1_id;
		end loop;
		return;	
	    ELSIF result.strahler=4 THEN
		RAISE NOTICE 'WSO_4 used';
		SELECT INTO resultcount count(*) from  ccm21.catchments c  where wso4_id =
			(select wso4_id from ccm21.riversegments r 
			join ccm21.catchments c on r.wso1_id=c.wso1_id 
			where r.gid=gid_);
		RAISE NOTICE 'number of primary catchments in the catchment =%',resultcount; 
		for pri in select wso1_id from ccm21.catchments c  where wso4_id =
			(select wso4_id from ccm21.riversegments r 
			join ccm21.catchments c on r.wso1_id=c.wso1_id 
			where r.gid=gid_) loop
		pri.wso1_id=CAST(pri.wso1_id AS int8);
		return next pri.wso1_id;
		end loop;
		return;			
	    ELSIF result.strahler=5 THEN
		RAISE NOTICE 'WSO_5 used';
		SELECT INTO resultcount count(*) from  ccm21.catchments c  where wso5_id =
			(select wso5_id from ccm21.riversegments r 
			join ccm21.catchments c on r.wso1_id=c.wso1_id 
			where r.gid=gid_);
		RAISE NOTICE 'number of primary catchments in the catchment =%',resultcount; 
		for pri in select wso1_id from ccm21.catchments c  where wso5_id =
			(select wso5_id from ccm21.riversegments r 
			join ccm21.catchments c on r.wso1_id=c.wso1_id 
			where r.gid=gid_) loop
		pri.wso1_id=CAST(pri.wso1_id AS int8);
		return next pri.wso1_id;
		end loop;
		return;	
	    ELSIF result.strahler=6 THEN
		RAISE NOTICE 'WSO_6 used';
		SELECT INTO resultcount count(*) from  ccm21.catchments c  where wso6_id =
			(select wso6_id from ccm21.riversegments r 
			join ccm21.catchments c on r.wso1_id=c.wso1_id 
			where r.gid=gid_);
		RAISE NOTICE 'number of primary catchments in the catchment =%',resultcount; 
		for pri in select wso1_id from ccm21.catchments c  where wso6_id =
			(select wso6_id from ccm21.riversegments r 
			join ccm21.catchments c on r.wso1_id=c.wso1_id 
			where r.gid=gid_) loop
		pri.wso1_id=CAST(pri.wso1_id AS int8);
		return next pri.wso1_id;
		end loop;
		return;		
	    ELSIF result.strahler=7 THEN
	    	RAISE NOTICE 'WSO_7 used';
		SELECT INTO resultcount count(*) from  ccm21.catchments c  where wso7_id =
			(select wso7_id from ccm21.riversegments r 
			join ccm21.catchments c on r.wso1_id=c.wso1_id 
			where r.gid=gid_);
		RAISE NOTICE 'number of primary catchments in the catchment =%',resultcount; 
		for pri in select wso1_id from ccm21.catchments c  where wso7_id =
			(select wso7_id from ccm21.riversegments r 
			join ccm21.catchments c on r.wso1_id=c.wso1_id 
			where r.gid=gid_) loop
		pri.wso1_id=CAST(pri.wso1_id AS int8);
		return next pri.wso1_id;
		end loop;
		return;	
	    ELSIF result.strahler=8 THEN
	    	RAISE NOTICE 'WSO_8 used';
		SELECT INTO resultcount count(*) from  ccm21.catchments c  where wso8_id =
			(select wso8_id from ccm21.riversegments r 
			join ccm21.catchments c on r.wso1_id=c.wso1_id 
			where r.gid=gid_);
		RAISE NOTICE 'number of primary catchments in the catchment =%',resultcount; 
		for pri in select wso1_id from ccm21.catchments c  where wso8_id =
			(select wso8_id from ccm21.riversegments r 
			join ccm21.catchments c on r.wso1_id=c.wso1_id 
			where r.gid=gid_) loop
		pri.wso1_id=CAST(pri.wso1_id AS int8);
		return next pri.wso1_id;
		end loop;
		return;		    	
	    ELSIF result.strahler=9 THEN
	    	RAISE NOTICE 'WSO_9 used';
		SELECT INTO resultcount count(*) from  ccm21.catchments c  where wso9_id =
			(select wso9_id from ccm21.riversegments r 
			join ccm21.catchments c on r.wso1_id=c.wso1_id 
			where r.gid=gid_);
		RAISE NOTICE 'number of primary catchments in the catchment =%',resultcount; 
		for pri in select wso1_id from ccm21.catchments c  where wso9_id =
			(select wso9_id from ccm21.riversegments r 
			join ccm21.catchments c on r.wso1_id=c.wso1_id 
			where r.gid=gid_) loop
		pri.wso1_id=CAST(pri.wso1_id AS int8);
		return next pri.wso1_id;
		end loop;
		return;		    	
	    ELSIF result.strahler=10 THEN
	    	RAISE NOTICE 'WSO_10 used';
		SELECT INTO resultcount count(*) from  ccm21.catchments c  where wso10_id =
			(select wso10_id from ccm21.riversegments r 
			join ccm21.catchments c on r.wso1_id=c.wso1_id 
			where r.gid=gid_);
		RAISE NOTICE 'number of primary catchments in the catchment =%',resultcount; 
		for pri in select wso1_id from ccm21.catchments c  where wso10_id =
			(select wso10_id from ccm21.riversegments r 
			join ccm21.catchments c on r.wso1_id=c.wso1_id 
			where r.gid=gid_) loop
		pri.wso1_id=CAST(pri.wso1_id AS int8);
		return next pri.wso1_id;
		end loop;
		return;		    		    		    		    	
	    ELSE RAISE NOTICE 'Strahler is larger than this function can handle !';
	    END IF;
	    RETURN;
	END;
	$$
LANGUAGE 'plpgsql' ;
COMMENT ON FUNCTION ccm21.pricatch_from_catch (gid_ numeric) IS 'Uses the gid from ccm21.riversegments, checks strahler, and select the catchment accordingly, then returns all primary catchments from this larger catchment';


--examples for use

SELECT  ccm21.pricatch_from_catch(234706); --la vilaine dans son ensemble
-- NOTICE:  result is (6)
-- NOTICE:  WSO_6 used
-- NOTICE:  number of primary catchments in the catchment =(759)
-- durée=2s
--
SELECT * from ccm21.catchments where wso1_id in (SELECT  ccm21.pricatch_from_catch(234706));

comment:4 Changed 15 years ago by cedric

And the second

.. still to be tested in Ggis

DROP TYPE IF EXISTS gid;
CREATE TYPE gid as (gid int);
DROP  FUNCTION IF EXISTS ccm21.upstream_segments(gid_ numeric);
CREATE OR REPLACE FUNCTION ccm21.upstream_segments(gid_ numeric) RETURNS setof int AS $$
	DECLARE
	gid_riversegments gid%rowtype;
	rec_pfafstette_length RECORD;
	pfafstette_chain_length int;
	pfafstette_value int8;
	BEGIN
	-- Valeur de la chaine pfafstette
	select into pfafstette_value CAST(round(pfafstette) AS int8) FROM ccm21.riversegments where gid=gid_;
	RAISE NOTICE 'pfafstette  is %', pfafstette_value;
	-- Longueur de la chaine pfafstette
	select into rec_pfafstette_length character_length(CAST(round(pfafstette) AS TEXT)) AS length 
		FROM ccm21.riversegments WHERE gid=gid_;
	pfafstette_chain_length=rec_pfafstette_length.length;	
	RAISE NOTICE 'pfafstette chain length is %', pfafstette_chain_length;
        -- loop
        for gid_riversegments in select gid from ccm21.riversegments
		where wso1_id in (SELECT  ccm21.pricatch_from_catch(gid_))
		and CAST(
			substring(CAST(round(pfafstette) AS TEXT),1,pfafstette_chain_length)
		as numeric)> pfafstette_value loop
	        gid_riversegments.gid=CAST(gid_riversegments.gid AS int8);
		return next gid_riversegments.gid;
		end loop;
		return;	
	END;
$$
LANGUAGE 'plpgsql' ;
COMMENT ON FUNCTION ccm21.upstream_segments(gid_ numeric) IS 'This function uses ccm21.pricatch_from_catch to get all the segments from a catchment and searches all segments with a pfafstetter higher than the segment in the bassin. The pfafstette are trucated to the level of the segment used as input in the function';
An example for use
select ccm21.upstream_segments(234706);

comment:5 Changed 15 years ago by cedric

  • Resolution set to fixed
  • Status changed from accepted to closed

comment:6 Changed 15 years ago by cedric

DEFECT =>
Should be >= instead of > and everything is OK

comment:7 Changed 15 years ago by cedric

Big problems linked with the fact that pfafstette lengths varies along one axis 191 the 18 19 and 18<191
Solution, writing a new table, with good lenght to pfafstette
DROP TYPE IF EXISTS gid;
CREATE TYPE gid as (gid int);
DROP FUNCTION IF EXISTS ccm21.upstream_segments(gid_ numeric);
CREATE OR REPLACE FUNCTION ccm21.upstream_segments(gid_ numeric) RETURNS setof int AS $$

DECLARE
gid_riversegments gid%ROWTYPE;
rec_pfafstette_length RECORD;
-- pfafstette_chain_length int;
pfafstette_max_chain_length int;
pfafstette_value int8;
BEGIN
-- filling a new table with the results from a catchment


DROP TABLE IF EXISTS ccm21.upstream_riversegments;
CREATE TABLE ccm21.upstream_riversegments AS SELECT * FROM ccm21.riversegments
--ALTER TABLE ccm21.upstream_riversegments ADD PRIMARY KEY..
where wso1_id in (SELECT ccm21.mycatchment(gid_));
-- pfafstette chain length
--select into pfafstette_chain_length character_length(CAST(round(pfafstette) AS TEXT)
-- FROM ccm21.upstream_riversegments WHERE gid=gid_;
-- RAISE NOTICE 'pfafstette chain length is %', pfafstette_chain_length;
-- pfafstette max chain length in the selected basin
select into pfafstette_max_chain_length max(character_length(CAST(round(pfafstette) AS TEXT)))

FROM ccm21.upstream_riversegments;

-- pfafstette value
-- the chain is lengthened to the pfafstette max chain length
Update ccm21.upstream_riversegments set pfafstette =
floor(pfafstette)*power(10,(pfafstette_max_chain_length-character_length(CAST(floor(pfafstette) AS TEXT) )));
select into pfafstette_value CAST(round(pfafstette) AS int8) FROM ccm21.riversegments where gid=gid_;
RAISE NOTICE 'pfafstette is %', pfafstette_value;

-- extracting the riversegments corresponding to the basin
FOR gid_riversegments in select gid from ccm21.upstream_riversegments where pfafstette >= pfafstette_value loop
gid_riversegments.gid=CAST(gid_riversegments.gid AS int8);
return next gid_riversegments.gid;
end loop;
return;
END;

$$
LANGUAGE 'plpgsql' ;
COMMENT ON FUNCTION ccm21.upstream_segments(gid_ numeric) IS 'This function uses ccm21.mycatchment to get all the segments from a catchment and searches all segments with a pfafstetter higher than the segment in the bassin. The pfafstette are trucated to the level of the segment used as input in the function';

comment:8 Changed 15 years ago by cedric

Still to do, add a primary key

comment:9 Changed 15 years ago by laurent

une fonction très pratique dans les contrib (tablefunc) de postgresql qui fait cela : connectby http://docs.postgresqlfr.org/8.4/tablefunc.html
Voir aussi la fonction parcourt graph dans l'aide de WITH RECURSIVE http://docs.postgresqlfr.org/8.4/queries-with.html

comment:10 Changed 15 years ago by laurent

  • Resolution fixed deleted
  • Status changed from closed to reopened

comment:11 Changed 7 years ago by cedric

  • Milestone Data integration deleted

Milestone Data integration deleted

comment:12 Changed 7 years ago by cedric

  • Resolution set to fixed
  • Status changed from reopened to closed
Note: See TracTickets for help on using tickets.