Python geos.WKBWriter类代码示例

def convert_to_2d(geom):
    """Convert a geometry from 3D to 2D"""
    from django.contrib.gis.geos import WKBWriter, WKBReader
    wkb_r = WKBReader()
    wkb_w = WKBWriter()
    wkb_w.outdim = 2
    return wkb_r.read(wkb_w.write(geom))

 def _handle_geom(self, geometry):
     """ Geometry processing (in place), depending on options """
     # Optional force 2D
     if self.options.get('force2d'):
         wkb_w = WKBWriter()
         wkb_w.outdim = 2
         geometry = GEOSGeometry(wkb_w.write(geometry), srid=geometry.srid)
     # Optional geometry simplification
     simplify = self.options.get('simplify')
     if simplify is not None:
         geometry = geometry.simplify(tolerance=simplify, preserve_topology=True)
     # Optional geometry reprojection
     if self.srid != geometry.srid:
         geometry.transform(self.srid)
     return geometry

 def _prepare_geometry(self, field, value):
     """
     Reduce geometry to two dimensions if models.GeometryField
     dim parameter does not request otherwise.
     """
     # TODO: DO we need to add 2D to 3D conversation filling z
     # with 0s or is this taken care of implicitely?
     # (cannot be easily tested since 3d is not supported by
     # spatialite backend in Django 1.6.)
     from django.contrib.gis.geos import WKBWriter, GEOSGeometry
     if isinstance(value, (str, unicode)):
         value = GEOSGeometry(value)
     wkb_writer = WKBWriter()
     if isinstance(value, GEOSGeometry):
         if value.hasz and field.dim == 2:
             value = GEOSGeometry(wkb_writer.write(value))
     return value

    def _handle_geom(self, value):
        """ Geometry processing (in place), depending on options """
        if value is None:
            geometry = None
        elif isinstance(value, dict) and 'type' in value:
            geometry = value
        else:
            if isinstance(value, GEOSGeometry):
                geometry = value
            else:
                try:
                    # this will handle string representations (e.g. ewkt, bwkt)
                    geometry = GEOSGeometry(value)
                except ValueError:
                    # if the geometry couldn't be parsed.
                    # we can't generate valid geojson
                    error_msg = 'The field ["%s", "%s"] could not be parsed as a valid geometry' % (
                        self.geometry_field, value
                    )
                    raise SerializationError(error_msg)

            # Optional force 2D
            if self.options.get('force2d'):
                wkb_w = WKBWriter()
                wkb_w.outdim = 2
                geometry = GEOSGeometry(wkb_w.write(geometry), srid=geometry.srid)
            # Optional geometry simplification
            simplify = self.options.get('simplify')
            if simplify is not None:
                geometry = geometry.simplify(tolerance=simplify, preserve_topology=True)
            # Optional geometry reprojection
            if geometry.srid and geometry.srid != self.srid:
                geometry.transform(self.srid)
            # Optional bbox
            if self.options.get('bbox_auto'):
                self._current['bbox'] = geometry.extent

        self._current['geometry'] = geometry

 def test_empty_polygon_wkb(self):
     p = Polygon(srid=4326)
     p_no_srid = Polygon()
     wkb_w = WKBWriter()
     wkb_w.srid = True
     for byteorder, hexes in enumerate([
         (b'000000000300000000', b'0020000003000010E600000000'),
         (b'010300000000000000', b'0103000020E610000000000000'),
     ]):
         wkb_w.byteorder = byteorder
         for srid, hex in enumerate(hexes):
             wkb_w.srid = srid
             self.assertEqual(wkb_w.write_hex(p), hex)
             self.assertEqual(GEOSGeometry(wkb_w.write_hex(p)), p if srid else p_no_srid)
             self.assertEqual(wkb_w.write(p), memoryview(binascii.a2b_hex(hex)))
             self.assertEqual(GEOSGeometry(wkb_w.write(p)), p if srid else p_no_srid)

