[Gdal-dev] Help in Transforming coordinates

bargi bhagwat.maimt at gmail.com
Mon Feb 25 00:05:12 EST 2008


Hi,

I have posted earlier a message regarding Coordinate Transformation.

I have tried a lot and was able to transformed the coordinate ,but not able
to show them....

I have used Point shapefile in which there was no ".prj" file for projection
,so I create a.prj file on my own which is as follow:

point.prj:

PROJCS["World_Mercator",GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137.0,298.257223563]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Mercator"],PARAMETER["False_Easting",2600000.0],PARAMETER["False_Northing",400000.0],PARAMETER["Central_Meridian",-77.0],PARAMETER["Standard_Parallel_1",15.0],UNIT["Meter",1.0]]


By using this projection information I create my own Target source and
transform the Coordinates.....

After transforming I used SetGeometryDirectly() of OGRFeauture to set the
changed geometry......

The details of code is as follow:

[code]
  OGRSpatialReference*    pSourceSRS = layer->GetSpatialRef();
         OGRSpatialReference*    pTargetSRS = pSourceSRS->Clone(); 
         pTargetSRS->SetProjCS("UTM");
         pTargetSRS->SetTM(0,0,1,0,0);
         
          OGRGeometry*  poGeometry_new;
          OGREnvelope envelope;
           double x,y;
          

        OGRCoordinateTransformation* poCT =
OGRCreateCoordinateTransformation(pSourceSRS,pTargetSRS); 
       
         OGRFeature *ogrFeature = 0;
          layer->ResetReading();   
     for ( int i = 0; i < nLayers; ++i ) {
         layer = ogrDataSource->GetLayer( i );
       
         while ( (ogrFeature = layer->GetNextFeature()) ) {
               k++;  
               poGeometry = ogrFeature->GetGeometryRef(); 
                if (poGeometry !=NULL)  { 
       
             if (wkbFlatten(poGeometry->getGeometryType()) ==wkbPoint) {
                 OGRPoint *poPoint = (OGRPoint *) poGeometry; 
                 x = poPoint->getX();
                 y = poPoint->getY();

               if (poCT == NULL || !poCT->Transform(1,&x,&y)) 
                      printf(" TRANSFORMATION FAILED");
                else {


                     poPoint->setX(x);    //Set coordinates

                     poPoint->setY(y);

                   poGeometry_new = (OGRGeometry*)poPoint;
                  
    ogrFeature->SetGeometryDirectly(poGeometry_new);

  
    OGRGeometry *geometry = ogrFeature->GetGeometryRef(); // For getting new
geometry
    
    geometry->getEnvelope( &envelope ); // At this line give segmentation
fault
    m_bounds.setCoords( envelope.MinX, envelope.MinY, envelope.MaxX,
envelope.MaxY );
     
    
[/code]

BEFORE TRANSFORM x=-92.0333        y=27.3167 

AFTER TRANSFORM x= -2.49294e+07   y= -4.0033e+06  


I have following Question...

1) Is the Projection I have created is correct for Points ?

2) Whether the target transformation is correct ?

3) Is the step for transformation are correct ?

4) Are the resultant values correct or they are very high then expectation
(as are seen ) ?

5) How can I got the information of projection from shapefile if there is no
.prj file ?

Please correct me where I am wrong and suggest the changes to do ??

Thanks



-- 
View this message in context: http://www.nabble.com/Help-in-Transforming-coordinates-tp15674107p15674107.html
Sent from the GDAL - Dev mailing list archive at Nabble.com.



More information about the gdal-dev mailing list