[geos-devel] Re: Computing Radials inside a Polygon
Jo
doublebyte at gmail.com
Sat Aug 29 17:53:41 EDT 2009
Thinking better, I could see what I was doing wrong.
I was rotating a point around (0,0), rather than rotating it around my
origin (the centroid). Therefore, I had to add the difference between
the point and centroid coordinates to the sin and cosine terms:
double rx=centroid->getX()+cos(fangle)*(x-centroid->getX())-sin(fangle)*(y-centroid->getY());
double ry=centroid->getY()+sin(fangle)*(x-centroid->getX())+cos(fangle)*(y-centroid->getY());
And calling consecutively this function, gives me the equally spaced
radials from a centroid of a polygon to its external ring; in this
example, using a 45 degrees angle:
http://ladybug.no-ip.org/files/radials_calc.png
cheers,
Jo
2009/8/29 <geos-devel-request at lists.osgeo.org>:
> Send geos-devel mailing list submissions to
> geos-devel at lists.osgeo.org
>
> To subscribe or unsubscribe via the World Wide Web, visit
> http://lists.osgeo.org/mailman/listinfo/geos-devel
> or, via email, send a message with subject or body 'help' to
> geos-devel-request at lists.osgeo.org
>
> You can reach the person managing the list at
> geos-devel-owner at lists.osgeo.org
>
> When replying, please edit your Subject line so it is more specific
> than "Re: Contents of geos-devel digest..."
>
>
> Today's Topics:
>
> 1. Re: Computing Radials inside a Polygon (Jo)
>
>
> ----------------------------------------------------------------------
>
> Message: 1
> Date: Sat, 29 Aug 2009 16:16:23 +0100
> From: Jo <doublebyte at gmail.com>
> Subject: [geos-devel] Re: Computing Radials inside a Polygon
> To: geos-devel at lists.osgeo.org
> Message-ID:
> <23ab5f0a0908290816j66928c44y58de481b26525657 at mail.gmail.com>
> Content-Type: text/plain; charset=ISO-8859-1
>
> Hi Jo,
>
> I think that each rotate angles should be calculated more simply.
> (I think rotation matrix is valid when centroid is equal to origine(0, 0).)
>
> The following code is pseudo code. (Sorry I haven't OGR environment...)
>
> ===================================
> // get minCircle radius. (If radius is not calculated yet)
> double radius=Distance(¢roid, &startPt); // I don't know that Distance
> function is implemented in OGR...
> // get start angle.
> double startAngle=atan2(&startPt.getY()-¢roid.getY(),
> &startPt.getX()-¢roid.getX());
>
> for (size_t i=0; i <n; ++i){
> OGRLineString dLine;
> dLine.addPoint(¢roid);
>
> double rx=centroid.getX() + radius*cos(startAngle+fangle*i);
> double ry=centroid.getY() + radius*sin(startAngle+fangle*i);
> :
> }
> ===================================
>
> Regards,
>
> Sanak.
>
> Hi Sanak,
> First of all, thanks for your reply!
> I think I didn't explain myself clearly though...
>
> The idea is to compute radials from the centroid of the polygon, to
> its external ring, at equal intervals. Something like this:
>
> http://ladybug.no-ip.org/files/radials_final.png
>
> I had the idea of starting with a vector that starts in the centroid
> and ends up in the external ring (at any point), and therefore its
> length equals the radius; to make things simple I choose a point that
> is in the same x-alignment as the centroid, and so the start angle in
> your code would be null. My start vector is something like this:
>
> http://ladybug.no-ip.org/files/start.png
>
> Than my idea, was to rotate this vector around the centroid, at a
> constant angle, and get the radials.
> Lets assume for simplicity that the angle is 90 degrees (PI/2 radians)
> and therefore I want to end up with four radials.
> As I compute the new locations of the rotated endpoint of the vector,
> I store this point on a stl container, so that in the end I can loop
> through all these locations and generate my vectors.
> Here is a GEOS version of my code:
>
> //////////GEOS VERSION /////
> std::vector<geos::geom::Point*> pts;
>
> geos::geom::GeometryFactory factory;
> geos::geom::Coordinate coordC(ogrCentroid.getX(),ogrCentroid.getY());
>
> //create centroid from given point
> geos::geom::Point* centroid = factory.createPoint(coordC);
>
> //create start pt as any point on the perimeter (for instance, a
> point aligned in x with the centroid)
> geos::geom::Coordinate coordS(ogrCentroid.getX()+radius,ogrCentroid.getY());
> geos::geom::Point* startPt = factory.createPoint(coordS);
>
> pts.push_back(startPt);
>
> //nr radials=4; angle=90 deg= PI/2
> double fangle=PI/2.0;
>
> for (size_t i=0; i <4; ++i){
>
> double x=pts.back()->getX();
> double y=pts.back()->getY();
>
> double rx=x*cos(fangle)-y*sin(fangle);
> double ry=x*sin(fangle)+y*cos(fangle);
>
> geos::geom::Coordinate coords(x+rx,y+ry);
> geos::geom::Point* aPoint = factory.createPoint(coords);
>
> pts.push_back(aPoint);
> }
>
> Let me mention that this code is *not* doing what I want, since:
> - I am getting the wrong angle between the start point and the first
> iteration point;
> - the angles do not appear to be 90 degrees;
> - I think also the angles are not constant;
>
> Here is a screenshot of my output:
>
> http://ladybug.no-ip.org/files/output.png
>
> So I guess there must be some assumption (or more than one, actually)
> that is wrong in my code :-)
> Again, I would be very grateful for any observations or suggestions...
> thanks in advance for your help!
>
> cheers,
> Jo
>
>
> ------------------------------
>
> _______________________________________________
> geos-devel mailing list
> geos-devel at lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/geos-devel
>
> End of geos-devel Digest, Vol 82, Issue 25
> ******************************************
>
--
"#define QUESTION ((bb) || !(bb))" (Shakespeare)
More information about the geos-devel
mailing list