I'm new to GDAL and am usingthe python bindings. I want to transform a 
newly created shapefile from WGS84 to a Lambert Conical Conformal so 
that I can calculate the length of the contained polylines. The script 
exits with no errors but when I look at the shapefile in arcCatalog I 
see that the projection is still WGS84. Any ideas? Also, is there a way 
to calculate the length directly from GDAL? Thanks.<br>
<br>import os<br>try:<br>    from osgeo import ogr<br>    from osgeo 
import osr<br><br>except:<br>    import ogr<br>    import osr<br><br>#~ 
start driver<br>driver = ogr.GetDriverByName(&#39;ESRI Shapefile&#39;)<br>#~ set
 workspace<br>
ws = &#39;C:\\Users\\deadpickle\\<div>Desktop\\case study\\verify\\&#39;<br>os.chdir(ws)<br><br>#~

 check if it exists<br>if os.path.exists(&#39;test.shp&#39;):<br>    
driver.DeleteDataSource(&#39;test.shp&#39;)<br><br>#~ create a new shp<br>
ds = driver.CreateDataSource(&#39;test.shp&#39;)<br><br>#~ create spatref 
for new shp<br>source = osr.SpatialReference()<br>source.SetFromUserInput(&#39;WGS84&#39;)<br><br>#~

 create a new data layer<br>layer = ds.CreateLayer(&#39;test&#39;,geom_type=ogr.wkbLineString,srs=source)<br>
<br>#~ add an id field<br>fieldDefn = ogr.FieldDefn(&#39;id&#39;, 
ogr.OFTInteger)<br>layer.CreateField(fieldDefn)<br><br>#~ create a 
geometry<br>line = ogr.Geometry(ogr.wkbLineString)<br><br>#~ get 
the def of the layer, scheme<br>
featureDefn = layer.GetLayerDefn()<br><br># create a new feature<br>feature

 = ogr.Feature(featureDefn)<br><br>#~ Read 2006track file<br>file = 
open(&quot;2006track.csv&quot;,&quot;r&quot;)<br>file.readline()<br>lnlist = 
file.readlines()<br>
for lines in lnlist:<br>    value = lines.split(&quot;,&quot;)<br>    if 
float(value[1]) == 2:<br>        print value<br>        
line.AddPoint(float(value[3]),float(value[2]))<br>        #~ set 
the field &#39;ID&#39;<br>        feature.SetField(&#39;id&#39;, float(value[1]))<br>
        <br>#~ set the geometry for this shp<br>feature.SetGeometry(line)<br><br>#
 add the feature to the output layer<br>layer.CreateFeature(feature)<br><br>#~

 feature = layer.GetFeature(0)<br>geom = feature.GetGeometryRef()<br>
<br>#~ change projection to lcc<br>target = osr.SpatialReference()<br>target.ImportFromProj4(&#39;+proj=lcc

 +lat_1= 33.000000 +lat_2=45.000000 +lat_0=39.000000 +lon_0=-96.000000 
+x_0=0.0 +y_0=0.0&#39;)<br>coordtrans = osr.CoordinateTransformation(source,

 target)<br>
<br>geom.Transform(coordtrans)<br><br># destroy the geometry and feature
 and close the data source<br><br>ds.Destroy()<br>line.Destroy()<br>feature.Destroy()<br>file.close()</div>