2

After reading this post How to georeference a vector layer with control points? I decided to give it a try to move/rotation/scale a polygon using Postgis.

I have two polygons, polygon_from and polygon_to

Below are their 4 vertices' longtitude-latitude. I intend to move polygon_from to polygon_to position.

 gid |      from_x      |      from_y      |       to_x       |       to_y
-----+------------------+------------------+------------------+------------------
   1 | 103.771789705736 | 1.29874218723039 | 103.771893494233 | 1.29880598084823
   2 | 103.7718553679   | 1.2987732903611  | 103.771964782989 | 1.29879195095702
   3 | 103.771887622999 | 1.29870762819626 | 103.771951679355 | 1.29871997730869
   4 | 103.771820808866 | 1.29867710104945 | 103.771879806787 | 1.29873515522057

I copy and pasted a function (below) from this website http://casoilresource.lawr.ucdavis.edu/drupal/book/export/html/532

-- Function: trans_param(IN sql text, OUT a double precision, OUT b double precision, OUT d double precision, OUT e double precision, OUT xoff double precision, OUT yoff double precision)

-- DROP FUNCTION trans_param(IN sql text, OUT a double precision, OUT b double precision, OUT d double precision, OUT e double precision, OUT xoff double precision, OUT yoff double precision);

CREATE OR REPLACE FUNCTION trans_param(IN sql text, OUT a double precision, OUT b double precision, OUT d double precision, OUT e double precision, OUT xoff double precision, OUT yoff double precision) AS $BODY$ DECLARE
        cc_row record;

        cc_det float;

        inv_cc00 float;
        inv_cc01 float;
        inv_cc02 float;
        inv_cc11 float;
        inv_cc12 float;
        inv_cc22 float;

BEGIN
        EXECUTE 'SELECT 
        count(a.to_x) as cc00,
        sum(a.to_x) as cc01,
        sum(a.to_y) as cc02,
        sum(a.to_x * b.to_x) as cc11,
        sum(a.to_x * b.to_y) as cc12,
        sum(a.to_y * b.to_y) as cc22,
        sum(a.from_y) as aa0,
        sum(a.from_y * b.to_x) as aa1,
        sum(a.from_y * b.to_y) as aa2,
        sum(a.from_x) as bb0,
        sum(a.from_x * b.to_x) as bb1,
        sum(a.from_x * b.to_y) as bb2
        from (' || sql || ') a,(' || sql|| ') b WHERE a.gid = b.gid' INTO cc_row;


        SELECT INTO cc_det 
                cc_row.cc00*cc_row.cc11*cc_row.cc22 + 
                cc_row.cc01*cc_row.cc12*cc_row.cc02 + 
                cc_row.cc02*cc_row.cc01*cc_row.cc12 - 
                cc_row.cc00*cc_row.cc12*cc_row.cc12 - 
                cc_row.cc01*cc_row.cc01*cc_row.cc22 - 
                cc_row.cc02*cc_row.cc11*cc_row.cc02;

        SELECT INTO inv_cc00 (cc_row.cc11*cc_row.cc22-cc_row.cc12*cc_row.cc12)/cc_det;
        SELECT INTO inv_cc01 (cc_row.cc12*cc_row.cc02-cc_row.cc01*cc_row.cc22)/cc_det;
        SELECT INTO inv_cc02 (cc_row.cc01*cc_row.cc12-cc_row.cc11*cc_row.cc02)/cc_det;
        SELECT INTO inv_cc11 (cc_row.cc00*cc_row.cc22-cc_row.cc02*cc_row.cc02)/cc_det;
        SELECT INTO inv_cc12 (cc_row.cc01*cc_row.cc02-cc_row.cc00*cc_row.cc12)/cc_det;
        SELECT INTO inv_cc22 (cc_row.cc00*cc_row.cc11-cc_row.cc01*cc_row.cc01)/cc_det;

        SELECT INTO a cc_row.bb0*inv_cc01+cc_row.bb1*inv_cc11+cc_row.bb2*inv_cc12;
        SELECT INTO b cc_row.bb0*inv_cc02+cc_row.bb1*inv_cc12+cc_row.bb2*inv_cc22;
        SELECT INTO d cc_row.aa0*inv_cc01+cc_row.aa1*inv_cc11+cc_row.aa2*inv_cc12;
        SELECT INTO e cc_row.aa0*inv_cc02+cc_row.aa1*inv_cc12+cc_row.aa2*inv_cc22;
        SELECT INTO xoff cc_row.bb0*inv_cc00+cc_row.bb1*inv_cc01+cc_row.bb2*inv_cc02;
        SELECT INTO yoff cc_row.aa0*inv_cc00+cc_row.aa1*inv_cc01+cc_row.aa2*inv_cc02;

END;

$BODY$   LANGUAGE 'plpgsql' VOLATILE; ALTER FUNCTION trans_param(IN sql text, OUT a double precision, OUT b double precision, OUT d double precision, OUT e double precision, OUT xoff double precision, OUT yoff double precision) OWNER TO postgres;

SELECT trans_param('select * from above_table');

After applying the PostGIS function, I obtained a,b,d,e,x_off,y_off parameters as below

(0.0824142911101262,47.5134245217341,0.00103114836995566,0.594639499427103,-9.93126928180573,-0.124261944539057)

However when I apply those parameters to the polygon_from, it goes to a very strange lat/long position and becomes an irregular shape (not what I expected).