def run():
    # Get data directory from environment
    datadir = os.environ.get('CORINE_DATA_DIRECTORY', '')
    if not datadir:
        print('Datadir not found, please specify CORINE_DATA_DIRECTORY env var.')
        return

    wkb_w = WKBWriter()
    wkb_w.outdim = 2

    sources = sorted(glob.glob(os.path.join(datadir, '*.sqlite')))

    print('Processing files', sources)

    for source in sources:
        # Detect file content either landcover or landcover change
        change = re.findall(r'^cha([^\_*\.sqlite]+)', os.path.basename(source))
        normal = re.findall(r'^clc([^\_*\.sqlite]+)', os.path.basename(source))

        if len(normal):
            # Select field mapping for landcover files
            mapping = const.FIELD_MAPPING

            # Get current year from regex match
            year = normal[0]

            # Set change flag
            change = False

        elif len(change):
            # Select field mapping for change files
            mapping = const.CHANGE_FIELD_MAPPING

            # Get current year from regex match
            year = change[0]

            # Get previous year based on base year
            previous = const.PREVIOUS_LOOKUP[year]
            code_previous_mapping = 'code_' + previous

            # Set change flag
            change = True

        else:
            raise ValueError('Could not interpret source.')

        # Mapping for the landcover code field source
        code_mapping = 'code_' + year

        # Convert regex match year to full year
        year = const.YEAR_MAPPING[year]

        Patch.objects.filter(year=year, change=change).delete()

        print('Processing {}data for year {}.'.format('change ' if change else '', year))

        # Get full nomenclature from nomenclature app. Convert to dict for speed.
        nomenclature = {int(x.code): x.id for x in Nomenclature.objects.all()}

        # Open datasource
        ds = DataSource(source)
        # Get layer from datasource
        lyr = ds[0]

        # Initiate counter and batch array
        counter = 0
        batch = []

        # Process features in layer
        for feat in lyr:
            counter += 1

            # Create patch instance without commiting
            patch = Patch(
                year=year,
                change=change,
                nomenclature_id=nomenclature[int(feat.get(code_mapping))],
            )

            try:
                # Make sure geom is a multi polygon
                multi = feat.geom.geos
                if multi.geom_type != 'MultiPolygon':
                    multi = MultiPolygon(multi)

                # If necessary, roundtrip through hex writer to drop z dim
                if multi.hasz:
                    multi = GEOSGeometry(wkb_w.write_hex(multi))

                patch.geom = multi
            except (GDALException, GEOSException):
                print(
                    'ERROR: Could not set geom for feature (objectid {}, id {}, counter {})'
                    .format(feat['OBJECTID'], feat['ID'], counter)
                )
                continue

            # Set previous landcover for change patches
            if change:
                patch.nomenclature_previous_id = nomenclature[int(feat.get(code_previous_mapping))]
