<div dir="ltr"><div>Hello,</div><div><br></div><div>I am fairly new to PostGIS and I am hoping some one can provide me a hand with an issue relating to scenes crossing the International Date Line (IDL) and finding the center point. </div><div><br></div><div>I currently maintain a PostgreSQL database that contains footprints for satellite scenes. The scenes have global coverage and range in size from ~1500Km squared down to ~100Km squared. The footprints are stored in a geography column as basic polygons as seen below. </div><div><br></div><div><font face="monospace, monospace">SRID=4326;POLYGON((179.255436389612 49.0944921724897,-179.572711844818 49.0656217040858,-179.619180111394 48.3778311277411,179.224815891643 48.4060150740941,179.255436389612 49.0944921724897))</font><br></div><div><br></div><div>I need to be able to determine the center point of each scene in order to provide some additional information to various applications using the data. Based on forums and various other web resources, it appears the way to calculate the center point is by converting geography to geometry and using the ST_CENTROID function as seen below.</div><div><br></div><div><font face="monospace, monospace">ST_ASEWKT(ST_CENTROID(ST_GEOMFROMEWKT(ST_ASEWKT(footprint))))</font><br></div><div><br></div><div>Using my polygon example from above, the function call would look like the following.</div><div><br></div><div><font face="monospace, monospace">ST_ASEWKT(ST_CENTROID(ST_GEOMFROMEWKT('SRID=4326;POLYGON((179.255436389612 49.0944921724897,-179.572711844818 49.0656217040858,-179.619180111394 48.3778311277411,179.224815891643 48.4060150740941,179.255436389612 49.0944921724897))')))</font><br></div><div><br></div><div>Here is where the problem lies. In the above example, the centroid returns a center point value of the following:</div><div><br></div><div><div><font face="monospace, monospace"> st_asewkt </font></div><div><font face="monospace, monospace"> ---------------------------------------------------- </font></div><div><font face="monospace, monospace"> SRID=4326;POINT(-0.148022213748242 48.7359898569253)</font></div></div><div><br></div><div>Note that the centroid longitude is -0.148022213748242. The footprint I provided crosses the IDL. The scene center longitude should be ~179.8181. </div><div><br></div><div>I was able to hack together an inelegant solution that will correct for the IDL on a Cartesian system up to a point. I have provided the code below. The code appears to work but will breakdown if the scene footprint is to large and I am sure I have not accounted for everything.</div><div><br></div><div><font face="monospace, monospace">CREATE FUNCTION test (in i_wkt varchar)<br>RETURNS TABLE</font></div><div><font face="monospace, monospace">( centerlat FLOAT</font></div><div><font face="monospace, monospace">, centerlon FLOAT</font></div><div><font face="monospace, monospace">);</font></div><div><font face="monospace, monospace">AS</font></div><div><font face="monospace, monospace">$$</font></div><div><font face="monospace, monospace">DECLARE </font></div><div><font face="monospace, monospace"> l_ctrlat FLOAT;</font></div><div><font face="monospace, monospace"> l_ctrlon FLOAT;</font></div><div><font face="monospace, monospace"> l_lonmin FLOAT;</font></div><div><font face="monospace, monospace"> l_lonmax FLOAT;</font></div><div><font face="monospace, monospace"> l_lontmp FLOAT;</font></div><div><font face="monospace, monospace"><br></font></div><div><font face="monospace, monospace">BEGIN</font></div><div><font face="monospace, monospace"> SELECT ST_Y(ST_CENTROID(ST_GEOMFROMEWKT(i_wkt)))</font></div><div><font face="monospace, monospace"> , ST_X(ST_CENTROID(ST_GEOMFROMEWKT(i_wkt)))</font></div><div><font face="monospace, monospace"> , ST_XMIN(ST_GEOMFROMEWKT(i_wkt))</font></div><div><font face="monospace, monospace"> , ST_XMAX(ST_GEOMFROMEWKT(i_wkt))</font></div><div><font face="monospace, monospace"> INTO l_ctrlat</font></div><div><font face="monospace, monospace"> , l_ctrlon</font></div><div><font face="monospace, monospace"> , l_lonmin</font></div><div><font face="monospace, monospace"> , l_lonmax;</font></div><div><font face="monospace, monospace"><br></font></div><div><font face="monospace, monospace"> -- See if the scene crossesthe IDL if so correct the negative value</font></div><div><font face="monospace, monospace"> IF (abs(l_lonmin) + abs(l_lonmax)) > 180 THEN</font></div><div><font face="monospace, monospace"> l_lontmp := 360 + l_lonmin;</font></div><div><font face="monospace, monospace"> </font><span style="font-family:monospace,monospace">l_ctrlon</span><font face="monospace, monospace"> := ((l_lonmax + l_lontmp)/2);</font></div><div><font face="monospace, monospace"> END IF;</font></div><div><font face="monospace, monospace"><br></font></div><div><font face="monospace, monospace"> -- Determine which side of the IDL the longitude belongs on</font></div><div><font face="monospace, monospace"> IF </font><span style="font-family:monospace,monospace">l_ctrlon</span><font face="monospace, monospace"> > 180 THEN</font></div><div><font face="monospace, monospace"> </font><span style="font-family:monospace,monospace">l_ctrlon</span><font face="monospace, monospace"> := l_lonctr - 360;</font></div><div><font face="monospace, monospace"> END IF;</font></div><div><font face="monospace, monospace"><br></font></div><div><font face="monospace, monospace"> RETURN QUERY</font></div><div><font face="monospace, monospace"> SELECT l_ctrlat AS centerlat</font></div><div><font face="monospace, monospace"> , </font><span style="font-family:monospace,monospace">l_ctrlon AS </span><font face="monospace, monospace">centerlon;</font><span style="font-family:monospace,monospace"> </span></div><div><font face="monospace, monospace"><br></font></div><div><font face="monospace, monospace">END;</font></div><div><font face="monospace, monospace">$$</font></div><div><font face="monospace, monospace">LANGUAGE plpgsql;</font></div><div><br></div><div><br></div><div>Is there a way to determine the center point for a geographic scene that crosses the IDL using current PostGIS functions that will remain true for any foot print size that does not require additional code? If so, what combination of functions do I need to use? </div><div><br></div><div>Do I need to re-project the scene prior to determining the center point with ST_CENTROID since I move the footprint to a Cartesian geometry type?</div><div><br></div><div>If no functions are available, is there a preferred method for determining center point of a polygon that crosses the IDL that will better allow for larger sizes of polygons in a geographic system? </div><div><br></div><div>I am currently running PostgreSQL 9.3.4 with PostGIS extension 2.1.0. <br></div><div><br></div><div>Any help you could provide would be appreciated.</div><div><br></div><div><br></div><div><br></div><br clear="all"><div><div class="gmail_signature"><div dir="ltr">Bradley J. Hazuka<br>Database Administrator<br>Innovate!, Inc.<br>Contractor to the U.S. Geological Survey (USGS)<br></div></div></div>
</div>