Update polygon_table Set the_geom=(Select st_affine(the_geom,0.0824142911101262,47.5134245217341,0.00103114836995566,0.594639499427103,-9.93126928180573,-0.124261944539057) Where gid=polygon_from_gid;

I am wondering how to use the above trans_param() function and those ST_Affine parameters to shift/rotate/scale a polygon?

  • I had a look at casoilresource.lawr.ucdavis.edu/drupal/book/export/html/532 (apparently by someone named Bruce) where it compares transformations from Grass, R, PostGIS and his own plpgsql trans_param() function. It seems that those transfors each use different naming or ordering of parameters. Thus, he, you or I may have got confused somewhere! Some parameters may need re-ordering (but I cannot say which/how). – Martin F Nov 23 '13 at 2:31
2

I had the same problem: the result of the trans_params function defined as we can read here https://casoilresource.lawr.ucdavis.edu/software/postgis-spatially-enabled-relational-database-sytem/affine-transformation-operations-postgis/user-defined-function-plpgsql-compute-transformation-parameters/ was incorrect and very different from the same result for example I obtained with the ogr2org command (with GPC parameter). I would apply that function to more than 50 control points in UTM metric coordinates. I discovered finally that the problem was due to the PRECISION of the variables declared in that function: just changing the type of variables from "float" to "numeric" (mainly for the "cc_det" variable) solved the problem.

I also simplified a bit the code, always taking in account the v.transform code from GRASS source code https://grass.osgeo.org/programming6/transform_8c_source.html, thanks to OPEN SOURCE !

finally the correct code is:

-- DROP FUNCTION public.trans_param(text);

CREATE OR REPLACE FUNCTION public.trans_param(
    IN sql text,
    OUT a numeric,
    OUT b numeric,
    OUT d numeric,
    OUT e numeric,
    OUT xoff numeric,
    OUT yoff numeric)
  RETURNS record AS
$BODY$
DECLARE
        cc_row record;

        cc_det numeric;

        inv_cc00 numeric;
        inv_cc01 numeric;
        inv_cc02 numeric;
        inv_cc11 numeric;
        inv_cc12 numeric;
        inv_cc22 numeric;

BEGIN
        EXECUTE 'SELECT
        count(a.to_x) as cc00,
        sum(a.to_x)::numeric as cc01,
        sum(a.to_y)::numeric as cc02,
        sum(a.to_x * a.to_x)::numeric as cc11,
        sum(a.to_x * a.to_y)::numeric as cc12,
        sum(a.to_y * a.to_y)::numeric as cc22,
        sum(a.from_y)::numeric as aa0,
        sum(a.from_y * a.to_x)::numeric as aa1,
        sum(a.from_y * a.to_y)::numeric as aa2,
        sum(a.from_x)::numeric as bb0,
        sum(a.from_x * a.to_x)::numeric as bb1,
        sum(a.from_x * a.to_y)::numeric as bb2
        from (' || sql || ') a' INTO cc_row;


        SELECT INTO cc_det
                (cc_row.cc00*cc_row.cc11*cc_row.cc22 +
                cc_row.cc01*cc_row.cc12*cc_row.cc02 +
                cc_row.cc02*cc_row.cc01*cc_row.cc12 -
                cc_row.cc00*cc_row.cc12*cc_row.cc12 -
                cc_row.cc01*cc_row.cc01*cc_row.cc22 -
                cc_row.cc02*cc_row.cc11*cc_row.cc02)::numeric;

        SELECT INTO inv_cc00 ((cc_row.cc11*cc_row.cc22-cc_row.cc12*cc_row.cc12)/cc_det)::numeric;
        SELECT INTO inv_cc01 ((cc_row.cc12*cc_row.cc02-cc_row.cc01*cc_row.cc22)/cc_det)::numeric;
        SELECT INTO inv_cc02 ((cc_row.cc01*cc_row.cc12-cc_row.cc11*cc_row.cc02)/cc_det)::numeric;
        SELECT INTO inv_cc11 ((cc_row.cc00*cc_row.cc22-cc_row.cc02*cc_row.cc02)/cc_det)::numeric;
        SELECT INTO inv_cc12 ((cc_row.cc01*cc_row.cc02-cc_row.cc00*cc_row.cc12)/cc_det)::numeric;
        SELECT INTO inv_cc22 ((cc_row.cc00*cc_row.cc11-cc_row.cc01*cc_row.cc01)/cc_det)::numeric;

        SELECT INTO a (cc_row.bb0*inv_cc01+cc_row.bb1*inv_cc11+cc_row.bb2*inv_cc12)::numeric;
        SELECT INTO b (cc_row.bb0*inv_cc02+cc_row.bb1*inv_cc12+cc_row.bb2*inv_cc22)::numeric;
        SELECT INTO d (cc_row.aa0*inv_cc01+cc_row.aa1*inv_cc11+cc_row.aa2*inv_cc12)::numeric;
        SELECT INTO e (cc_row.aa0*inv_cc02+cc_row.aa1*inv_cc12+cc_row.aa2*inv_cc22)::numeric;
        SELECT INTO xoff (cc_row.bb0*inv_cc00+cc_row.bb1*inv_cc01+cc_row.bb2*inv_cc02)::numeric;
        SELECT INTO yoff (cc_row.aa0*inv_cc00+cc_row.aa1*inv_cc01+cc_row.aa2*inv_cc02)::numeric;

END;

$BODY$
  LANGUAGE plpgsql VOLATILE
  COST 100;

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.