/alaspatial/src/main/java/org/ala/spatial/util/ComplexRegion.java
Java | 622 lines | 413 code | 66 blank | 143 comment | 201 complexity | 630d210101a2e30699ccfd159334021a MD5 | raw file
1package org.ala.spatial.util; 2 3import java.util.ArrayList; 4 5/** 6 * ComplexRegion is a collection of SimpleRegion, expect POLYGONs for now. 7 * 8 * treat as a shape file, overlapping regions cancel out presence. 9 * 10 * TODO: clockwise/anticlockwise identification 11 * 12 * @author Adam Collins 13 */ 14public class ComplexRegion extends SimpleRegion { 15 16 static public SimpleRegion parseComplexRegion(String[] polygons) { 17 ComplexRegion cr = new ComplexRegion(); 18 19 for (String s : polygons) { 20 cr.addPolygon(SimpleRegion.parseSimpleRegion(s)); 21 } 22 23 cr.useMask(-1, -1, -1); 24 25 return cr; 26 } 27 /** 28 * list of SimpleRegion members 29 */ 30 ArrayList<SimpleRegion> simpleregions; 31 /** 32 * bounding box for all, see SimpleRegion boundingbox. 33 */ 34 double[][] boundingbox_all; 35 /** 36 * value assigned 37 */ 38 int value; 39 /** 40 * array for speeding up isWithin 41 */ 42 byte[][] mask; 43 Object[][] maskDepth; 44 /** 45 * mask height 46 */ 47 int mask_height; 48 /** 49 * mask width 50 */ 51 int mask_width; 52 /** 53 * mask multiplier for longitude inputs 54 */ 55 double mask_long_multiplier; 56 /** 57 * mask mulitplier for latitude inputs 58 */ 59 double mask_lat_multiplier; 60 /** 61 * maintain mapping for simpleregions belonging to the same polygon 62 */ 63 ArrayList<Integer> polygons; 64 65 /** 66 * Constructor for empty ComplexRegion 67 */ 68 public ComplexRegion() { 69 super(); 70 71 simpleregions = new ArrayList(); 72 boundingbox_all = new double[2][2]; 73 value = -1; 74 mask = null; 75 polygons = new ArrayList(); 76 } 77 78 /** 79 * sets integer value stored 80 * @param value_ as int 81 */ 82 public void setValue(int value_) { 83 value = value_; 84 } 85 86 /** 87 * gets integer value stored 88 * @return int 89 */ 90 public int getValue() { 91 return value; 92 } 93 94 /** 95 * gets the bounding box for shapes in this ComplexRegion 96 * 97 * @return bounding box for ComplexRegion as double [][] 98 */ 99 @Override 100 public double[][] getBoundingBox() { 101 return boundingbox_all; 102 } 103 104 /** 105 * adds a new polygon 106 * 107 * note: if a mask is in use must call <code>useMask</code> again 108 * @param points_ points = double[n][2] 109 * where 110 * n is number of points 111 * [][0] is longitude 112 * [][1] is latitude 113 */ 114 public void addPolygon(SimpleRegion sr) { 115 simpleregions.add(sr); 116 117 /* update boundingbox_all */ 118 double[][] bb = sr.getBoundingBox(); 119 if (simpleregions.size() == 1 || boundingbox_all[0][0] > bb[0][0]) { 120 boundingbox_all[0][0] = bb[0][0]; 121 } 122 if (simpleregions.size() == 1 || boundingbox_all[1][0] < bb[1][0]) { 123 boundingbox_all[1][0] = bb[1][0]; 124 } 125 if (simpleregions.size() == 1 || boundingbox_all[0][1] > bb[0][1]) { 126 boundingbox_all[0][1] = bb[0][1]; 127 } 128 if (simpleregions.size() == 1 || boundingbox_all[1][1] < bb[1][1]) { 129 boundingbox_all[1][1] = bb[1][1]; 130 } 131 132 bounding_box = boundingbox_all; 133 } 134 135 /** 136 * returns true when the point provided is within the ComplexRegion 137 * 138 * uses <code>mask</code> when available 139 * 140 * note: type UNDEFINED implies no boundary, always returns true. 141 * 142 * @param longitude 143 * @param latitude 144 * @return true iff point is within or on the edge of this ComplexRegion 145 */ 146 @Override 147 public boolean isWithin(double longitude, double latitude) { 148 if (simpleregions.size() == 1) { 149 return simpleregions.get(0).isWithin(longitude, latitude); 150 } 151 if (boundingbox_all[0][0] > longitude || boundingbox_all[1][0] < longitude 152 || boundingbox_all[0][1] > latitude || boundingbox_all[1][1] < latitude) { 153 return false; 154 } 155 156 short[] countsIn = new short[polygons.get(polygons.size() - 1) + 1]; 157 //int count_in = 0; //count of regions overlapping the point 158 if (mask != null) { 159 /* use mask if exists */ 160 int long1 = (int) Math.floor((longitude - boundingbox_all[0][0]) * mask_long_multiplier); 161 int lat1 = (int) Math.floor((latitude - boundingbox_all[0][1]) * mask_lat_multiplier); 162 163 if (long1 == mask[0].length) { 164 long1--; 165 } 166 if (lat1 == mask.length) { 167 lat1--; 168 } 169 170 if (mask[lat1][long1] == SimpleRegion.GI_FULLY_PRESENT) { 171 return true; 172 } else if (mask[lat1][long1] == SimpleRegion.GI_UNDEFINED 173 || mask[lat1][long1] == SimpleRegion.GI_ABSENCE) { 174 return false; 175 } 176 //partial, try maskDepth and sum overlaps 177 if (maskDepth != null && maskDepth[lat1][long1] != null) { 178 int[] d = (int[]) maskDepth[lat1][long1]; 179 for (int i = 0; i < d.length; i++) { 180 if (simpleregions.get(d[i]).isWithin(longitude, latitude)) { 181 countsIn[polygons.get(d[i])]++; 182 } 183 } 184 /* true iif within an odd number of regions for any polygon*/ 185 for (int i = 0; i < countsIn.length; i++) { 186 if (countsIn[i] % 2 == 1) { 187 return true; 188 } 189 } 190 } 191 192 } 193 194 /* check for all SimpleRegions */ 195 for (int i = 0; i < simpleregions.size(); i++) { 196 if (simpleregions.get(i).isWithin(longitude, latitude)) { 197 countsIn[polygons.get(i)]++; 198 } 199 } 200 201 /* true iif within an odd number of regions for any polygon*/ 202 for (int i = 0; i < countsIn.length; i++) { 203 if (countsIn[i] % 2 == 1) { 204 return true; 205 } 206 } 207 return false; 208 } 209 210 /** 211 * builds a grid (mask) to speed up isWithin. 212 * 213 * TODO: split out shapes with large numbers of points in GI_PARTIALLY_PRESENT grid cells. 214 * 215 * TODO: automatic(best) height/width specification 216 * 217 * @param width 218 * @param height 219 */ 220 public void useMask(int width, int height, int depthThreashold) { 221 //calculate defaults for -1 inputs 222 double[][] bb = getBoundingBox(); 223 int length = 0; 224 for (SimpleRegion sr : simpleregions) { 225 length += sr.getNumberOfPoints(); 226 } 227 int w = (int) ((bb[1][0] - bb[0][0]) * 3); 228 int h = (int) ((bb[1][1] - bb[0][1]) * 3); 229 if (length > 5000) { 230 w = 200; 231 h = 200; 232 } 233 if (w > 200) { 234 w = 200; 235 } 236 if (h > 200) { 237 h = 200; 238 } 239 if (width == -1) { 240 width = w; 241 } 242 if (height == -1) { 243 height = h; 244 } 245 if (depthThreashold == -1) { 246 depthThreashold = 100; 247 } 248 if (width < 3 || height < 3) { 249 return; 250 } 251 252 int i, j; 253 254 /* class variables assignment */ 255 mask_width = width; 256 mask_height = height; 257 mask_long_multiplier = 258 mask_width / (double) (boundingbox_all[1][0] - boundingbox_all[0][0]); 259 mask_lat_multiplier = 260 mask_height / (double) (boundingbox_all[1][1] - boundingbox_all[0][1]); 261 262 /* end result mask */ 263 mask = new byte[height][width]; 264 ArrayList<Integer>[][] md = null; 265 if (simpleregions.size() > depthThreashold) { 266 //use mask depth as well 267 md = new ArrayList[height][width]; 268 } 269 270 /* temp mask for current Polygon */ 271 byte[][] shapemask = new byte[height][width]; 272 273 /* temp mask for current SimpleRegion */ 274 byte[][] shapemaskregion = new byte[height][width]; 275 276 int k = 0; 277 while (k < simpleregions.size()) { 278 int p = k; 279 for (; k < simpleregions.size() 280 && (p == k || polygons.get(k - 1) == polygons.get(k)); k++) { 281 282 SimpleRegion sr = simpleregions.get(k); 283 sr.getOverlapGridCells(boundingbox_all[0][0], boundingbox_all[0][1], boundingbox_all[1][0], boundingbox_all[1][1], width, height, shapemaskregion); 284 285 //shapemaskregion into shapemask 286 for (i = 0; i < height; i++) { 287 for (j = 0; j < width; j++) { 288 if (shapemaskregion[i][j] == SimpleRegion.GI_PARTIALLY_PRESENT 289 || shapemask[i][j] == SimpleRegion.GI_PARTIALLY_PRESENT) { 290 shapemask[i][j] = SimpleRegion.GI_PARTIALLY_PRESENT; //partially inside 291 if (md != null) { 292 if (md[i][j] == null) { 293 md[i][j] = new ArrayList<Integer>(); 294 } 295 md[i][j].add(k); 296 } 297 } else if (shapemaskregion[i][j] == SimpleRegion.GI_FULLY_PRESENT) { 298 if (shapemask[i][j] == SimpleRegion.GI_FULLY_PRESENT) { 299 shapemask[i][j] = SimpleRegion.GI_ABSENCE; //completely inside 300 } else { 301 shapemask[i][j] = SimpleRegion.GI_FULLY_PRESENT;//completely outside (inside of a cutout region) 302 } 303 if (md != null) { 304 if (md[i][j] == null) { 305 md[i][j] = new ArrayList<Integer>(); 306 } 307 md[i][j].add(k); 308 } 309 } 310 311 /* reset shapemaskregion for next part */ 312 shapemaskregion[i][j] = 0; 313 } 314 } 315 } 316 317 //shapemask into mask 318 for (i = 0; i < height; i++) { 319 for (j = 0; j < width; j++) { 320 if (shapemask[i][j] == SimpleRegion.GI_FULLY_PRESENT 321 || mask[i][j] == SimpleRegion.GI_FULLY_PRESENT) { 322 mask[i][j] = SimpleRegion.GI_FULLY_PRESENT; //partially inside 323 } else if (shapemask[i][j] == SimpleRegion.GI_PARTIALLY_PRESENT) { 324 mask[i][j] = SimpleRegion.GI_PARTIALLY_PRESENT; 325 } 326 327 /* reset shapemask for next part */ 328 shapemask[i][j] = 0; 329 } 330 } 331 } 332 333 //maskDepth to int[] 334 if (md != null) { 335 maskDepth = new Object[md.length][md[0].length]; 336 for (i = 0; i < height; i++) { 337 for (j = 0; j < width; j++) { 338 if (md[i][j] != null && mask[i][j] == SimpleRegion.GI_PARTIALLY_PRESENT) { 339 int[] d = new int[md[i][j].size()]; 340 for (k = 0; k < d.length; k++) { 341 d[k] = md[i][j].get(k); 342 } 343 maskDepth[i][j] = d; 344 } 345 } 346 } 347 } 348 } 349 350 /** 351 * determines overlap with a grid 352 * 353 * for type POLYGON 354 * when <code>three_state_map</code> is not null populate it with one of: 355 * GI_UNDEFINED 356 * GI_PARTIALLY_PRESENT 357 * GI_FULLY_PRESENT 358 * GI_ABSENCE 359 * 360 * @param longitude1 361 * @param latitude1 362 * @param longitude2 363 * @param latitude2 364 * @param xres number of longitude segements as int 365 * @param yres number of latitude segments as int 366 * @return (x,y) as double [][2] for each grid cell at least partially falling 367 * within the specified region of the specified resolution beginning at 0,0 368 * for minimum longitude and latitude through to xres,yres for maximums * 369 */ 370 @Override 371 public int[][] getOverlapGridCells(double longitude1, double latitude1, double longitude2, double latitude2, int width, int height, byte[][] three_state_map, boolean noCellsReturned) { 372 int i, j; 373 // if no threestate map exists, create one 374 if (mask == null) { 375 mask = new byte[height][width]; 376 three_state_map = mask; 377 } 378 379 byte[][] tmpMask = new byte[height][width]; 380 381 int k = 0; 382 while (k < simpleregions.size()) { 383 int p = k; 384 for (; k < simpleregions.size() 385 && (p == k || polygons.get(k - 1).equals(polygons.get(k))); k++) { 386 387 SimpleRegion sr = simpleregions.get(k); 388 sr.getOverlapGridCells_Acc(longitude1, latitude1, longitude2, latitude2, width, height, tmpMask); 389 } 390 391 fillAccMask(p, k - 1, longitude1, latitude1, longitude2, latitude2, width, height, tmpMask, true); 392 393 //tmpMask into mask 394 for (i = 0; i < height; i++) { 395 for (j = 0; j < width; j++) { 396 if (tmpMask[i][j] == 2 || mask[i][j] == 2) { 397 mask[i][j] = 2; 398 } else if (tmpMask[i][j] == 1) { 399 mask[i][j] = 1; 400 } 401 402 /* reset shapemask for next part */ 403 tmpMask[i][j] = 0; 404 } 405 } 406 } 407 408 boolean cellsReturned = !noCellsReturned; 409 if (cellsReturned) { 410 int[][] data = new int[width * height][2]; 411 int p = 0; 412 for (i = 0; i < height; i++) { 413 for (j = 0; j < width; j++) { 414 if (mask[i][j] != GI_UNDEFINED) { //undefined == absence 415 data[p][0] = j; 416 data[p][1] = i; 417 p++; 418 } 419 } 420 } 421 data = java.util.Arrays.copyOf(data, p); 422 } 423 424 425 return null; 426 } 427 428 @Override 429 public boolean isWithin_EPSG900913(double longitude, double latitude) { 430 short[] countsIn = new short[polygons.get(polygons.size() - 1) + 1]; 431 /* check for all SimpleRegions */ 432 for (int i = 0; i < simpleregions.size(); i++) { 433 if (simpleregions.get(i).isWithin_EPSG900913(longitude, latitude)) { 434 countsIn[polygons.get(i)]++; 435 } 436 } 437 438 /* true iif within an odd number of regions for any polygon*/ 439 for (int i = 0; i < countsIn.length; i++) { 440 if (countsIn[i] % 2 == 1) { 441 return true; 442 } 443 } 444 445 return false; 446 } 447 448 @Override 449 public int[][] getOverlapGridCells_EPSG900913(double longitude1, double latitude1, double longitude2, double latitude2, int width, int height, byte[][] three_state_map, boolean noCellsReturned) { 450 int i, j; 451 452 int[][] output = null; 453 byte[][] mask = three_state_map; 454 455 // if no threestate map exists, create one 456 if (mask == null) { 457 mask = new byte[height][width]; 458 three_state_map = mask; 459 } 460 461 byte[][] tmpMask = new byte[height][width]; 462 463 int k = 0; 464 while (k < simpleregions.size()) { 465 int p = k; 466 for (; k < simpleregions.size() 467 && (p == k || polygons.get(k - 1).equals(polygons.get(k))); k++) { 468 469 SimpleRegion sr = simpleregions.get(k); 470 sr.getOverlapGridCells_Acc_EPSG900913(longitude1, latitude1, longitude2, latitude2, width, height, tmpMask); 471 } 472 473 fillAccMask_EPSG900913(k, p - 1, longitude1, latitude1, longitude2, latitude2, width, height, tmpMask, true); 474 475 //tmpMask into mask 476 for (i = 0; i < height; i++) { 477 for (j = 0; j < width; j++) { 478 if (tmpMask[i][j] == 2 || mask[i][j] == 2) { 479 mask[i][j] = 2; 480 } else if (tmpMask[i][j] == 1) { 481 mask[i][j] = 1; 482 } 483 484 /* reset shapemask for next part */ 485 tmpMask[i][j] = 0; 486 } 487 } 488 } 489 490 boolean cellsReturned = !noCellsReturned; 491 if (cellsReturned) { 492 int[][] data = new int[width * height][2]; 493 int p = 0; 494 for (i = 0; i < height; i++) { 495 for (j = 0; j < width; j++) { 496 if (mask[i][j] != GI_UNDEFINED) { //undefined == absence 497 data[p][0] = j; 498 data[p][1] = i; 499 p++; 500 } 501 } 502 } 503 data = java.util.Arrays.copyOf(data, p); 504 return data; 505 } 506 507 return null; 508 } 509 510 void addSet(ArrayList<SimpleRegion> simpleRegions) { 511 int nextSetNumber = (polygons.size() > 0) ? polygons.get(polygons.size() - 1) + 1 : 0; 512 513 for (int i = 0; i < simpleRegions.size(); i++) { 514 addPolygon(simpleRegions.get(i)); 515 polygons.add(nextSetNumber); 516 } 517 } 518 519 public int[][] fillAccMask(int startPolygon, int endPolygon, double longitude1, double latitude1, double longitude2, double latitude2, int width, int height, byte[][] three_state_map, boolean noCellsReturned) { 520 double divx = (longitude2 - longitude1) / width; 521 double divy = (latitude2 - latitude1) / height; 522 523 int i, j; 524 //do raster check 525 int[][] data = null; 526 boolean cellsReturned = !noCellsReturned; 527 if (cellsReturned) { 528 data = new int[width * height][2]; 529 } 530 int p = 0; 531 for (j = 0; j < three_state_map[0].length; j++) { 532 for (i = 0; i < three_state_map.length; i++) { 533 if (three_state_map[i][j] == GI_PARTIALLY_PRESENT) { 534 //if it is partially present, do nothing 535 } else if ((j == 0 || three_state_map[i][j - 1] == GI_PARTIALLY_PRESENT)) { 536 if (i > 0 537 && (three_state_map[i - 1][j] == GI_FULLY_PRESENT 538 || three_state_map[i - 1][j] == GI_ABSENCE)) { 539 //use same as LHS 540 three_state_map[i][j] = three_state_map[i - 1][j]; 541 } else { 542 int count = 0; 543 for (int k = startPolygon; k <= endPolygon; k++) { 544 count += isWithin(j * divx + divx / 2 + longitude1, i * divy + divy / 2 + latitude1) ? 1 : 0; 545 } 546 if (count % 2 == 1) { 547 //if the previous was partially present, test 548 three_state_map[i][j] = GI_FULLY_PRESENT; 549 } //else absent 550 } //else absent 551 } else { 552 //if the previous was fully present, repeat 553 //if the previous was absent, repeat 554 three_state_map[i][j] = three_state_map[i][j - 1]; 555 } 556 557 //apply to cells; 558 if (cellsReturned && three_state_map[i][j] != GI_UNDEFINED) { //undefined == absence 559 data[p][0] = j; 560 data[p][1] = i; 561 p++; 562 } 563 } 564 } 565 if (data != null) { 566 data = java.util.Arrays.copyOf(data, p); 567 } 568 return data; 569 } 570 571 public int[][] fillAccMask_EPSG900913(int startPolygon, int endPolygon, double longitude1, double latitude1, double longitude2, double latitude2, int width, int height, byte[][] three_state_map, boolean noCellsReturned) { 572 double divx = (longitude2 - longitude1) / width; 573 double divy = (latitude2 - latitude1) / height; 574 575 int i, j; 576 //do raster check 577 int[][] data = null; 578 boolean cellsReturned = !noCellsReturned; 579 if (cellsReturned) { 580 data = new int[width * height][2]; 581 } 582 int p = 0; 583 for (j = 0; j < three_state_map[0].length; j++) { 584 for (i = 0; i < three_state_map.length; i++) { 585 if (three_state_map[i][j] == GI_PARTIALLY_PRESENT) { 586 //if it is partially present, do nothing 587 } else if ((j == 0 || three_state_map[i][j - 1] == GI_PARTIALLY_PRESENT)) { 588 if (i > 0 589 && (three_state_map[i - 1][j] == GI_FULLY_PRESENT 590 || three_state_map[i - 1][j] == GI_ABSENCE)) { 591 //use same as LHS 592 three_state_map[i][j] = three_state_map[i - 1][j]; 593 } else { 594 int count = 0; 595 for (int k = startPolygon; k <= endPolygon; k++) { 596 count += isWithin_EPSG900913(j * divx + divx / 2 + longitude1, i * divy + divy / 2 + latitude1) ? 1 : 0; 597 } 598 if (count % 2 == 1) { 599 //if the previous was partially present, test 600 three_state_map[i][j] = GI_FULLY_PRESENT; 601 } //else absent 602 } //else absent 603 } else { 604 //if the previous was fully present, repeat 605 //if the previous was absent, repeat 606 three_state_map[i][j] = three_state_map[i][j - 1]; 607 } 608 609 //apply to cells; 610 if (cellsReturned && three_state_map[i][j] != GI_UNDEFINED) { //undefined == absence 611 data[p][0] = j; 612 data[p][1] = i; 613 p++; 614 } 615 } 616 } 617 if (data != null) { 618 data = java.util.Arrays.copyOf(data, p); 619 } 620 return data; 621 } 622}