ArcView shape files and GeoBUGS

I think this is it:
http://www.biostat.umn.edu/~brad/yuecui/index.html
http://www.biostat.umn.edu/~brad/yuecui/convert.r

#/09/23/03 The only different btw this convert.r and convert.s is that
#in convert.s, we us AsciiToInt("character")-AsciiToInt("0") to change
#character into numeric, in convert.r, we use type.convert("character")
#to do this.

#/09/23/03
#This version of convert.s corrects the problems found in convert2&3.s

#Tried every combination of options for scan, the following doesn't work.
#checking by comparing the number of observations in fortest and in .cgm file
#fortest<-scan("test.cgm",what=list(name=""),widths=c(90),multi.line=F,sep="\n")
#source("convert.r")

convert<-function(cgmfile)
  {
#outfile <- "test.txt"
 outfile<-paste(cgmfile,".txt",sep="")

#This one works
#fortest<-scan("test.cgm",what=list(name=""),sep="\n")
 fortest<-scan(paste(cgmfile,".cgm",sep=""),what=list(name=""),sep="\n")

fortest<-fortest$name
fortest<-fortest[grep("VIS",fortest)]
fortest<-as.matrix(fortest)

totpolyn<-length(fortest[grep("POLYGON_SET",fortest)])

polyn<-0
count<-0;
indicator<-0;
#First tried rep(0,1000), but there is a polygon with 167 rows, so we get
#167*3*2=1002 coordinates in one polygon. Error occurs with NA output
coord<-rep(0,5000)

write(paste("map:",totpolyn,"\n"),outfile)
for (i in 1:totpolyn) {
    write(paste( i, paste("grid",i,sep="")),outfile,append=T)
  }

for (i in 1:length(fortest)){
  for(j in (1:nchar(fortest[i]))){
      letter<-substring(fortest[i],j,j)
      #At first try to use AsciiToInt,stupid
      #if (AsciiToInt(letter)<=AsciiToInt("9")&
      #    AsciiToInt(letter)>=AsciiToInt("0"))
      if(letter<="9"&letter>="0"){
                         if (indicator==0) {count<-count+1;
                                            indicator<- 1
                                         }
                            coord[count]<-coord[count]*10+
                                 type.convert(letter)
                                 }
      #add the following else if statement because if and only when
      #a POLYGON_SET is encounted, polyn is increased by 1
      else if (letter=="P"){ polyn<-polyn+1;indicator<-0}

      #a CLOSEVIS is encounted, output the coordinates set, this is
      #done independently with increase of polyn since a POLYGON_SET may
      #consist of several small polygons
      else if(letter=="C"){
        coordmat<-cbind(rep(paste("grid",polyn,sep=""),count/2),
                        coord[2*(1:(count/2))-1],
                        coord[2*(1:(count/2))])

        if (polyn == 1) { write("",outfile,append=T)}
        else            { write(c(NA,NA,NA),outfile,append=T)}
        write(t(coordmat),outfile,append=T,ncol=3)
        #polyn<-polyn+1
        #the above statement is not right because 2 CLOSEVIS may lie in one
        #same polygon, but we output a polygon coordinates set whenever a
        #CLOSEVIS is encountered, although it may have the same label,
        #denoted by "grid&polyn" with previous or next coordinate sets.

        count<-0
        indicator<-0
        coord<-rep(0,5000)
        #if(grep("\r",substring(fortest[i],j+1,nchar(fortest[i]))))
        #break
       }

      else        indicator<-0
   }
}
write("END",outfile,append=T)
}

ArcView shape files and GeoBUGS

Thread (2 messages)