t.create--> t.register--> t.vect.observer.strds in GRASS7 ? Looks really nifty, could be useful with the Landsat Cube data sets <a href="view-source:http://landsat.usgs.gov/documents/Oct27_29_2009_huang_LST_boston.ppt" class="attribute-value">http://landsat.usgs.gov/documents/Oct27_29_2009_huang_LST_boston.ppt</a><div>
<br></div><div><br><div><br></div><div>Doug</div><div><br><br><div class="gmail_quote">On Fri, Feb 22, 2013 at 2:53 PM, Michael Barton <span dir="ltr"><<a href="mailto:Michael.Barton@asu.edu" target="_blank">Michael.Barton@asu.edu</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div style="word-wrap:break-word">
Thanks Doug,
<div><br>
</div>
<div>I knew I could do it with a script and v.what.rast. </div>
<div><br>
</div>
<div>What I was hoping was that there is a shortcut already usable in GRASS modules. Looks like the new temporal GIS tools may be able to do it.</div>
<div><br>
</div>
<div><div class="im">Michael<br>
<div><span style="text-indent:0px;letter-spacing:normal;font-variant:normal;text-align:auto;font-style:normal;font-weight:normal;line-height:normal;border-collapse:separate;text-transform:none;font-size:medium;white-space:normal;font-family:'Lucida Grande';word-spacing:0px"><span style="text-indent:0px;letter-spacing:normal;font-variant:normal;font-style:normal;font-weight:normal;line-height:normal;border-collapse:separate;text-transform:none;font-size:medium;white-space:normal;font-family:'Lucida Grande';word-spacing:0px">
<div style="word-wrap:break-word">
<span style="text-indent:0px;letter-spacing:normal;font-variant:normal;font-style:normal;font-weight:normal;line-height:normal;border-collapse:separate;text-transform:none;font-size:medium;white-space:normal;font-family:'Lucida Grande';word-spacing:0px">
<div style="word-wrap:break-word">
<span style="text-indent:0px;letter-spacing:normal;font-variant:normal;font-style:normal;font-weight:normal;line-height:normal;border-collapse:separate;text-transform:none;font-size:medium;white-space:normal;font-family:'Lucida Grande';word-spacing:0px">
<div style="word-wrap:break-word">
<span style="text-indent:0px;letter-spacing:normal;font-variant:normal;font-style:normal;font-weight:normal;line-height:normal;border-collapse:separate;text-transform:none;font-size:medium;white-space:normal;font-family:'Lucida Grande';word-spacing:0px">
<div style="word-wrap:break-word">
<span style="text-indent:0px;letter-spacing:normal;font-variant:normal;font-style:normal;font-weight:normal;line-height:normal;border-collapse:separate;text-transform:none;font-size:medium;white-space:normal;font-family:'Lucida Grande';word-spacing:0px">
<div style="word-wrap:break-word">
<span style="text-indent:0px;letter-spacing:normal;font-variant:normal;font-style:normal;font-weight:normal;line-height:normal;border-collapse:separate;text-transform:none;font-size:medium;white-space:normal;font-family:'Lucida Grande';word-spacing:0px">
<div style="word-wrap:break-word">
<span style="text-indent:0px;letter-spacing:normal;font-variant:normal;font-style:normal;font-weight:normal;line-height:normal;border-collapse:separate;text-transform:none;font-size:13px;white-space:normal;font-family:'Lucida Grande';word-spacing:0px">
<div style="word-wrap:break-word">
<span style="text-indent:0px;letter-spacing:normal;font-variant:normal;font-style:normal;font-weight:normal;line-height:normal;border-collapse:separate;text-transform:none;font-size:13px;white-space:normal;font-family:'Lucida Grande';word-spacing:0px">
<div style="word-wrap:break-word">
<span style="text-indent:0px;letter-spacing:normal;font-variant:normal;font-style:normal;font-weight:normal;line-height:normal;border-collapse:separate;text-transform:none;font-size:13px;white-space:normal;font-family:'Lucida Grande';word-spacing:0px">
<div style="word-wrap:break-word">
<div><font face="Helvetica" size="3"><span style="font-size:12px">______________________________</span></font></div>
<div><font face="Helvetica" size="3"><span style="font-size:12px">C. Michael Barton </span></font></div>
<div><font face="Helvetica" size="3"><span style="font-size:12px">Director, Center for Social Dynamics & Complexity</span></font></div>
<div><font face="Helvetica" size="3"><span style="font-size:12px">Professor of Anthropology, School of Human Evolution & Social Change</span></font></div>
<div><span style="font-size:12px;font-family:Helvetica">Arizona State University</span></div>
<div><font face="Helvetica" size="3"><span style="font-size:12px">Tempe, AZ 85287-2402</span></font></div>
<div><font face="Helvetica" size="3"><span style="font-size:12px">USA</span></font></div>
<div><font face="Helvetica" size="3"><span style="font-size:12px"><br>
</span></font></div>
<div><font face="Helvetica" size="3"><span style="font-size:12px">voice:<span> </span><span style="white-space:pre-wrap">
</span>480-965-6262 (SHESC), 480-727-9746 (CSDC)</span></font></div>
<div><font face="Helvetica" size="3"><span style="font-size:12px">fax: 480-965-7671(SHESC), 480-727-0709 (CSDC)</span></font></div>
<div><font face="Helvetica" size="3"><span style="font-size:12px">www: <span style="white-space:pre-wrap">
</span><a href="http://csdc.asu.edu" target="_blank">http://csdc.asu.edu</a>, <a href="http://shesc.asu.edu" target="_blank">
http://shesc.asu.edu</a></span></font></div>
<div><font face="Helvetica" size="3"><span style="font-size:12px"><span style="white-space:pre-wrap"></span><a href="http://www.public.asu.edu/~cmbarton" target="_blank">http://www.public.asu.edu/~cmbarton</a></span></font></div>
</div>
</span></div>
</span></div>
</span></div>
</span></div>
</span></div>
</span></div>
</span></div>
</span></div>
</span></span></div>
<br>
</div><div>
<div>On Feb 22, 2013, at 12:31 PM, "Newcomb, Doug" <<a href="mailto:doug_newcomb@fws.gov" target="_blank">doug_newcomb@fws.gov</a>></div>
<div> wrote:</div><div><div class="h5">
<br>
<blockquote type="cite">
<div>Michael,</div>
<div>You could use v.what.rast in a python script , iterating the raster layers with an sqlite database back end. I think that would go up to 2000 columns for each point you have.</div>
<div><br>
</div>
<div>Alternatively, you could use a bit of python with gdal. I was trying to do something similar in GRASS, to change the Z value of each point in a text LiDAR file, from absolute above sea level to relative the elevation of the ground. For 25.5 billion text
points and a Statewide 20 ft elevation grid (6.5 ? billion cells) , it was a bit slow using r.what. So I converted the LiDAR data to (7) 3.3 billion point LAS format files and exported the GRASS layer to an Erdas imagine format file and wrote the following
ugly python script: </div>
<div><br>
</div>
<div>#!/usr/bin/python</div>
<div>import os,string,glob,re,gdal</div>
<div>from liblas import file</div>
<div>from liblas import header</div>
<div>from liblas import point</div>
<div>from gdalconst import *</div>
<div>h=header.Header()</div>
<div>#enter the LAS point file name</div>
<div>infile=raw_input("Enter the input lidar data points file: ")</div>
<div>#Hardcoded edras imagine file, you will have to use an array for the different data layer names</div>
<div>imgfile="/gisdata2/raster/allnc_20ft_el.img"</div>
<div>#print "suggest /gisdata2/raster for output dir\n"</div>
<div>inarr=infile.split('.')</div>
<div>#This sets the LAS output file, substitute your output text file instead</div>
<div>outfil=inarr[0]+"_norm.las"</div>
<div>#outfil=raw_input("Enter output text file name: ")</div>
<div>#This part reads the LAS file, if you </div>
<div>l=file.File(infile,mode='r')</div>
<div>#Outputs the LAS file</div>
<div>lout=file.File(outfil,mode='w',header=h)</div>
<div># register all of the drivers, hopefully your gdal speaks GRASS</div>
<div>gdal.AllRegister()</div>
<div>#opening and closing the image layers might take some time if you are reading thousands of images</div>
<div>ds=gdal.Open(imgfile,GA_ReadOnly)</div>
<div>if ds is None:</div>
<div> print 'Could not open image'</div>
<div> sys.exit(1)</div>
<div># get image size</div>
<div>rows = ds.RasterYSize</div>
<div>cols = ds.RasterXSize</div>
<div>bands = ds.RasterCount</div>
<div># get georeference info, not sure how this would work for GRASS data layers</div>
<div>transform = ds.GetGeoTransform()</div>
<div>xOrigin = transform[0]</div>
<div>yOrigin = transformAsArray(xOffset, yOffset, 1, 1)</div>
<div>pixelWidth = transform[1]</div>
<div>pixelHeight = transform[5]</div>
<div>
<div>for p in l:</div>
<div> x=float(p.x)</div>
<div> y=float(p.y)</div>
<div> z=float(p.z)</div>
<div> # compute pixel offset</div>
<div> xOffset = int((x - xOrigin) / pixelWidth)</div>
<div> yOffset = int((y - yOrigin) / pixelHeight)</div>
<div> band = ds.GetRasterBand(1) # 1-based index 0? 1?</div>
<div> data = band.Readr(value) :continue</div>
<div> value = data[0,0]</div>
<div> #print value,"11","\n"</div>
<div> if "nan" in st[3]</div>
<div> znorm = z-value</div>
<div> #print znorm,"\n"</div>
<div> pt=point.Point()</div>
<div> pt.x=p.x</div>
<div> pt.y=p.y</div>
<div> pt.z=znorm</div>
<div> lout.write(pt)</div>
<div><br>
</div>
<div>l.close()</div>
<div>lout.close()</div>
<div>#25561312019 points in allreturns</div>
</div>
<div><br>
</div>
<div>This managed to process everything over a weekend( in 7 parallel threads), which was fast enough for me at the time. Approaching your problem , if your data layers are all n GRASS and your version of gdal can read GRASS data layers, I would grab the
list of GRASS layers via glob and iterate the layers , writing the name of the raster layer and the value of the raster layer at the point coordinates out to a text file as you state below. </div>
<div><br>
</div>
<div>Hope this helps,</div>
<div>Doug</div>
<div> </div>
<br>
<div class="gmail_quote">On Fri, Feb 22, 2013 at 1:23 PM, Michael Barton <span dir="ltr">
<<a href="mailto:Michael.Barton@asu.edu" target="_blank">Michael.Barton@asu.edu</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
But I want to do it with a time series of hundreds or thousands of maps.<br>
<div><br>
Michael<br>
______________________________<br>
C. Michael Barton<br>
Director, Center for Social Dynamics & Complexity<br>
Professor of Anthropology, School of Human Evolution & Social Change<br>
Arizona State University<br>
Tempe, AZ 85287-2402<br>
USA<br>
<br>
voice: 480-965-6262 (SHESC), 480-727-9746 (CSDC)<br>
fax: 480-965-7671(SHESC), 480-727-0709 (CSDC)<br>
www: <a href="http://csdc.asu.edu/" target="_blank">http://csdc.asu.edu</a>, <a href="http://shesc.asu.edu/" target="_blank">
http://shesc.asu.edu</a><br>
<a href="http://www.public.asu.edu/~cmbarton" target="_blank">http://www.public.asu.edu/~cmbarton</a><br>
<br>
</div>
On Feb 22, 2013, at 10:53 AM, Markus Neteler <<a href="mailto:neteler@osgeo.org" target="_blank">neteler@osgeo.org</a>><br>
<div>
<div> wrote:<br>
<br>
> On Fri, Feb 22, 2013 at 6:41 PM, Michael Barton <<a href="mailto:Michael.Barton@asu.edu" target="_blank">Michael.Barton@asu.edu</a>> wrote:<br>
>> Is there tool somewhere, including in the new temporal GIS modules, to<br>
>> sample the value of a raster series at one cell location? I'd like to get<br>
>> text something like this for a cell that I specify with xy coordinates or a<br>
>> cat:<br>
>><br>
>> map1, value 1<br>
>> map2, value 2<br>
>> map3, value 3<br>
>> map4, value 4<br>
>> map5, value 5<br>
>> …<br>
>><br>
><br>
> I do such queries with r.what which accepts multiple input:<br>
><br>
> <a href="http://grass.osgeo.org/grass64/manuals/r.what.html" target="_blank">http://grass.osgeo.org/grass64/manuals/r.what.html</a><br>
><br>
> Markus<br>
<br>
_______________________________________________<br>
grass-dev mailing list<br>
<a href="mailto:grass-dev@lists.osgeo.org" target="_blank">grass-dev@lists.osgeo.org</a><br>
<a href="http://lists.osgeo.org/mailman/listinfo/grass-dev" target="_blank">http://lists.osgeo.org/mailman/listinfo/grass-dev</a><br>
</div>
</div>
</blockquote>
</div>
<br>
<br clear="all">
<div><br>
</div>
-- <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>
</blockquote>
</div></div></div>
<br>
</div>
</div>
</blockquote></div><br><br clear="all"><div><br></div>-- <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>
</div></div>