extra/shapefiles/shapefiles.factor

   1 ! Copyright (C) 2018 John Benediktsson
   2 ! See http://factorcode.org/license.txt for BSD license
   3
   4 USING: accessors byte-arrays classes combinators endian io
   5 io.encodings.binary io.files io.streams.byte-array kernel locals
   6 math math.order math.statistics sequences sequences.extras sets
   7
   8 IN: shapefiles
   9
  10 SINGLETON: null-shape
  11 TUPLE: point x y ;
  12 TUPLE: multipoint box points ;
  13 TUPLE: polygon box parts points ;
  14 TUPLE: polyline box parts points ;
  15 TUPLE: point-z < point z m ;
  16 TUPLE: polyline-z < polyline z-range z-array m-range m-array ;
  17 TUPLE: polygon-z < polygon z-range z-array m-range m-array ;
  18 TUPLE: multipoint-z < multipoint z-range z-array m-range m-array ;
  19 TUPLE: point-m < point m ;
  20 TUPLE: polyline-m < polyline m-range m-array ;
  21 TUPLE: polygon-m < polygon m-range m-array ;
  22 TUPLE: multipoint-m < multipoint m-range m-array ;
  23 TUPLE: multipatch box parts part-types points z-range z-array m-range m-array ;
  24
  25 <PRIVATE
  26
  27 : read-int ( -- n )
  28     4 read le> ;
  29
  30 : read-ints ( n -- parts )
  31     [ read-int ] replicate ;
  32
  33 : read-double ( -- n )
  34     8 read le> bits>double ;
  35
  36 : read-doubles ( n -- array )
  37     [ read-double ] replicate ;
  38
  39 : read-box ( -- box )
  40     4 read-doubles ;
  41
  42 : read-range ( -- range )
  43     2 read-doubles ;
  44
  45 : read-point ( -- point )
  46     read-double read-double point boa ;
  47
  48 : read-points ( n -- points )
  49     [ read-point ] replicate ;
  50
  51 : (read-multipoint) ( -- box points )
  52     read-box read-int read-points ;
  53
  54 : read-multipoint ( -- multipoint )
  55     (read-multipoint) multipoint boa ;
  56
  57 : read-poly ( -- box parts points )
  58     read-box read-int read-int [ read-ints ] dip read-points ;
  59
  60 : read-point-z ( -- point-z )
  61     read-double read-double read-double read-double point-z boa ;
  62
  63 : read-poly-z ( -- box parts points z-range z-array m-range m-array )
  64     read-poly read-range over length
  65     [ read-doubles read-range ] [ read-doubles ] bi ;
  66
  67 : read-multipoint-z ( -- multipoint-z )
  68     (read-multipoint) read-range over length
  69     [ read-doubles read-range ] [ read-doubles ] bi
  70     multipoint-z boa ;
  71
  72 : read-point-m ( -- point-m )
  73     read-double read-double read-double point-m boa ;
  74
  75 : read-poly-m ( -- box parts points m-range m-array )
  76     read-poly read-range over length read-doubles ;
  77
  78 : read-multipoint-m ( -- multipoint-m )
  79     (read-multipoint) read-range over length read-doubles
  80     multipoint-m boa ;
  81
  82 : read-multipatch ( -- multipatch )
  83     read-box read-int read-int
  84     [ [ read-ints ] [ read-ints ] bi ] dip
  85     [ read-points read-range ]
  86     [ read-doubles read-range ]
  87     [ read-doubles ] tri multipatch boa ;
  88
  89 : read-shape ( -- shape )
  90     read-int {
  91         { 0 [ null-shape ] }
  92         { 1 [ read-point ] }
  93         { 3 [ read-poly polyline boa ] }
  94         { 5 [ read-poly polygon boa ] }
  95         { 8 [ read-multipoint ] }
  96         { 11 [ read-point-z ] }
  97         { 13 [ read-poly-z polyline-z boa ] }
  98         { 15 [ read-poly-z polygon-z boa ] }
  99         { 18 [ read-multipoint-z ] }
 100         { 21 [ read-point-m ] }
 101         { 23 [ read-poly-m polyline-m boa ] }
 102         { 25 [ read-poly-m polygon-m boa ] }
 103         { 28 [ read-multipoint-m ] }
 104         { 31 [ read-multipatch ] }
 105     } case ;
 106
 107 TUPLE: header { file-code initial: 9994 } file-length
 108 { version initial: 1000 } shape-type x-min y-min x-max y-max
 109 z-min z-max m-min m-max ;
 110
 111 : read-header ( -- header )
 112     4 read be> dup 9994 assert=
 113     20 read drop ! unused
 114     4 read be>
 115     read-int dup 1000 assert=
 116     read-int
 117     read-double
 118     read-double
 119     read-double
 120     read-double
 121     read-double
 122     read-double
 123     read-double
 124     read-double
 125     header boa ;
 126
 127 TUPLE: record number content-length shape ;
 128
 129 : read-record ( -- record/f )
 130     4 read [ be> 4 read be> read-shape record boa ] [ f ] if* ;
 131
 132 : read-records ( -- records )
 133     [ read-record ] loop>array ;
 134
 135 : read-shp ( -- header records )
 136     read-header read-records ;
 137
 138 : file>shp ( path -- header records )
 139     binary [ read-shp ] with-file-reader ;
 140
 141 TUPLE: index offset content-length ;
 142
 143 : read-index ( -- index/f )
 144     4 read [ be> 4 read be> index boa ] [ f ] if* ;
 145
 146 : read-indices ( -- indices )
 147     [ read-index ] loop>array ;
 148
 149 : read-shx ( -- header indices )
 150     read-header read-indices ;
 151
 152 : file>shx ( path -- header indices )
 153     binary [ read-shx ] with-file-reader ;
 154
 155 : num-records ( path -- n )
 156     ".shx" append binary [
 157         read-header file-length>> 2 * 100 - 8 /
 158     ] with-file-reader ;
 159
 160 : nth-index ( n path -- index )
 161     ".shx" append binary [
 162         8 * 100 + seek-absolute seek-input read-index
 163     ] with-file-reader ;
 164
 165 : nth-record ( n path -- record )
 166     [ nth-index offset>> ] keep ".shp" append binary [
 167         2 * seek-absolute seek-input read-record
 168     ] with-file-reader ;
 169
 170 : write-int ( n -- )
 171     4 >le write ;
 172
 173 : write-double ( n -- )
 174     double>bits 8 >le write ;
 175
 176 : write-point ( point -- )
 177     [ x>> ] [ y>> ] bi [ write-double ] bi@ ;
 178
 179 :: update-box ( header shape -- header )
 180     header shape points>> :> points
 181     points [ x>> ] map minmax :> ( x-min x-max )
 182     points [ y>> ] map minmax :> ( y-min y-max )
 183     [ x-min [ or ] keep min ] change-x-min
 184     [ x-max [ or ] keep max ] change-x-max
 185     [ y-min [ or ] keep min ] change-y-min
 186     [ y-max [ or ] keep max ] change-y-max
 187     { x-min y-min x-max y-max } shape box<< ;
 188
 189 :: update-z-range ( header shape -- header )
 190     header shape z-array>> minmax :> ( z-min z-max )
 191     [ z-min [ or ] keep min ] change-z-min
 192     [ z-max [ or ] keep max ] change-z-max
 193     { z-min z-max } shape z-range<< ;
 194
 195 :: update-m-range ( header shape -- header )
 196     header shape m-array>> minmax :> ( m-min m-max )
 197     [ m-min [ or ] keep min ] change-m-min
 198     [ m-max [ or ] keep max ] change-m-max
 199     { m-min m-max } shape m-range<< ;
 200
 201 GENERIC: update-bounds ( header shape -- header )
 202
 203 M: object update-bounds drop ;
 204
 205 M: polyline update-bounds update-box ;
 206
 207 M: polygon update-bounds update-box ;
 208
 209 M: multipoint update-bounds update-box ;
 210
 211 M: polyline-z update-bounds
 212     [ call-next-method ] [ update-z-range ] [ update-m-range ] tri ;
 213
 214 M: polygon-z update-bounds
 215     [ call-next-method ] [ update-z-range ] [ update-m-range ] tri ;
 216
 217 M: multipoint-z update-bounds
 218     [ call-next-method ] [ update-z-range ] [ update-m-range ] tri ;
 219
 220 M: polyline-m update-bounds
 221     [ call-next-method ] [ update-m-range ] bi ;
 222
 223 M: polygon-m update-bounds
 224     [ call-next-method ] [ update-m-range ] bi ;
 225
 226 M: multipoint-m update-bounds
 227     [ call-next-method ] [ update-m-range ] bi ;
 228
 229 M: multipatch update-bounds
 230     [ update-box ] [ update-z-range ] [ update-m-range ] tri ;
 231
 232 GENERIC: (write-shape) ( shape -- )
 233
 234 M: null-shape (write-shape) drop ;
 235
 236 M: point (write-shape) write-point ;
 237
 238 : write-poly ( poly -- )
 239     {
 240         [ box>> [ write-double ] each ]
 241         [ parts>> length write-int ]
 242         [ points>> length write-int ]
 243         [ parts>> [ write-int ] each ]
 244         [ points>> [ write-point ] each ]
 245     } cleave ; inline
 246
 247 M: polyline (write-shape) write-poly ;
 248
 249 M: polygon (write-shape) write-poly ;
 250
 251 M: multipoint (write-shape)
 252     {
 253         [ box>> [ write-double ] each ]
 254         [ points>> length write-int ]
 255         [ points>> [ write-point ] each ]
 256     } cleave ;
 257
 258 M: point-z (write-shape)
 259     [ call-next-method ] [ z>> ] [ m>> ] tri [ write-double ] bi@ ;
 260
 261 : write-z ( shape -- )
 262     [ z-range>> ] [ z-array>> ] bi [ [ write-double ] each ] bi@ ; inline
 263
 264 : write-m ( shape -- )
 265     [ m-range>> ] [ m-array>> ] bi [ [ write-double ] each ] bi@ ; inline
 266
 267 : write-poly-z ( poly -- )
 268     [ write-poly ] [ write-z ] [ write-m ] tri ; inline
 269
 270 M: polyline-z (write-shape) write-poly-z ;
 271
 272 M: polygon-z (write-shape) write-poly-z ;
 273
 274 M: multipoint-z (write-shape)
 275     [ call-next-method ] [ write-z ] [ write-m ] tri ;
 276
 277 M: point-m (write-shape)
 278     [ call-next-method ] [ m>> write-double ] bi ;
 279
 280 : write-poly-m ( poly -- )
 281     [ write-poly ] [ write-m ] bi ; inline
 282
 283 M: polyline-m (write-shape) write-poly-m ;
 284
 285 M: polygon-m (write-shape) write-poly-m ;
 286
 287 M: multipoint-m (write-shape)
 288     [ call-next-method ] [ write-m ] bi ;
 289
 290 M: multipatch (write-shape)
 291     {
 292         [ box>> [ write-double ] each ]
 293         [ parts>> length write-int ]
 294         [ points>> length write-int ]
 295         [ parts>> [ write-int ] each ]
 296         [ part-types>> [ write-int ] each ]
 297         [ points>> [ write-point ] each ]
 298         [ write-z ]
 299         [ write-m ]
 300     } cleave ;
 301
 302 GENERIC: shape-type ( shape -- shape-type )
 303 M: null-shape shape-type drop 0 ;
 304 M: point shape-type drop 1 ;
 305 M: polyline shape-type drop 3 ;
 306 M: polygon shape-type drop 5 ;
 307 M: multipoint shape-type drop 8 ;
 308 M: point-z shape-type drop 11 ;
 309 M: polyline-z shape-type drop 13 ;
 310 M: polygon-z shape-type drop 15 ;
 311 M: multipoint-z shape-type drop 18 ;
 312 M: point-m shape-type drop 21 ;
 313 M: polyline-m shape-type drop 23 ;
 314 M: polygon-m shape-type drop 25 ;
 315 M: multipoint-m shape-type drop 28 ;
 316 M: multipatch shape-type drop 31 ;
 317
 318 : write-shape ( shape -- )
 319     [ shape-type write-int ] [ (write-shape) ] bi ;
 320
 321 : write-header ( header -- )
 322     {
 323         [ file-code>> 4 >be write ]
 324         [ drop 20 <byte-array> write ] ! unused
 325         [ file-length>> 4 >be write ]
 326         [ version>> write-int ]
 327         [ shape-type>> write-int ]
 328         [ x-min>> 0.0 or write-double ]
 329         [ y-min>> 0.0 or write-double ]
 330         [ x-max>> 0.0 or write-double ]
 331         [ y-max>> 0.0 or write-double ]
 332         [ z-min>> 0.0 or write-double ]
 333         [ z-max>> 0.0 or write-double ]
 334         [ m-min>> 0.0 or write-double ]
 335         [ m-max>> 0.0 or write-double ]
 336     } cleave ;
 337
 338 : write-record ( shape index -- )
 339     1 + 4 >be write
 340     binary [ write-shape ] with-byte-writer
 341     [ length 2/ 4 >be write ] [ write ] bi ;
 342
 343 ERROR: non-null-shapes-must-be-same-type shape-types ;
 344
 345 : non-null-shape-types ( shapes -- shape-types )
 346     [ null-shape? ] reject [ class-of ] map members ;
 347
 348 : check-shape-types ( shapes -- )
 349     non-null-shape-types dup length 1 >
 350     [ non-null-shapes-must-be-same-type ] [ drop ] if ;
 351
 352 : write-shp ( shapes -- header indices )
 353     [ header new ] dip {
 354         [ check-shape-types ]
 355         [ first shape-type >>shape-type ]
 356         [ [ update-bounds ] each ]
 357         [ ]
 358     } cleave binary [
 359         [
 360             [ tell-output 100 + 2/ ] 2dip write-record
 361             tell-output 100 + 8 - 2/ over - index boa
 362         ] map-index
 363     ] with-byte-writer swap [
 364         [ length 100 + 2/ >>file-length [ write-header ] keep ]
 365         [ write ] bi
 366     ] dip ;
 367
 368 : write-index ( index -- )
 369     [ offset>> ] [ content-length>> ] bi [ 4 >be write ] bi@ ;
 370
 371 : write-shx ( header indices -- )
 372     [ length 8 * 100 + 2/ >>file-length write-header ]
 373     [ [ write-index ] each ] bi ;
 374
 375 PRIVATE>
 376
 377 : load-shapes ( path -- shapes )
 378     ".shp" append file>shp nip [ shape>> ] map ;
 379
 380 : save-shapes ( shapes path -- )
 381     [ ".shp" append binary [ write-shp ] with-file-writer ]
 382     [ ".shx" append binary [ write-shx ] with-file-writer ] bi ;