pointcompression.go (10631B)
1 // Copyright 2017 Google Inc. All rights reserved. 2 // 3 // Licensed under the Apache License, Version 2.0 (the "License"); 4 // you may not use this file except in compliance with the License. 5 // You may obtain a copy of the License at 6 // 7 // http://www.apache.org/licenses/LICENSE-2.0 8 // 9 // Unless required by applicable law or agreed to in writing, software 10 // distributed under the License is distributed on an "AS IS" BASIS, 11 // WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 12 // See the License for the specific language governing permissions and 13 // limitations under the License. 14 15 package s2 16 17 import ( 18 "errors" 19 "fmt" 20 21 "github.com/golang/geo/r3" 22 ) 23 24 // maxEncodedVertices is the maximum number of vertices, in a row, to be encoded or decoded. 25 // On decode, this defends against malicious encodings that try and have us exceed RAM. 26 const maxEncodedVertices = 50000000 27 28 // xyzFaceSiTi represents the The XYZ and face,si,ti coordinates of a Point 29 // and, if this point is equal to the center of a Cell, the level of this cell 30 // (-1 otherwise). This is used for Loops and Polygons to store data in a more 31 // compressed format. 32 type xyzFaceSiTi struct { 33 xyz Point 34 face int 35 si, ti uint32 36 level int 37 } 38 39 const derivativeEncodingOrder = 2 40 41 func appendFace(faces []faceRun, face int) []faceRun { 42 if len(faces) == 0 || faces[len(faces)-1].face != face { 43 return append(faces, faceRun{face, 1}) 44 } 45 faces[len(faces)-1].count++ 46 return faces 47 } 48 49 // encodePointsCompressed uses an optimized compressed format to encode the given values. 50 func encodePointsCompressed(e *encoder, vertices []xyzFaceSiTi, level int) { 51 var faces []faceRun 52 for _, v := range vertices { 53 faces = appendFace(faces, v.face) 54 } 55 encodeFaces(e, faces) 56 57 type piQi struct { 58 pi, qi uint32 59 } 60 verticesPiQi := make([]piQi, len(vertices)) 61 for i, v := range vertices { 62 verticesPiQi[i] = piQi{siTitoPiQi(v.si, level), siTitoPiQi(v.ti, level)} 63 } 64 piCoder, qiCoder := newNthDerivativeCoder(derivativeEncodingOrder), newNthDerivativeCoder(derivativeEncodingOrder) 65 for i, v := range verticesPiQi { 66 f := encodePointCompressed 67 if i == 0 { 68 // The first point will be just the (pi, qi) coordinates 69 // of the Point. NthDerivativeCoder will not save anything 70 // in that case, so we encode in fixed format rather than varint 71 // to avoid the varint overhead. 72 f = encodeFirstPointFixedLength 73 } 74 f(e, v.pi, v.qi, level, piCoder, qiCoder) 75 } 76 77 var offCenter []int 78 for i, v := range vertices { 79 if v.level != level { 80 offCenter = append(offCenter, i) 81 } 82 } 83 e.writeUvarint(uint64(len(offCenter))) 84 for _, idx := range offCenter { 85 e.writeUvarint(uint64(idx)) 86 e.writeFloat64(vertices[idx].xyz.X) 87 e.writeFloat64(vertices[idx].xyz.Y) 88 e.writeFloat64(vertices[idx].xyz.Z) 89 } 90 } 91 92 func encodeFirstPointFixedLength(e *encoder, pi, qi uint32, level int, piCoder, qiCoder *nthDerivativeCoder) { 93 // Do not ZigZagEncode the first point, since it cannot be negative. 94 codedPi, codedQi := piCoder.encode(int32(pi)), qiCoder.encode(int32(qi)) 95 // Interleave to reduce overhead from two partial bytes to one. 96 interleaved := interleaveUint32(uint32(codedPi), uint32(codedQi)) 97 98 // Write as little endian. 99 bytesRequired := (level + 7) / 8 * 2 100 for i := 0; i < bytesRequired; i++ { 101 e.writeUint8(uint8(interleaved)) 102 interleaved >>= 8 103 } 104 } 105 106 // encodePointCompressed encodes points into e. 107 // Given a sequence of Points assumed to be the center of level-k cells, 108 // compresses it into a stream using the following method: 109 // - decompose the points into (face, si, ti) tuples. 110 // - run-length encode the faces, combining face number and count into a 111 // varint32. See the faceRun struct. 112 // - right shift the (si, ti) to remove the part that's constant for all cells 113 // of level-k. The result is called the (pi, qi) space. 114 // - 2nd derivative encode the pi and qi sequences (linear prediction) 115 // - zig-zag encode all derivative values but the first, which cannot be 116 // negative 117 // - interleave the zig-zag encoded values 118 // - encode the first interleaved value in a fixed length encoding 119 // (varint would make this value larger) 120 // - encode the remaining interleaved values as varint64s, as the 121 // derivative encoding should make the values small. 122 // In addition, provides a lossless method to compress a sequence of points even 123 // if some points are not the center of level-k cells. These points are stored 124 // exactly, using 3 double precision values, after the above encoded string, 125 // together with their index in the sequence (this leads to some redundancy - it 126 // is expected that only a small fraction of the points are not cell centers). 127 // 128 // To encode leaf cells, this requires 8 bytes for the first vertex plus 129 // an average of 3.8 bytes for each additional vertex, when computed on 130 // Google's geographic repository. 131 func encodePointCompressed(e *encoder, pi, qi uint32, level int, piCoder, qiCoder *nthDerivativeCoder) { 132 // ZigZagEncode, as varint requires the maximum number of bytes for 133 // negative numbers. 134 zzPi := zigzagEncode(piCoder.encode(int32(pi))) 135 zzQi := zigzagEncode(qiCoder.encode(int32(qi))) 136 // Interleave to reduce overhead from two partial bytes to one. 137 interleaved := interleaveUint32(zzPi, zzQi) 138 e.writeUvarint(interleaved) 139 } 140 141 type faceRun struct { 142 face, count int 143 } 144 145 func decodeFaceRun(d *decoder) faceRun { 146 faceAndCount := d.readUvarint() 147 ret := faceRun{ 148 face: int(faceAndCount % numFaces), 149 count: int(faceAndCount / numFaces), 150 } 151 if ret.count <= 0 && d.err == nil { 152 d.err = errors.New("non-positive count for face run") 153 } 154 return ret 155 } 156 157 func decodeFaces(numVertices int, d *decoder) []faceRun { 158 var frs []faceRun 159 for nparsed := 0; nparsed < numVertices; { 160 fr := decodeFaceRun(d) 161 if d.err != nil { 162 return nil 163 } 164 frs = append(frs, fr) 165 nparsed += fr.count 166 } 167 return frs 168 } 169 170 // encodeFaceRun encodes each faceRun as a varint64 with value numFaces * count + face. 171 func encodeFaceRun(e *encoder, fr faceRun) { 172 // It isn't necessary to encode the number of faces left for the last run, 173 // but since this would only help if there were more than 21 faces, it will 174 // be a small overall savings, much smaller than the bound encoding. 175 coded := numFaces*uint64(fr.count) + uint64(fr.face) 176 e.writeUvarint(coded) 177 } 178 179 func encodeFaces(e *encoder, frs []faceRun) { 180 for _, fr := range frs { 181 encodeFaceRun(e, fr) 182 } 183 } 184 185 type facesIterator struct { 186 faces []faceRun 187 // How often have we yet shown the current face? 188 numCurrentFaceShown int 189 curFace int 190 } 191 192 func (fi *facesIterator) next() (ok bool) { 193 if len(fi.faces) == 0 { 194 return false 195 } 196 fi.curFace = fi.faces[0].face 197 fi.numCurrentFaceShown++ 198 199 // Advance fs if needed. 200 if fi.faces[0].count <= fi.numCurrentFaceShown { 201 fi.faces = fi.faces[1:] 202 fi.numCurrentFaceShown = 0 203 } 204 205 return true 206 } 207 208 func decodePointsCompressed(d *decoder, level int, target []Point) { 209 faces := decodeFaces(len(target), d) 210 211 piCoder := newNthDerivativeCoder(derivativeEncodingOrder) 212 qiCoder := newNthDerivativeCoder(derivativeEncodingOrder) 213 214 iter := facesIterator{faces: faces} 215 for i := range target { 216 decodeFn := decodePointCompressed 217 if i == 0 { 218 decodeFn = decodeFirstPointFixedLength 219 } 220 pi, qi := decodeFn(d, level, piCoder, qiCoder) 221 if ok := iter.next(); !ok && d.err == nil { 222 d.err = fmt.Errorf("ran out of faces at target %d", i) 223 return 224 } 225 target[i] = Point{facePiQitoXYZ(iter.curFace, pi, qi, level)} 226 } 227 228 numOffCenter := int(d.readUvarint()) 229 if d.err != nil { 230 return 231 } 232 if numOffCenter > len(target) { 233 d.err = fmt.Errorf("numOffCenter = %d, should be at most len(target) = %d", numOffCenter, len(target)) 234 return 235 } 236 for i := 0; i < numOffCenter; i++ { 237 idx := int(d.readUvarint()) 238 if d.err != nil { 239 return 240 } 241 if idx >= len(target) { 242 d.err = fmt.Errorf("off center index = %d, should be < len(target) = %d", idx, len(target)) 243 return 244 } 245 target[idx].X = d.readFloat64() 246 target[idx].Y = d.readFloat64() 247 target[idx].Z = d.readFloat64() 248 } 249 } 250 251 func decodeFirstPointFixedLength(d *decoder, level int, piCoder, qiCoder *nthDerivativeCoder) (pi, qi uint32) { 252 bytesToRead := (level + 7) / 8 * 2 253 var interleaved uint64 254 for i := 0; i < bytesToRead; i++ { 255 rr := d.readUint8() 256 interleaved |= (uint64(rr) << uint(i*8)) 257 } 258 259 piCoded, qiCoded := deinterleaveUint32(interleaved) 260 261 return uint32(piCoder.decode(int32(piCoded))), uint32(qiCoder.decode(int32(qiCoded))) 262 } 263 264 func zigzagEncode(x int32) uint32 { 265 return (uint32(x) << 1) ^ uint32(x>>31) 266 } 267 268 func zigzagDecode(x uint32) int32 { 269 return int32((x >> 1) ^ uint32((int32(x&1)<<31)>>31)) 270 } 271 272 func decodePointCompressed(d *decoder, level int, piCoder, qiCoder *nthDerivativeCoder) (pi, qi uint32) { 273 interleavedZigZagEncodedDerivPiQi := d.readUvarint() 274 piZigzag, qiZigzag := deinterleaveUint32(interleavedZigZagEncodedDerivPiQi) 275 return uint32(piCoder.decode(zigzagDecode(piZigzag))), uint32(qiCoder.decode(zigzagDecode(qiZigzag))) 276 } 277 278 // We introduce a new coordinate system (pi, qi), which is (si, ti) 279 // with the bits that are constant for cells of that level shifted 280 // off to the right. 281 // si = round(s * 2^31) 282 // pi = si >> (31 - level) 283 // = floor(s * 2^level) 284 // If the point has been snapped to the level, the bits that are 285 // shifted off will be a 1 in the msb, then 0s after that, so the 286 // fractional part discarded by the cast is (close to) 0.5. 287 288 // stToPiQi returns the value transformed to the PiQi coordinate space. 289 func stToPiQi(s float64, level uint) uint32 { 290 return uint32(s * float64(int(1)<<level)) 291 } 292 293 // siTiToPiQi returns the value transformed into the PiQi coordinate spade. 294 // encodeFirstPointFixedLength encodes the return value using level bits, 295 // so we clamp si to the range [0, 2**level - 1] before trying to encode 296 // it. This is okay because if si == maxSiTi, then it is not a cell center 297 // anyway and will be encoded separately as an off-center point. 298 func siTitoPiQi(siTi uint32, level int) uint32 { 299 s := uint(siTi) 300 const max = maxSiTi - 1 301 if s > max { 302 s = max 303 } 304 305 return uint32(s >> (maxLevel + 1 - uint(level))) 306 } 307 308 // piQiToST returns the value transformed to ST space. 309 func piQiToST(pi uint32, level int) float64 { 310 // We want to recover the position at the center of the cell. If the point 311 // was snapped to the center of the cell, then math.Modf(s * 2^level) == 0.5. 312 // Inverting STtoPiQi gives: 313 // s = (pi + 0.5) / 2^level. 314 return (float64(pi) + 0.5) / float64(int(1)<<uint(level)) 315 } 316 317 func facePiQitoXYZ(face int, pi, qi uint32, level int) r3.Vector { 318 return faceUVToXYZ(face, stToUV(piQiToST(pi, level)), stToUV(piQiToST(qi, level))).Normalize() 319 }