Skip to content

reading select values from a binary file

18 messages · Rick Reeves, Barry Rowlingson, Karl Ove Hufthammer +6 more

#
On Wed, Dec 9, 2009 at 4:58 PM, Jonathan Thayn <jthayn at ilstu.edu> wrote:
You need to seek help. Or help(seek). Yes really.

Barry
#
Hello All:
Faced with a similar challenge, and NOT wanting to resort to writing a C 
language function
employing fseek(), I just used two readBin calls: One to read (and 
implicitly discard) data up
to the spot I actually wanted to read, and a second call to read the 
desired information.
here is a code sample:

                      fileCon = file(sTsFileName,"rb")
#
# Offset has been calculated above: number of NumericByte-sized elements 
preceeding
# the data of interest.
#
                      if (iOffset > 0)
                      
{                                                                                          

                         TestVec = readBin(fileCon,numeric(),size = 
iNumericBytes,n=iOffset)
                      
}                                                        
                      dTimeSeries = readBin(fileCon,numeric(),size = 
iNumericBytes,n=iNumElements) 
                      close(fileCon)

Crude, but effective. I WIS that we could add a seekBin() function that 
positions the file pointer
to the desired spot.

Hope this helps,
Rick R
Barry Rowlingson wrote:

  
    
#
On Wed, Dec 9, 2009 at 6:20 PM, rick reeves <reeves at nceas.ucsb.edu> wrote:
The R seek function does this. I hadn't tested it when I posted,
hoping the original poster would work it from my message. So I just
now tested it:

 I created a file with ascii a to z in, and then:

 > con = file("file.tst","rb")
 # jump to 13th letter:
 > seek(con,13)
 [1] 0
 # read it:
 > readBin(con,what="raw",n=1)
 [1] 6e     # thats ascii 'n'
 # jump back:
 > seek(con,1)
 [1] 14 # returns where we were, we've now moved...
 > readBin(con,what="raw",n=1)
 [1] 62  # ascii 'b'. Things start at zero....

 Is that what you are trying to do?

Barry
#
If you can read the whole file at once it might be better to push the
results into a matrix for indexing in R. This could also be useful for
reading chunks of the file, and processing in parts.  seek is fine but
if you are moving around for three values each read it will be slow.
If the arrangement complex,  R's "[" indexing can save you from making
simple errors.

BTW, what is the format of the file? If GDAL can read it the rgdal
package provides a simpler route.

Cheers, Mike.

On Thu, Dec 10, 2009 at 5:33 AM, Barry Rowlingson
<b.rowlingson at lancaster.ac.uk> wrote:
#
On Wed, 9 Dec 2009 10:58:01 -0600 Jonathan Thayn <jthayn at ilstu.edu> 
wrote:
You can use the 'readBinFragments' function in the 'R.utils' package.

Also, if your file is an a 'standard' format, the 'raster' package may 
be of use: http://raster.r-forge.r-project.org/
6 days later
#
I tried using the seek() function, but it did not work as I thought it would. Perhaps I have done something wrong:
[1] 1
[1] 2
[1] 1792 #this is not correct
[1] 15
[1] 1 #the seek(cc,0) seems to have reset it
[1] 2
[1] 4 #this is not the sixth number
[1] 5
[1] 10
[1] 7 #this is not the twelfth number
[1] 14
[1] 256
[1] 3
[1] 1
[1] 2
[1] 10


I tried the readBinFragments() function in the R.utils package and it works, so long as I use a seek(cc,0) between calls:

readBinFragments(cc,"integer",idxs=c(2,6,7,13,19),size=2,endian="big")
[1]  2  6  7 13 19
[1] 0 0 0
[1] 0 0 0
[1] 52
[1]  4  9 12


Am I using seek() incorrectly?



Jonathan B. Thayn
On Dec 9, 2009, at 12:33 PM, Barry Rowlingson wrote:

            
#
?seek

The "where" argument is an offset in bytes, so you need to offet 2
bytes for each of your integers:

writeBin(1:20,"test.bip",size=2,endian="big")
cc <- file("test.bip","rb")
readBin(cc,"integer",n=1,size=2,endian="big")
seek(cc, 26)
readBin(cc,"integer",n=1,size=2,endian="big")
## [1] 14
close(cc)

Use seek(cc) to find the current offset on the connection.

HTH
On Thu, Dec 17, 2009 at 5:12 AM, Jonathan Thayn <jthayn at ilstu.edu> wrote:
#
Hi,

I am trying to plot a shape file  in R. I can do  it easily using plot 
or spplot function.  But I want plot this map with  lattice xyplot 
function. Any one has any idea? I do not want use "maps" package.

Help will be appreciated...

Thanks

Zia

library(maptools)
library(lattice)
bound<-read.shape("Boundary.shp")
#OR
bound<-readShapePoly("Boundary.shp")
plot(bound)
spplot(bound)

# Or following way

bound<-Map2poly(bound)
xybd<-attr(bound,"maplim")

plot(xybd$x,xybd$y, type="n",asp=1, xlab="", ylab="")
for (i in 1:length(bound)) lines(bound[[i]]) # How do I apply this code 
in following  panel.function of  xyplot?


windows(width=5, height=5)
xyplot(y~x, as.data.frame(tala),pch=1, cex=.8,col="black",
       xlab=list("Easting (m)",cex=.9),ylab=list("Northing (m)",cex=.9),
       aspect = "iso",
       panel = function(...) {
              panel.xyplot(...)
              panel.abline(h = 0:4*5000 + 545000, v= 0:4*5000 + 2650000,
              col = "light grey")
              panel.points(coordinates(tala.nr),cex=.95,pch=20, col="black")
              panel.text(2671000, 563000, cex=1, " (a)")
              },
)
#
What's wrong with using spplot? It is based on xyplot anyway.

