PostGIS – How to Convert LineStringZ to LineStringM Using Alter Column

linestringpostgispostgresqlsql

I am attempting to import linestring data with or without m values, then test whether the data has m values, then add m values to the geometry if they are missing. So far, I have working tests for the m values, and I can create the new geometry with m values included when using a select statement. However, I am trying to use this test in an alter table/alter column statement, which doesn't allow subqueries. I needed the subquery that I am using in order to convert a linestringz to a linestringm. I thought there would be a postgis function to swap coordinates between z and m, like st_force3dm(geom, st_z(geom)), which obviously doesn't work because st_z requires a point. Is there a way to swap a linestringz to linestringm using the z values as the new m values, and within an alter table statement?

This is the code I have so far. It tests whether m values exists, whether they are all 0, whether z values exist and whether they are all 0. If it finds m values that are not all 0, it returns the original geometry. If it finds no z or m values, or z and m values are all 0, it uses the geometry length to create m values. If it finds z values that are not all 0, it creates a linestringm using the z values as the m values. This last part doesn't work in an alter table statement because I have used a subquery to assign the z values to m values. I'm open to suggestions about better ways to achieve this task.

Doesn't work due to subquery in alter table:

alter table mt_altis_statewide_routes alter column geom type geometry(multilinestringm,4326) using 
    case when ((st_m(st_startpoint(geom)) = 0 and st_m(st_endpoint(geom)) = 0) or st_m(st_startpoint(geom)) is null) 
        and ((st_z(st_startpoint(geom)) = 0 and st_z(st_endpoint(geom)) = 0) or st_z(st_endpoint(geom)) is null) 
        then st_addmeasure(geom,0,length_miles(geom)) 
    when ((st_m(st_startpoint(geom)) = 0 and st_m(st_endpoint(geom)) = 0) or st_m(st_startpoint(geom)) is null) 
        and (st_z(st_startpoint(geom)) != 0 or st_z(st_endpoint(geom)) != 0) and st_z(st_startpoint(geom)) is not null 
        then (select st_setsrid(st_makeline(st_makepointm(st_x(geom),st_y(geom),st_z(geom)) order by st_z(geom)),4326) from (select (st_dumppoints(geom)).geom) as f) 
    when ((st_z(st_startpoint(geom)) = 0 and st_z(st_endpoint(geom)) = 0) or st_z(st_startpoint(geom)) is null) 
        and (st_m(st_startpoint(geom)) != 0 or st_m(st_endpoint(geom)) != 0) and st_m(st_startpoint(geom)) is not null 
        then geom end;

Works as select statement:

select 
    case when ((st_m(st_startpoint(geom)) = 0 and st_m(st_endpoint(geom)) = 0) or st_m(st_startpoint(geom)) is null) 
        and ((st_z(st_startpoint(geom)) = 0 and st_z(st_endpoint(geom)) = 0) or st_z(st_endpoint(geom)) is null) 
        then st_addmeasure(geom,0,length_miles(geom)) 
    when ((st_m(st_startpoint(geom)) = 0 and st_m(st_endpoint(geom)) = 0) or st_m(st_startpoint(geom)) is null) 
        and (st_z(st_startpoint(geom)) != 0 or st_z(st_endpoint(geom)) != 0) and st_z(st_startpoint(geom)) is not null 
        then (select st_setsrid(st_makeline(st_makepointm(st_x(geom),st_y(geom),st_z(geom)) order by st_z(geom)),4326) from (select (st_dumppoints(geom)).geom) as f) 
    when ((st_z(st_startpoint(geom)) = 0 and st_z(st_endpoint(geom)) = 0) or st_z(st_startpoint(geom)) is null) 
        and (st_m(st_startpoint(geom)) != 0 or st_m(st_endpoint(geom)) != 0) and st_m(st_startpoint(geom)) is not null 
        then geom end geom 
from mt_state_maintained_routes msmr ;

Best Answer

As said in the comments it would be necessary put your conversion code in a function like this:

CREATE OR REPLACE FUNCTION make_multilinestringm(geometry) RETURNS geometry(multilinestringm,4326) AS $$
DECLARE 
    geom geometry;  