#.........这里部分代码省略.........

    def test04_wkbwriter(self):
        wkb_w = WKBWriter()

        # Representations of 'POINT (5 23)' in hex -- one normal and
        # the other with the byte order changed.
        g = GEOSGeometry('POINT (5 23)')
        hex1 = b'010100000000000000000014400000000000003740'
        wkb1 = memoryview(binascii.a2b_hex(hex1))
        hex2 = b'000000000140140000000000004037000000000000'
        wkb2 = memoryview(binascii.a2b_hex(hex2))

        self.assertEqual(hex1, wkb_w.write_hex(g))
        self.assertEqual(wkb1, wkb_w.write(g))

        # Ensuring bad byteorders are not accepted.
        for bad_byteorder in (-1, 2, 523, 'foo', None):
            # Equivalent of `wkb_w.byteorder = bad_byteorder`
            self.assertRaises(ValueError, wkb_w._set_byteorder, bad_byteorder)

        # Setting the byteorder to 0 (for Big Endian)
        wkb_w.byteorder = 0
        self.assertEqual(hex2, wkb_w.write_hex(g))
        self.assertEqual(wkb2, wkb_w.write(g))

        # Back to Little Endian
        wkb_w.byteorder = 1

        # Now, trying out the 3D and SRID flags.
        g = GEOSGeometry('POINT (5 23 17)')
        g.srid = 4326

        hex3d = b'0101000080000000000000144000000000000037400000000000003140'
        wkb3d = memoryview(binascii.a2b_hex(hex3d))
        hex3d_srid = b'01010000A0E6100000000000000000144000000000000037400000000000003140'
        wkb3d_srid = memoryview(binascii.a2b_hex(hex3d_srid))

        # Ensuring bad output dimensions are not accepted
        for bad_outdim in (-1, 0, 1, 4, 423, 'foo', None):
            # Equivalent of `wkb_w.outdim = bad_outdim`
            self.assertRaises(ValueError, wkb_w._set_outdim, bad_outdim)

        # Now setting the output dimensions to be 3
        wkb_w.outdim = 3

        self.assertEqual(hex3d, wkb_w.write_hex(g))
        self.assertEqual(wkb3d, wkb_w.write(g))

        # Telling the WKBWriter to include the srid in the representation.
        wkb_w.srid = True
        self.assertEqual(hex3d_srid, wkb_w.write_hex(g))
        self.assertEqual(wkb3d_srid, wkb_w.write(g))

    def test_empty_point_wkb(self):
        p = Point(srid=4326)
        wkb_w = WKBWriter()

        wkb_w.srid = False
        with self.assertRaisesMessage(ValueError, 'Empty point is not representable in WKB.'):
            wkb_w.write(p)
        with self.assertRaisesMessage(ValueError, 'Empty point is not representable in WKB.'):
            wkb_w.write_hex(p)

        wkb_w.srid = True
        for byteorder, hex in enumerate([
            b'0020000001000010E67FF80000000000007FF8000000000000',
            b'0101000020E6100000000000000000F87F000000000000F87F',
        ]):
            wkb_w.byteorder = byteorder
            self.assertEqual(wkb_w.write_hex(p), hex)
            self.assertEqual(GEOSGeometry(wkb_w.write_hex(p)), p)
            self.assertEqual(wkb_w.write(p), memoryview(binascii.a2b_hex(hex)))
            self.assertEqual(GEOSGeometry(wkb_w.write(p)), p)

    def test04_wkbwriter(self):
        wkb_w = WKBWriter()

        # Representations of 'POINT (5 23)' in hex -- one normal and
        # the other with the byte order changed.
        g = GEOSGeometry('POINT (5 23)')
        hex1 = '010100000000000000000014400000000000003740'
        wkb1 = buffer(binascii.a2b_hex(hex1))
        hex2 = '000000000140140000000000004037000000000000'
        wkb2 = buffer(binascii.a2b_hex(hex2))

        self.assertEqual(hex1, wkb_w.write_hex(g))
        self.assertEqual(wkb1, wkb_w.write(g))

        # Ensuring bad byteorders are not accepted.
        for bad_byteorder in (-1, 2, 523, 'foo', None):
            # Equivalent of `wkb_w.byteorder = bad_byteorder`
            self.assertRaises(ValueError, wkb_w._set_byteorder, bad_byteorder)
            
        # Setting the byteorder to 0 (for Big Endian)
        wkb_w.byteorder = 0
        self.assertEqual(hex2, wkb_w.write_hex(g))
        self.assertEqual(wkb2, wkb_w.write(g))

        # Back to Little Endian
        wkb_w.byteorder = 1

        # Now, trying out the 3D and SRID flags.
        g = GEOSGeometry('POINT (5 23 17)')
        g.srid = 4326
        
        hex3d = '0101000080000000000000144000000000000037400000000000003140'
        wkb3d = buffer(binascii.a2b_hex(hex3d))
        hex3d_srid = '01010000A0E6100000000000000000144000000000000037400000000000003140'
        wkb3d_srid = buffer(binascii.a2b_hex(hex3d_srid))

        # Ensuring bad output dimensions are not accepted
        for bad_outdim in (-1, 0, 1, 4, 423, 'foo', None):
            # Equivalent of `wkb_w.outdim = bad_outdim`
            self.assertRaises(ValueError, wkb_w._set_outdim, bad_outdim)

        # These tests will fail on 3.0.0 because of a bug that was fixed in 3.1:
        # http://trac.osgeo.org/geos/ticket/216
        if not geos_version_info()['version'].startswith('3.0.'):
            # Now setting the output dimensions to be 3
            wkb_w.outdim = 3

            self.assertEqual(hex3d, wkb_w.write_hex(g))
            self.assertEqual(wkb3d, wkb_w.write(g))

            # Telling the WKBWriter to inlcude the srid in the representation.
            wkb_w.srid = True
            self.assertEqual(hex3d_srid, wkb_w.write_hex(g))
            self.assertEqual(wkb3d_srid, wkb_w.write(g))

dry_run = len(argv) == 2 and argv[1] == "dry-run"
print "dry_run:", dry_run

print "deleting all places",
for table in ["appfd_alias", "appfd_alternatename", "appfd_place"]:
    call("sudo -u postgres psql -c 'TRUNCATE " + table + " CASCADE;' dbfd", shell=True)
print "done"

delimiter = "\t"
null = "null"

output_file = open("/tmp/allCountriesCleaned.txt", "wb")
writer = csv.writer(output_file, delimiter=delimiter)

wkb_w = WKBWriter()
wkb_w.srid = True
print "writing output file"
with open("/tmp/allCountries.txt", "r") as f:
    counter = 0
    places_to_create = []
    for line in f:
      try:
        counter += 1
        geonameid, name, asciiname, alternatenames, latitude, longitude, feature_class, feature_code, country_code, cc2, admin1_code, admin2_code, admin3_code, admin4_code, population, elevation, dem, timezone, modification_date = line.split("\t")
        if feature_code == "ADM1": admin_level = "1"
        elif feature_code == "ADM2": admin_level = "2"
        elif feature_code == "ADM3": admin_level = "3"
        elif feature_code == "ADM4": admin_level = "4"
        elif feature_code == "ADM5": admin_level = "5"
        else: admin_level = "null"