If you want to use xyplot, you can use the panel function that sp defines:
sp.polygons(shp)

You can read in the shapefile with the RGDAL package:
shp <- readOGR(...)


2009/12/17 Zia Ahmed <zua3 at cornell.edu>:
-- 
Felix Andrews / ???
Postdoctoral Fellow
Integrated Catchment Assessment and Management (iCAM) Centre
Fenner School of Environment and Society [Bldg 48a]
The Australian National University
Canberra ACT 0200 Australia
M: +61 410 400 963
T: + 61 2 6125 4670
E: felix.andrews at anu.edu.au
CRICOS Provider No. 00120C
#
Hello all,

I have a data frame with lat and lon coordinates for some sample areas 
(not a grid) .

I was to create an nb object that I can use with spdep but am unsure as 
to how to convert it.

Does anybody have any ideas?


Thanks,

Pete
#
Thanks Andrews. It works fine. Nothing wrong with sp plot.  I need to 
plot boundary around of sampling locations. R-script I have  written, I 
used xyplot and levelplot.

Thanks  again

Zia
Felix Andrews wrote:
#
Hello all,

I have a data frame with lat and lon coordinates for some sample areas 
(not a grid) .

I was to create an nb object that I can use with spdep but am unsure as 
to how to convert it.

Does anybody have any ideas?

Thanks,


Pete
#
On Wed, 16 Dec 2009, Zia Ahmed wrote:

            
This is no longer possible in current maptools - read.shape() is no longer 
in the package namespace.
Nor is this possible now - after years of deprecation, these S3 classes 
are now defunct and their methods not in the namespace of maptools. Only 
S4 sp classes are available. Retaining the antiquated older versions, and 
trying to help users who instist on not reading warnings about 
deprecation, has ceased to be worthwhile.

Everyone looking for advise in the archives will see messages mentioning 
read.shape() etc., this statement should suffice to indicate what is 
going on.

Roger

  
    
#
On Wed, 16 Dec 2009, Zia Ahmed wrote:

            
This look like a two-step question, assuming that you need a line 
representation of possible polygon boundaries from the shapefile. If it is 
a SpatialLines object, omit the first step (it will work on 
SpatialPolygons* objects with appropriate slot names, but the intention 
is to use lines, so make sure that they are):

library(sp)
library(maptools)
xx <- readShapeSpatial(system.file("shapes/sids.shp",
   package="maptools")[1])
proj4string(xx) <- CRS("+proj=longlat +ellps=clrk66")
# assigning geographical coordinates for this data set
class(xx)

lxx <- as(xx, "SpatialLines")
class(lxx)

lns <- slot(lxx, "lines")
olns <- lapply(lns, slot, "Lines")
ocrds <- matrix(nrow=0, ncol=2)
for (i in seq(along=olns)) {
   for (j in seq(along=olns[[i]])) {
     crds <- rbind(slot(olns[[i]][[j]], "coords"), c(NA, NA))
     ocrds <- rbind(ocrds, crds)
   }
}

plot(ocrds, type="l")
# sanity check

pts <- coordinates(spsample(xx, n=400, type="random"))

library(lattice)
xyplot(y~x, as.data.frame(pts), pch=1, cex=.8, col="black",
   xlab=list("Easting (m)",cex=.9), ylab=list("Northing (m)",cex=.9),
   aspect = mapasp(xx), # setting aspect for geographical coordinates
   panel = function(...) {
     panel.xyplot(...)
     panel.lines(ocrds)
   }
)

In summary: make an S-style coordinates matrix (lines separated by NA 
rows), then pass through panel.lines().

Hope this helps,

Roger

  
    
#
On Wed, 16 Dec 2009, Pete Larson wrote:

            
Create a matrix of coordinates, if not using sp classes, for example:

xy <- cbind(df$x, df$y)

or if your data frame is a SpatialPointsDataFrame, just:

xy <- coordinates(obj)

Then your choice of dnearneigh(), knn2nb(knearneigh()), graph2nb() with 
three graph variants, tri2nb() for a Delaunay triangulation. The first two 
of these also take a longlat= argument in case your coordinates are 
geographical, not projected. Triangulation and graph-based methods need 
planar coordinates, so you would have to project them first.

Hope this helps,

Roger

  
    
#
Roger,

OK, I did what you indicated here:

xy<-coordinates(fish)
test<-knn2nb(knearneigh(xy))

with a summary:
summary(test)
Neighbour list object:
Number of regions: 341
Number of nonzero links: 341
Percentage nonzero weights: 0.2932551
Average number of links: 1
Non-symmetric neighbours list

Then I want to run a spatial regression like this:

spautolm(N.fish.spp ~ TEMPannual + Elevation, data=fish, 
family="SAR",listw=test, method="Matrix")

and I get this error:

Error in spautolm(N.fish.spp ~ TEMPannual + Elevation, data = fish, 
family = "SAR",  :
  No neighbourhood list

I thought that the "test" object WAS a neighborhood list?

Any help would be appreciated.

Pete
Roger Bivand wrote:
#
On Thu, 17 Dec 2009, Pete Larson wrote:

            
class(test)
?nb2listw

Neighbourhood "nb" just says that there is a link, not the weight to put 
on it. Use nb2listw() to make the weights "listw" object, after 
considering the options available. If the actual process doesn't match 
your guess, you may not be able to model it adequately. Is the choice of 
first nearest neighbour considered, as you see, the neighbour 
relationships are not symmetric? Are your points in geographical or 
projected coordinates - it does make a difference? Marine fish data are 
often in geographical coordinates (so use longlat=TRUE), but yours are 
fresh-water?

Roger