BEGIN
    geom = $1;
   case  
       when (( st_m(st_startpoint(geom)) = 0 and st_m(st_endpoint(geom)) = 0) or st_m(st_startpoint(geom)) is null) 
            and ((st_z(st_startpoint(geom)) = 0 and st_z(st_endpoint(geom)) = 0) or st_z(st_endpoint(geom)) is null) 
            then return st_addmeasure(geom,0,length_miles(geom));        
       when ((st_m(st_startpoint(geom)) = 0 and st_m(st_endpoint(geom)) = 0) or st_m(st_startpoint(geom)) is null) 
            and (st_z(st_startpoint(geom)) != 0 or st_z(st_endpoint(geom)) != 0) and st_z(st_startpoint(geom)) is not null 
            then return (select st_setsrid(st_makeline(st_makepointm(st_x(geom),st_y(geom),st_z(geom)) order by st_z(geom)),4326) from (select (st_dumppoints(geom)).geom) as f); 
       when ((st_z(st_startpoint(geom)) = 0 and st_z(st_endpoint(geom)) = 0) or st_z(st_startpoint(geom)) is null) 
            and (st_m(st_startpoint(geom)) != 0 or st_m(st_endpoint(geom)) != 0) and st_m(st_startpoint(geom)) is not null 
            then return geom; 
    end case;    
END;
$$
STRICT
LANGUAGE plpgsql IMMUTABLE;

and then execute alter sentence using it

alter table mt_altis_statewide_routes alter column geom type geometry(multilinestringm,4326) using make_multilinestringm(geom);

Related Solutions

[GIS] ny character limit in st_union and array_agg

The maximum size of any Postgres object, including the results of array_agg, ST_Collect, and other array-producing aggregates, is 1 GB.

If you hit this limit, these functions won't silently return null. They'll fail with a message like this:

ERROR: array size exceeds the maximum allowed (1073741823).

For what it's worth, this object size limit roughly corresponds to a 2D LineString with 67 million vertices. I wrote a post about this topic a while ago with a few more details.

[GIS] How to make projected geometry and update it to postgresql

As user30184 outlined: It's a very common task, and there's plenty of documentation: ST_Transform.

To apply this, you need to figure out the EPSG codes of your projections.

UTM 35N probably is EPSG:32635, and your lat/lon coordinates could be anything. One of the more common ones is WGS84 as used in GPS with code EPSG:4326

The Postgis query then would be

SELECT ST_Transform(geom, 4326,32635)

If your postgis table already has the proper source SRID, you can do a simpler

SELECT ST_Transform(geom, 32635)

EDIT after you updated your question a lot:

The issue is not with ST_Transform, the issue is your query. You explicitly set all geometries to be the same. This is how UPDATE in combination with subqueries works. I really recommend you dig through the documentation a bit more to grasp this.

The solution to this issue can be found in the docs as well, by using a WHERE clause at the end. While not standard SQL, PostgreSQL allows UPDATE WHERE in order to apply subqueries to specific rows they match. In most cases, one has a unique identifier that you can use (just hand it down from the sub queries). In your case, depending on your data, you could use the timestamp or the lon/lat/alt itself to do so, as all these should match the same anyway.:

UPDATE info SET geom = sq.geom 
FROM (
    SELECT lon, lat,alt, ST_Transform(geom,32635) as geom
    FROM (
        SELECT lon,lat,alt, ST_setsrid(geom_v,4326) as geom 
        from (
            select lon,lat,alt, ST_MakePoint(lon,lat,alt) as geom_v 
            from (
                select latitude as lat, longitude as lon, altitude as alt from info order by gid asc
            ) as fpp
        )as ftt
    ) as tr
) as sq
WHERE info.lon = sq.lon AND info.lat = sq.lat AND info.alt = sq.alt

However, your entire query is unnecessarily convoluted (extremely so!). You could (and should) just use this:

UPDATE test SET geom = ST_Transform(ST_SetSRID(ST_MakePoint(lon,lat, alt),4326),32635)

to achieve the same result with a much faster and easier-to-read query.

Best Answer

Related Solutions

[GIS] ny character limit in st_union and array_agg

[GIS] How to make projected geometry and update it to postgresql

Related Question