<div dir="ltr"><div><div><div><div><div><div><div><div><div><div><div><div>Dear Jan and Doug,<br><br></div>Thanks for all these ideas. Finally I report back. <br></div>Finally I find a point shapefile data equivalent to the previous ploygon ones.<br>
</div>So I take the easy approach (which must be computationally unoptimized):<br></div>1. read all the attribute data to a csv or pandas dataframe object (<a href="http://pandas.pydata.org/">http://pandas.pydata.org/</a>) <br>
</div>2. search each (latitude,longtitude) pair, then assign to the column "geoindex" the tuple which is used to <br></div> index the final numpy array which is going to store the statiscal information data.<br>
</div>3. use pandas.DataFrame.groupby('geoindex').agg(np.men, np.std, np.sum) etc, to generate the <br></div> generate the statiscal information I want. <br></div>4. then create an empty numpy array which conforms to the geospatial extent of the area but with 0.5 degree<br>
</div> resolution, and assign the array element values by uisng the geoindex information in the dataframe.<br><br></div><div>if anybody is interested in this simple way, I put below a function used to extract the statistical information.<br>
</div><br>def dataframe_extract_statistic_info(df,target_field=None,groupby_field='geoindex'):<br> """ <br> Extract information for target_field in dataframe "df" grouped by <br>
groupby_field. <br> <br> Arguments: <br>
---------- <br> 1. target_field: the field for statistical information. <br> 2. groupby_field: the field used to group data. <br>
<br> Notes: <br> ------ <br>
1. extracted statstical information are: <br> sum, mean, std, 95percentile, 50percentile(median), 5percentile <br> 2. the groupby method will drop the null row in the groupby_field <br>
automatically. <br> """ <br> df_size_grp = df.groupby(groupby_field) <br>
dft = df_size_grp[target_field].agg([np.sum,np.mean,np.std]) <br> df_number = df_size_grp.size() <br> dft['number'] = df_number <br>
<br> dflist = [] <br> statnamelist = ['per95','median','per05'] <br>
perlist = [95,50,5] <br> func_dict = {} <br> for name,per in zip(statnamelist,perlist): <br>
dflist.append(df_size_grp[target_field].agg({name:lambda x: np.percentile(x.values,per)}))<br> dflist.extend([dft]) <br> dfall = pa.concat(dflist,axis=1) <br>
return dfall <br><br><br></div>cheers,<br><br>Chao<br></div><div class="gmail_extra"><br><br><div class="gmail_quote">On Thu, Mar 28, 2013 at 1:39 PM, Newcomb, Doug <span dir="ltr"><<a href="mailto:doug_newcomb@fws.gov" target="_blank">doug_newcomb@fws.gov</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">Chao,<div>Depending on how detailed you want your stats to be, it might be easier run the zonal statistics plugin in qgis or run v.rast.stats in GRASS 70. <br>
<br><br>Doug</div><div><div><div class="h5"><br><div class="gmail_quote">On Wed, Mar 27, 2013 at 6:38 PM, Jan Heckman <span dir="ltr"><<a href="mailto:jan.heckman@gmail.com" target="_blank">jan.heckman@gmail.com</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">Hello Chao,<div><br><div>I may not understand your wish quite well, but I think you need some datastructure which conforms to the elements of the grid and receives addtional information derived from a number of shapes with attributes.</div>
<div>Your output could be another shapefile, with the summed results as attributes, like the number of fires.</div><div>In the case that the only spatial relationship between shapes and gridelements is whether or not a shape falls inside, you can indeed simplify the polygon shape to pointshapes (centerpoints) which can be converted to a raster of 1's and 0's. Then you could use a dictionary to put the shapeattributes in the right place. I see a problem when several shapes fall into the same grid-element, having different attribute values.</div>
<div>So you would be better off with a way of rasterbuilding, where each attribute is summed per gridelement over all the shapes that fall inside. Such routines exist at least in several programs.</div><div>Alternatively, you can make a polygon shapefile which conforms to the rasterelements, then do a spatial join between your (centerpoint) shapefile and and the rastershape (effectively adding a field to the pointshapes which is the index into the grid), then add the several attribute data per gridelement by e.g. making a pivot table.</div>
<div>If you have to work with shapes that are larger than the gridelements, you would need another solution; there would be a choice between sampling the shapes or distributing shape attributes proportionally to calculated overlap with (polygonal) gridelements.</div>
<div><br></div><div>So, more simply still, assuming the grid is regular and it is possible to compute the gridindex from pointcoordinates:</div><div><br></div><div>Add an index field to the centerpoint shape which has the attributes;</div>
<div>calculate the gridindex from the coordinates and put this into the index field;</div><div>bring the table e.g. into excel;</div><div>do a pivot table for the attributes (count or sum).</div><div>If necessary, add the centerpoint coordinates to the resulting table;</div>
<div>convert the table into a new point shape with the desired attributes.</div><div><br></div><div>Hope I got you right and this helps,</div><div>Jan<br><br><div class="gmail_quote"><div><div>On Wed, Mar 27, 2013 at 9:52 PM, Chao YUE <span dir="ltr"><<a href="mailto:chaoyuejoy@gmail.com" target="_blank">chaoyuejoy@gmail.com</a>></span> wrote:<br>
</div></div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div><div>Dear all,<br><br>I am using the python binding of GDAL for reading polygon shapefiles.<br>
<br>I have a polygon shapefile which contains the fire polygon and burned area, with typical fires burning on 10-100ha. <br>
my model simulation result is based on 0.5degree X 0.5degree. So I need to extract from this shapefile on each gridcell<br>
of the 0.5d X 0.5 grid the information about like number of fires and mean fire size etc. <br><br>The approach I am thinking currently is to calculate the center of each of fire polygon, and check which 0.5d X 0.5d gridcell<br>
it falls in, then build a huge dictionary with key as the numpy index of the 0.5dX0.5d array, with the dictionary values as <br>the list(or array or pandas dataframe) which include the information of all the polygons falling within the the same grid cell.<br>
<br>I don't know if there is a better approach than this? I appreciate any comments.<br clear="all"><br>best regards,<br><br>Chao<span><font color="#888888"><br><br>-- <br><div>***********************************************************************************</div>
<div>Chao YUE<br>Laboratoire des Sciences du Climat et de l'Environnement (LSCE-IPSL)<br>UMR 1572 CEA-CNRS-UVSQ<br>Batiment 712 - Pe 119<br>91191 GIF Sur YVETTE Cedex</div>
<div>Tel: (33) 01 69 08 29 02; Fax:01.69.08.77.16<br></div>
<div>************************************************************************************</div>
</font></span><br></div></div>_______________________________________________<br>
gdal-dev mailing list<br>
<a href="mailto:gdal-dev@lists.osgeo.org" target="_blank">gdal-dev@lists.osgeo.org</a><br>
<a href="http://lists.osgeo.org/mailman/listinfo/gdal-dev" target="_blank">http://lists.osgeo.org/mailman/listinfo/gdal-dev</a><br></blockquote></div><br></div></div>
<br>_______________________________________________<br>
gdal-dev mailing list<br>
<a href="mailto:gdal-dev@lists.osgeo.org" target="_blank">gdal-dev@lists.osgeo.org</a><br>
<a href="http://lists.osgeo.org/mailman/listinfo/gdal-dev" target="_blank">http://lists.osgeo.org/mailman/listinfo/gdal-dev</a><br></blockquote></div><br><br clear="all"><div><br></div></div></div><span class="HOEnZb"><font color="#888888">-- <br>
<div>Doug Newcomb</div><div>
USFWS</div><div>Raleigh, NC</div><div>919-856-4520 ext. 14 <a href="mailto:doug_newcomb@fws.gov" target="_blank">doug_newcomb@fws.gov</a></div><div>---------------------------------------------------------------------------------------------------------</div>
<div>The opinions I express are my own and are not representative of the official policy of the U.S.Fish and Wildlife Service or Dept. of the Interior. Life is too short for undocumented, proprietary data formats.</div>
</font></span></div>
</blockquote></div><br><br clear="all"><br>-- <br><div>***********************************************************************************</div>
<div>Chao YUE<br>Laboratoire des Sciences du Climat et de l'Environnement (LSCE-IPSL)<br>UMR 1572 CEA-CNRS-UVSQ<br>Batiment 712 - Pe 119<br>91191 GIF Sur YVETTE Cedex</div>
<div>Tel: (33) 01 69 08 29 02; Fax:01.69.08.77.16<br></div>
<div>************************************************************************************</div>
</div>