Python sanity_check.check_layer函数代码示例

def smart_clip(layer_to_clip, mask_layer, callback=None):
    """Smart clip a vector layer with another.

    Issue https://github.com/inasafe/inasafe/issues/3186

    :param layer_to_clip: The vector layer to clip.
    :type layer_to_clip: QgsVectorLayer

    :param mask_layer: The vector layer to use for clipping.
    :type mask_layer: QgsVectorLayer

    :param callback: A function to all to indicate progress. The function
        should accept params 'current' (int), 'maximum' (int) and 'step' (str).
        Defaults to None.
    :type callback: function

    :return: The clip vector layer.
    :rtype: QgsVectorLayer

    .. versionadded:: 4.0
    """
    output_layer_name = smart_clip_steps['output_layer_name']
    processing_step = smart_clip_steps['step_name']

    writer = create_memory_layer(
        output_layer_name,
        layer_to_clip.geometryType(),
        layer_to_clip.crs(),
        layer_to_clip.fields()
    )
    writer.startEditing()

    # first build up a list of clip geometries
    request = QgsFeatureRequest().setSubsetOfAttributes([])
    iterator = mask_layer.getFeatures(request)
    feature = next(iterator)
    geometries = QgsGeometry(feature.geometry())

    # use prepared geometries for faster intersection tests
    # noinspection PyArgumentList
    engine = QgsGeometry.createGeometryEngine(geometries.geometry())
    engine.prepareGeometry()

    extent = mask_layer.extent()

    for feature in layer_to_clip.getFeatures(QgsFeatureRequest(extent)):

        if engine.intersects(feature.geometry().geometry()):
            out_feat = QgsFeature()
            out_feat.setGeometry(feature.geometry())
            out_feat.setAttributes(feature.attributes())
            writer.addFeature(out_feat)

    writer.commitChanges()

    writer.keywords = layer_to_clip.keywords.copy()
    writer.keywords['title'] = output_layer_name
    check_layer(writer)
    return writer

def union(union_a, union_b):
    """Union of two vector layers.

    Issue https://github.com/inasafe/inasafe/issues/3186

    :param union_a: The vector layer for the union.
    :type union_a: QgsVectorLayer

    :param union_b: The vector layer for the union.
    :type union_b: QgsVectorLayer

    :return: The clip vector layer.
    :rtype: QgsVectorLayer

    .. versionadded:: 4.0
    """
    output_layer_name = union_steps['output_layer_name']
    output_layer_name = output_layer_name % (
        union_a.keywords['layer_purpose'],
        union_b.keywords['layer_purpose']
    )

    keywords_union_1 = union_a.keywords
    keywords_union_2 = union_b.keywords
    inasafe_fields_union_1 = keywords_union_1['inasafe_fields']
    inasafe_fields_union_2 = keywords_union_2['inasafe_fields']
    inasafe_fields = inasafe_fields_union_1
    inasafe_fields.update(inasafe_fields_union_2)

    parameters = {'INPUT': union_a,
                  'OVERLAY': union_b,
                  'OUTPUT': 'memory:'}

    # TODO implement callback through QgsProcessingFeedback object

    initialize_processing()

    feedback = create_processing_feedback()
    context = create_processing_context(feedback=feedback)
    result = processing.run('native:union', parameters, context=context)
    if result is None:
        raise ProcessingInstallationError

    union_layer = result['OUTPUT']
    union_layer.setName(output_layer_name)

    # use to avoid modifying original source
    union_layer.keywords = dict(union_a.keywords)
    union_layer.keywords['inasafe_fields'] = inasafe_fields
    union_layer.keywords['title'] = output_layer_name
    union_layer.keywords['layer_purpose'] = 'aggregate_hazard'
    union_layer.keywords['hazard_keywords'] = keywords_union_1.copy()
    union_layer.keywords['aggregation_keywords'] = keywords_union_2.copy()

    fill_hazard_class(union_layer)

    check_layer(union_layer)
    return union_layer

def reproject(layer, output_crs, callback=None):
    """Reproject a vector layer to a specific CRS.

    Issue https://github.com/inasafe/inasafe/issues/3183

    :param layer: The layer to reproject.
    :type layer: QgsVectorLayer

    :param output_crs: The destination CRS.
    :type output_crs: QgsCoordinateReferenceSystem

    :param callback: A function to all to indicate progress. The function
        should accept params 'current' (int), 'maximum' (int) and 'step' (str).
        Defaults to None.
    :type callback: function

    :return: Reprojected memory layer.
    :rtype: QgsVectorLayer

    .. versionadded:: 4.0
    """
    output_layer_name = reproject_steps['output_layer_name']
    output_layer_name = output_layer_name % layer.keywords['layer_purpose']
    processing_step = reproject_steps['step_name']

    input_crs = layer.crs()
    input_fields = layer.fields()
    feature_count = layer.featureCount()

    reprojected = create_memory_layer(
        output_layer_name, layer.geometryType(), output_crs, input_fields)
    reprojected.startEditing()

    crs_transform = QgsCoordinateTransform(input_crs, output_crs)

    out_feature = QgsFeature()

    for i, feature in enumerate(layer.getFeatures()):
        geom = feature.geometry()
        geom.transform(crs_transform)
        out_feature.setGeometry(geom)
        out_feature.setAttributes(feature.attributes())
        reprojected.addFeature(out_feature)

        if callback:
            callback(current=i, maximum=feature_count, step=processing_step)

    reprojected.commitChanges()

    # We transfer keywords to the output.
    # We don't need to update keywords as the CRS is dynamic.
    reprojected.keywords = layer.keywords
    reprojected.keywords['title'] = output_layer_name
    check_layer(reprojected)
    return reprojected

def intersection(source, mask):
    """Intersect two layers.

    Issue https://github.com/inasafe/inasafe/issues/3186

    :param source: The vector layer to clip.
    :type source: QgsVectorLayer

    :param mask: The vector layer to use for clipping.
    :type mask: QgsVectorLayer

    :return: The clip vector layer.
    :rtype: QgsVectorLayer

    .. versionadded:: 4.0
    """
    output_layer_name = intersection_steps['output_layer_name']
    output_layer_name = output_layer_name % (
        source.keywords['layer_purpose'])

    parameters = {'INPUT': source,
                  'OVERLAY': mask,
                  'OUTPUT': 'memory:'}

    # TODO implement callback through QgsProcessingFeedback object

    initialize_processing()

    feedback = create_processing_feedback()
    context = create_processing_context(feedback=feedback)
    result = processing.run('native:intersection', parameters, context=context)
    if result is None:
        raise ProcessingInstallationError

    intersect = result['OUTPUT']
    intersect.setName(output_layer_name)
    intersect.keywords = dict(source.keywords)
    intersect.keywords['title'] = output_layer_name
    intersect.keywords['layer_purpose'] = \
        layer_purpose_exposure_summary['key']
    intersect.keywords['inasafe_fields'] = \
        dict(source.keywords['inasafe_fields'])
    intersect.keywords['inasafe_fields'].update(
        mask.keywords['inasafe_fields'])
    intersect.keywords['hazard_keywords'] = \
        dict(mask.keywords['hazard_keywords'])
    intersect.keywords['exposure_keywords'] = dict(source.keywords)
    intersect.keywords['aggregation_keywords'] = dict(
        mask.keywords['aggregation_keywords'])

    check_layer(intersect)
    return intersect

def clean_layer(layer):
    """Clean a vector layer.

    :param layer: The vector layer.
    :type layer: QgsVectorLayer

    :return: The buffered vector layer.
    :rtype: QgsVectorLayer
    """
    output_layer_name = clean_geometry_steps['output_layer_name']
    output_layer_name = output_layer_name % layer.keywords['layer_purpose']

    # start editing
    layer.startEditing()
    count = 0

    # iterate through all features
    request = QgsFeatureRequest().setSubsetOfAttributes([])
    for feature in layer.getFeatures(request):
        geom = feature.geometry()
        was_valid, geometry_cleaned = geometry_checker(geom)

        if was_valid:
            # Do nothing if it was valid
            pass
        elif not was_valid and geometry_cleaned:
            # Update the geometry if it was not valid, and clean now
            layer.changeGeometry(feature.id(), geometry_cleaned, True)
        else:
            # Delete if it was not valid and not able to be cleaned
            count += 1
            layer.deleteFeature(feature.id())

    if count:
        LOGGER.critical(
            '%s features have been removed from %s because of invalid '
            'geometries.' % (count, layer.name()))
    else:
        LOGGER.info(
            'No feature has been removed from the layer: %s' % layer.name())

    # save changes
    layer.commitChanges()
    layer.keywords['title'] = output_layer_name

    check_layer(layer)
    return layer

def clean_layer(layer, callback=None):
    """Clean a vector layer.

    :param layer: The vector layer.
    :type layer: QgsVectorLayer

    :param callback: A function to all to indicate progress. The function
        editing should accept params 'current' (int), 'maximum' (int) and
        'step' (str). Defaults to None.
    :type callback: function

    :return: The buffered vector layer.
    :rtype: QgsVectorLayer
    """
    output_layer_name = clean_geometry_steps['output_layer_name']
    processing_step = clean_geometry_steps['step_name']
    output_layer_name = output_layer_name % layer.keywords['layer_purpose']

    # start editing
    layer.startEditing()

    # iterate through all features
    for feature in layer.getFeatures():
        geom = feature.geometry()
        geometry_cleaned = geometry_checker(geom)
        if geometry_cleaned:
            feature.setGeometry(geometry_cleaned)
        else:
            layer.deleteFeature(feature.id())

    # save changes
    layer.commitChanges()

    layer.keywords['title'] = output_layer_name

    check_layer(layer)
    return layer

def clip_by_extent(layer, extent):
    """Clip a raster using a bounding box using processing.

    Issue https://github.com/inasafe/inasafe/issues/3183

    :param layer: The layer to clip.
    :type layer: QgsRasterLayer

    :param extent: The extent.
    :type extent: QgsRectangle

    :return: Clipped layer.
    :rtype: QgsRasterLayer

    .. versionadded:: 4.0
    """
    parameters = dict()
    # noinspection PyBroadException
    try:
        output_layer_name = quick_clip_steps['output_layer_name']
        output_layer_name = output_layer_name % layer.keywords['layer_purpose']

        output_raster = unique_filename(suffix='.tif', dir=temp_dir())

        # We make one pixel size buffer on the extent to cover every pixels.
        # See https://github.com/inasafe/inasafe/issues/3655
        pixel_size_x = layer.rasterUnitsPerPixelX()
        pixel_size_y = layer.rasterUnitsPerPixelY()
        buffer_size = max(pixel_size_x, pixel_size_y)
        extent = extent.buffered(buffer_size)

        if is_raster_y_inverted(layer):
            # The raster is Y inverted. We need to switch Y min and Y max.
            bbox = [
                str(extent.xMinimum()),
                str(extent.xMaximum()),
                str(extent.yMaximum()),
                str(extent.yMinimum())
            ]
        else:
            # The raster is normal.
            bbox = [
                str(extent.xMinimum()),
                str(extent.xMaximum()),
                str(extent.yMinimum()),
                str(extent.yMaximum())
            ]

        # These values are all from the processing algorithm.
        # https://github.com/qgis/QGIS/blob/master/python/plugins/processing/
        # algs/gdal/ClipByExtent.py
        # Please read the file to know these parameters.
        parameters['INPUT'] = layer.source()
        parameters['NO_DATA'] = ''
        parameters['PROJWIN'] = ','.join(bbox)
        parameters['DATA_TYPE'] = 5
        parameters['COMPRESS'] = 4
        parameters['JPEGCOMPRESSION'] = 75
        parameters['ZLEVEL'] = 6
        parameters['PREDICTOR'] = 1
        parameters['TILED'] = False
        parameters['BIGTIFF'] = 0
        parameters['TFW'] = False
        parameters['EXTRA'] = ''
        parameters['OUTPUT'] = output_raster

        initialize_processing()
        feedback = create_processing_feedback()
        context = create_processing_context(feedback=feedback)

        result = processing.run(
            "gdal:cliprasterbyextent",
            parameters,
            context=context)

        if result is None:
            raise ProcessingInstallationError

        clipped = QgsRasterLayer(result['OUTPUT'], output_layer_name)

        # We transfer keywords to the output.
        clipped.keywords = layer.keywords.copy()
        clipped.keywords['title'] = output_layer_name

        check_layer(clipped)

    except Exception as e:
        # This step clip_raster_by_extent was nice to speedup the analysis.
        # As we got an exception because the layer is invalid, we are not going
        # to stop the analysis. We will return the original raster layer.
        # It will take more processing time until we clip the vector layer.
        # Check https://github.com/inasafe/inasafe/issues/4026 why we got some
        # exceptions with this step.
        LOGGER.exception(parameters)
        LOGGER.exception(
            'Error from QGIS clip raster by extent. Please check the QGIS '
            'logs too !')
        LOGGER.info(
            'Even if we got an exception, we are continuing the analysis. The '
            'layer was not clipped.')
#.........这里部分代码省略.........

def reclassify(layer, exposure_key=None):
    """Reclassify a continuous vector layer.

    This function will modify the input.

    For instance if you want to reclassify like this table :
            Original Value     |   Class
            - ∞ < val <= 0     |     1
            0   < val <= 0.5   |     2
            0.5 < val <= 5     |     3
            5   < val <  + ∞   |     6

    You need a dictionary :
        ranges = OrderedDict()
        ranges[1] = [None, 0]
        ranges[2] = [0.0, 0.5]
        ranges[3] = [0.5, 5]
        ranges[6] = [5, None]

    :param layer: The raster layer.
    :type layer: QgsRasterLayer

    :param exposure_key: The exposure key.
    :type exposure_key: str

    :return: The classified vector layer.
    :rtype: QgsVectorLayer

    .. versionadded:: 4.0
    """
    output_layer_name = reclassify_vector_steps['output_layer_name']
    output_layer_name = output_layer_name % layer.keywords['title']

    # This layer should have this keyword, or it's a mistake from the dev.
    inasafe_fields = layer.keywords['inasafe_fields']
    continuous_column = inasafe_fields[hazard_value_field['key']]

    if exposure_key:
        classification_key = active_classification(
            layer.keywords, exposure_key)
        thresholds = active_thresholds_value_maps(layer.keywords, exposure_key)
        layer.keywords['thresholds'] = thresholds
        layer.keywords['classification'] = classification_key
    else:
        classification_key = layer.keywords.get('classification')
        thresholds = layer.keywords.get('thresholds')

    if not thresholds:
        raise InvalidKeywordsForProcessingAlgorithm(
            'thresholds are missing from the layer %s'
            % layer.keywords['layer_purpose'])

    continuous_index = layer.fields().lookupField(continuous_column)

    classified_field = QgsField()
    classified_field.setType(hazard_class_field['type'])
    classified_field.setName(hazard_class_field['field_name'])
    classified_field.setLength(hazard_class_field['length'])
    classified_field.setPrecision(hazard_class_field['precision'])

    layer.startEditing()
    layer.addAttribute(classified_field)

    classified_field_index = layer.fields(). \
        lookupField(classified_field.name())

    for feature in layer.getFeatures():
        attributes = feature.attributes()
        source_value = attributes[continuous_index]
        classified_value = reclassify_value(source_value, thresholds)
        if (classified_value is None
                or (hasattr(classified_value, 'isNull')
                    and classified_value.isNull())):
            layer.deleteFeature(feature.id())
        else:
            layer.changeAttributeValue(
                feature.id(), classified_field_index, classified_value)

    layer.commitChanges()
    layer.updateFields()

    # We transfer keywords to the output.
    inasafe_fields[hazard_class_field['key']] = (
        hazard_class_field['field_name'])

    value_map = {}

    hazard_classes = definition(classification_key)['classes']
    for hazard_class in reversed(hazard_classes):
        value_map[hazard_class['key']] = [hazard_class['value']]

    layer.keywords['value_map'] = value_map
    layer.keywords['title'] = output_layer_name

    check_layer(layer)
    return layer

def from_counts_to_ratios(layer):
    """Transform counts to ratios.

    Formula: ratio = subset count / total count

    :param layer: The vector layer.
    :type layer: QgsVectorLayer

    :return: The layer with new ratios.
    :rtype: QgsVectorLayer

    .. versionadded:: 4.0
    """
    output_layer_name = recompute_counts_steps['output_layer_name']

    exposure = definition(layer.keywords['exposure'])
    inasafe_fields = layer.keywords['inasafe_fields']

    layer.keywords['title'] = output_layer_name

    if not population_count_field['key'] in inasafe_fields:
        # There is not a population count field. Let's skip this layer.
        LOGGER.info(
            'Population count field {population_count_field} is not detected '
            'in the exposure. We will not compute a ratio from this field '
            'because the formula needs Population count field. Formula: '
            'ratio = subset count / total count.'.format(
                population_count_field=population_count_field['key']))
        return layer

    layer.startEditing()

    mapping = {}
    non_compulsory_fields = get_non_compulsory_fields(
        layer_purpose_exposure['key'], exposure['key'])
    for count_field in non_compulsory_fields:
        exists = count_field['key'] in inasafe_fields
        if count_field['key'] in list(count_ratio_mapping.keys()) and exists:
            ratio_field = definition(count_ratio_mapping[count_field['key']])

            field = create_field_from_definition(ratio_field)
            layer.addAttribute(field)
            name = ratio_field['field_name']
            layer.keywords['inasafe_fields'][ratio_field['key']] = name
            mapping[count_field['field_name']
                    ] = layer.fields().lookupField(name)
            LOGGER.info(
                'Count field {count_field} detected in the exposure, we are '
                'going to create a equivalent field {ratio_field} in the '
                'exposure layer.'.format(
                    count_field=count_field['key'],
                    ratio_field=ratio_field['key']))
        else:
            LOGGER.info(
                'Count field {count_field} not detected in the exposure. We '
                'will not compute a ratio from this field.'.format(
                    count_field=count_field['key']))

    if len(mapping) == 0:
        # There is not a subset count field. Let's skip this layer.
        layer.commitChanges()
        return layer

    for feature in layer.getFeatures():
        total_count = feature[inasafe_fields[population_count_field['key']]]

        for count_field, index in list(mapping.items()):
            count = feature[count_field]
            try:
                # For #4669, fix always get 0
                new_value = count / float(total_count)
            except TypeError:
                new_value = ''
            except ZeroDivisionError:
                new_value = 0
            layer.changeAttributeValue(feature.id(), index, new_value)

    layer.commitChanges()
    check_layer(layer)
    return layer

#.........这里部分代码省略.........
        clip_geometries.append(QgsGeometry(mask_feature.geometry()))

    # are we clipping against a single feature? if so,
    # we can show finer progress reports
    if len(clip_geometries) > 1:
        # noinspection PyTypeChecker,PyCallByClass,PyArgumentList
        combined_clip_geom = QgsGeometry.unaryUnion(clip_geometries)
        single_clip_feature = False
    else:
        combined_clip_geom = clip_geometries[0]
        single_clip_feature = True

    # use prepared geometries for faster intersection tests
    # noinspection PyArgumentList
    engine = QgsGeometry.createGeometryEngine(combined_clip_geom.geometry())
    engine.prepareGeometry()

    tested_feature_ids = set()

    for i, clip_geom in enumerate(clip_geometries):
        request = QgsFeatureRequest().setFilterRect(clip_geom.boundingBox())
        input_features = [f for f in layer_to_clip.getFeatures(request)]

        if not input_features:
            continue

        if single_clip_feature:
            total = 100.0 / len(input_features)
        else:
            total = 0

        for current, in_feat in enumerate(input_features):
            if not in_feat.geometry():
                continue

            if in_feat.id() in tested_feature_ids:
                # don't retest a feature we have already checked
                continue

            tested_feature_ids.add(in_feat.id())

            if not engine.intersects(in_feat.geometry().geometry()):
                continue

            if not engine.contains(in_feat.geometry().geometry()):
                cur_geom = in_feat.geometry()
                new_geom = combined_clip_geom.intersection(cur_geom)
                if new_geom.wkbType() == QgsWKBTypes.Unknown \
                        or QgsWKBTypes.flatType(
                            new_geom.geometry().wkbType()) == \
                                QgsWKBTypes.GeometryCollection:
                    int_com = in_feat.geometry().combine(new_geom)
                    int_sym = in_feat.geometry().symDifference(new_geom)
                    if not int_com or not int_sym:
                        # LOGGER.debug(
                        #     tr('GEOS geoprocessing error: One or more input '
                        #        'features have invalid geometry.'))
                        pass
                    else:
                        new_geom = int_com.difference(int_sym)
                        if new_geom.isGeosEmpty()\
                                or not new_geom.isGeosValid():
                            # LOGGER.debug(
                            #     tr('GEOS geoprocessing error: One or more '
                            #        'input features have invalid geometry.'))
                            pass
            else:
                # clip geometry totally contains feature geometry,
                # so no need to perform intersection
                new_geom = in_feat.geometry()

            try:
                out_feat = QgsFeature()
                out_feat.setGeometry(new_geom)
                out_feat.setAttributes(in_feat.attributes())
                if new_geom.type() == layer_to_clip.geometryType():
                    writer.addFeature(out_feat)
            except:
                LOGGER.debug(
                    tr('Feature geometry error: One or more output features '
                       'ignored due to invalid geometry.'))
                continue

            # TODO implement callback
            if single_clip_feature:
                # progress.setPercentage(int(current * total))
                pass

        if not single_clip_feature:
            # coarse progress report for multiple clip geometries
            # progress.setPercentage(100.0 * i / len(clip_geoms))
            pass

    # End copy/paste from Processing plugin.
    writer.commitChanges()

    writer.keywords = layer_to_clip.keywords.copy()
    writer.keywords['title'] = output_layer_name
    check_layer(writer)
    return writer

#.........这里部分代码省略.........

    .. versionadded:: 4.0
    """
    output_layer_name = zonal_stats_steps['output_layer_name']

    exposure = raster.keywords['exposure']
    if raster.crs().authid() != vector.crs().authid():
        layer = reproject(vector, raster.crs())

        # We prepare the copy
        output_layer = create_memory_layer(
            output_layer_name,
            vector.geometryType(),
            vector.crs(),
            vector.fields()
        )
        copy_layer(vector, output_layer)
    else:
        layer = create_memory_layer(
            output_layer_name,
            vector.geometryType(),
            vector.crs(),
            vector.fields()
        )
        copy_layer(vector, layer)

    input_band = layer.keywords.get('active_band', 1)
    analysis = QgsZonalStatistics(
        layer,
        raster,
        'exposure_',
        input_band,
        QgsZonalStatistics.Sum)
    result = analysis.calculateStatistics(None)
    LOGGER.debug(tr('Zonal stats on %s : %s' % (raster.source(), result)))

    output_field = exposure_count_field['field_name'] % exposure
    if raster.crs().authid() != vector.crs().authid():
        output_layer.startEditing()
        field = create_field_from_definition(
            exposure_count_field, exposure)

        output_layer.addAttribute(field)
        new_index = output_layer.fields().lookupField(field.name())
        old_index = layer.fields().lookupField('exposure_sum')
        for feature_input, feature_output in zip(
                layer.getFeatures(), output_layer.getFeatures()):
            output_layer.changeAttributeValue(
                feature_input.id(), new_index, feature_input[old_index])
        output_layer.commitChanges()
        layer = output_layer
    else:
        fields_to_rename = {
            'exposure_sum': output_field
        }
        if qgis_version() >= 21600:
            rename_fields(layer, fields_to_rename)
        else:
            copy_fields(layer, fields_to_rename)
            remove_fields(layer, list(fields_to_rename.keys()))
        layer.commitChanges()

    # The zonal stats is producing some None values. We need to fill these
    # with 0. See issue : #3778
    # We should start a new editing session as previous fields need to be
    # committed first.
    layer.startEditing()
    request = QgsFeatureRequest()
    expression = '\"%s\" is None' % output_field
    request.setFilterExpression(expression)
    request.setFlags(QgsFeatureRequest.NoGeometry)
    index = layer.fields().lookupField(output_field)
    for feature in layer.getFeatures():
        if feature[output_field] is None:
            layer.changeAttributeValue(feature.id(), index, 0)
    layer.commitChanges()

    layer.keywords = raster.keywords.copy()
    layer.keywords['inasafe_fields'] = vector.keywords['inasafe_fields'].copy()
    layer.keywords['inasafe_default_values'] = (
        raster.keywords['inasafe_default_values'].copy())

    key = exposure_count_field['key'] % raster.keywords['exposure']

    # Special case here, one field is the exposure count and the total.
    layer.keywords['inasafe_fields'][key] = output_field
    layer.keywords['inasafe_fields'][total_field['key']] = output_field

    layer.keywords['exposure_keywords'] = raster.keywords.copy()
    layer.keywords['hazard_keywords'] = vector.keywords[
        'hazard_keywords'].copy()
    layer.keywords['aggregation_keywords'] = (
        vector.keywords['aggregation_keywords'])
    layer.keywords['layer_purpose'] = (
        layer_purpose_aggregate_hazard_impacted['key'])

    layer.keywords['title'] = output_layer_name

    check_layer(layer)
    return layer

#.........这里部分代码省略.........
        input_field = summary_rule['input_field']
        case_field = summary_rule['case_field']
        if aggregate_hazard.fields().lookupField(input_field['field_name']) \
                == -1:
            continue
        if aggregate_hazard.fields().lookupField(case_field['field_name']) \
                == -1:
            continue

        summary_value = 0
        for area in aggregate_hazard.getFeatures():
            case_value = area[case_field['field_name']]
            if case_value in summary_rule['case_values']:
                summary_value += area[input_field['field_name']]

        summary_values[key] = summary_value

    for area in analysis.getFeatures(request):
        total = 0
        for i, val in enumerate(unique_hazard):
            if (val == ''
                    or val is None
                    or (hasattr(val, 'isNull')
                        and val.isNull())):
                val = 'NULL'
            sum = flat_table.get_value(hazard_class=val)
            total += sum
            analysis.changeAttributeValue(area.id(), shift + i, sum)

            affected = post_processor_affected_function(
                exposure=exposure,
                hazard=hazard,
                classification=classification,
                hazard_class=val)
            if affected == not_exposed_class['key']:
                not_exposed_sum += sum
            elif affected:
                affected_sum += sum
            else:
                not_affected_sum += sum

        # Total Affected field
        analysis.changeAttributeValue(
            area.id(), shift + len(unique_hazard), affected_sum)

        # Total Not affected field
        analysis.changeAttributeValue(
            area.id(), shift + len(unique_hazard) + 1, not_affected_sum)

        # Total Exposed field
        analysis.changeAttributeValue(
            area.id(), shift + len(unique_hazard) + 2, total - not_exposed_sum)

        # Total Not exposed field
        analysis.changeAttributeValue(
            area.id(), shift + len(unique_hazard) + 3, not_exposed_sum)

        # Total field
        analysis.changeAttributeValue(
            area.id(), shift + len(unique_hazard) + 4, total)

        # Any absolute postprocessors
        for i, field in enumerate(absolute_values.values()):
            value = field[0].get_value(
                all='all'
            )
            analysis.changeAttributeValue(
                area.id(), shift + len(unique_hazard) + 5 + i, value)

        # Summarizer of custom attributes
        for key, summary_value in list(summary_values.items()):
            summary_field = summary_rules[key]['summary_field']
            field = create_field_from_definition(summary_field)
            analysis.addAttribute(field)
            field_index = analysis.fields().lookupField(field.name())
            # noinspection PyTypeChecker
            analysis.keywords['inasafe_fields'][summary_field['key']] = (
                summary_field['field_name'])

            analysis.changeAttributeValue(
                area.id(), field_index, summary_value)

    # Sanity check ± 1 to the result. Disabled for now as it seems ± 1 is not
    # enough. ET 13/02/17
    # total_computed = (
    #     affected_sum + not_affected_sum + not_exposed_sum)
    # if not -1 < (total_computed - total) < 1:
    #     raise ComputationError

    analysis.commitChanges()

    analysis.keywords['title'] = layer_purpose_analysis_impacted['name']
    if qgis_version() >= 21600:
        analysis.setName(analysis.keywords['title'])
    else:
        analysis.setLayerName(analysis.keywords['title'])
    analysis.keywords['layer_purpose'] = layer_purpose_analysis_impacted['key']

    check_layer(analysis)
    return analysis

def polygonize(layer, callback=None):
    """Polygonize a raster layer into a vector layer using GDAL.

    Issue https://github.com/inasafe/inasafe/issues/3183

    :param layer: The layer to reproject.
    :type layer: QgsRasterLayer

    :param callback: A function to all to indicate progress. The function
        should accept params 'current' (int) and 'maximum' (int). Defaults to
        None.
    :type callback: function

    :return: Reprojected memory layer.
    :rtype: QgsRasterLayer

    .. versionadded:: 4.0
    """
    output_layer_name = polygonize_steps['output_layer_name']
    processing_step = polygonize_steps['step_name']
    output_layer_name = output_layer_name % layer.keywords['layer_purpose']
    gdal_layer_name = polygonize_steps['gdal_layer_name']

    if layer.keywords.get('layer_purpose') == 'exposure':
        output_field = exposure_type_field
    else:
        output_field = hazard_value_field

    input_raster = gdal.Open(layer.source(), gdal.GA_ReadOnly)

    srs = osr.SpatialReference()
    srs.ImportFromWkt(input_raster.GetProjectionRef())

    temporary_dir = temp_dir(sub_dir='pre-process')
    out_shapefile = unique_filename(
        suffix='-%s.shp' % output_layer_name, dir=temporary_dir)

    driver = ogr.GetDriverByName("ESRI Shapefile")
    destination = driver.CreateDataSource(out_shapefile)

    output_layer = destination.CreateLayer(gdal_layer_name, srs)

    # We have no other way to use a shapefile. We need only the first 10 chars.
    field_name = output_field['field_name'][0:10]
    fd = ogr.FieldDefn(field_name, ogr.OFTInteger)
    output_layer.CreateField(fd)

    input_band = input_raster.GetRasterBand(1)
    # Fixme : add our own callback to Polygonize
    gdal.Polygonize(input_band, None, output_layer, 0, [], callback=None)
    destination.Destroy()

    vector_layer = QgsVectorLayer(out_shapefile, output_layer_name, 'ogr')

    # Let's remove polygons which were no data
    request = QgsFeatureRequest()
    expression = '"%s" = %s' % (field_name, no_data_value)
    request.setFilterExpression(expression)
    vector_layer.startEditing()
    for feature in vector_layer.getFeatures(request):
        vector_layer.deleteFeature(feature.id())
    vector_layer.commitChanges()

    # We transfer keywords to the output.
    vector_layer.keywords = layer.keywords.copy()
    vector_layer.keywords[
        layer_geometry['key']] = layer_geometry_polygon['key']

    vector_layer.keywords['title'] = output_layer_name
    # We just polygonized the raster layer. inasafe_fields do not exist.
    vector_layer.keywords['inasafe_fields'] = {
        output_field['key']: field_name
    }

    check_layer(vector_layer)
    return vector_layer

#.........这里部分代码省略.........
        Defaults to None.
    :type callback: function

    :return: The classified raster layer.
    :rtype: QgsRasterLayer

    .. versionadded:: 4.0
    """
    output_layer_name = reclassify_raster_steps['output_layer_name']
    processing_step = reclassify_raster_steps['step_name']
    output_layer_name = output_layer_name % layer.keywords['layer_purpose']

    if exposure_key:
        classification_key = active_classification(
            layer.keywords, exposure_key)
        thresholds = active_thresholds_value_maps(layer.keywords, exposure_key)
        layer.keywords['thresholds'] = thresholds
        layer.keywords['classification'] = classification_key
    else:
        classification_key = layer.keywords.get('classification')
        thresholds = layer.keywords.get('thresholds')
    if not thresholds:
        raise InvalidKeywordsForProcessingAlgorithm(
            'thresholds are missing from the layer %s'
            % layer.keywords['layer_purpose'])

    if not classification_key:
        raise InvalidKeywordsForProcessingAlgorithm(
            'classification is missing from the layer %s'
            % layer.keywords['layer_purpose'])

    ranges = {}
    value_map = {}
    hazard_classes = definition(classification_key)['classes']
    for hazard_class in hazard_classes:
        ranges[hazard_class['value']] = thresholds[hazard_class['key']]
        value_map[hazard_class['key']] = [hazard_class['value']]

    if overwrite_input:
        output_raster = layer.source()
    else:
        output_raster = unique_filename(suffix='.tiff', dir=temp_dir())

    driver = gdal.GetDriverByName('GTiff')

    raster_file = gdal.Open(layer.source())
    band = raster_file.GetRasterBand(1)
    no_data = band.GetNoDataValue()
    source = band.ReadAsArray()
    destination = source.copy()

    for value, interval in ranges.iteritems():
        v_min = interval[0]
        v_max = interval[1]

        if v_min is None:
            destination[np.where(source <= v_max)] = value

        if v_max is None:
            destination[np.where(source > v_min)] = value

        if v_min < v_max:
            destination[np.where((v_min < source) & (source <= v_max))] = value

    # Tag no data cells
    destination[np.where(source == no_data)] = no_data_value

    # Create the new file.
    output_file = driver.Create(
        output_raster, raster_file.RasterXSize, raster_file.RasterYSize, 1)
    output_file.GetRasterBand(1).WriteArray(destination)
    output_file.GetRasterBand(1).SetNoDataValue(no_data_value)

    # CRS
    output_file.SetProjection(raster_file.GetProjection())
    output_file.SetGeoTransform(raster_file.GetGeoTransform())
    output_file.FlushCache()

    del output_file

    if not isfile(output_raster):
        raise FileNotFoundError

    reclassified = QgsRasterLayer(output_raster, output_layer_name)

    # We transfer keywords to the output.
    reclassified.keywords = layer.keywords.copy()
    reclassified.keywords['layer_mode'] = 'classified'

    value_map = {}

    hazard_classes = definition(classification_key)['classes']
    for hazard_class in reversed(hazard_classes):
        value_map[hazard_class['key']] = [hazard_class['value']]

    reclassified.keywords['value_map'] = value_map
    reclassified.keywords['title'] = output_layer_name

    check_layer(reclassified)
    return reclassified

#.........这里部分代码省略.........

    aggregate_hazard.startEditing()

    shift = aggregate_hazard.fields().count()
    add_fields(
        aggregate_hazard,
        absolute_values,
        [affected_field, total_field],
        unique_exposure,
        exposure_count_field
    )

    flat_table = FlatTable('aggregation_id', 'hazard_id', 'exposure_class')

    request = QgsFeatureRequest()
    request.setFlags(QgsFeatureRequest.NoGeometry)
    LOGGER.debug('Computing the aggregate hazard summary.')
    for feature in impact.getFeatures(request):
        # Field_index can be equal to 0.
        if field_index is not None:
            value = feature[field_index]
        else:
            value = 1

        aggregation_value = feature[aggregation_id]
        hazard_value = feature[hazard_id]
        if not hazard_value or isinstance(hazard_value, QPyNullVariant):
            hazard_value = not_exposed_class['key']
        exposure_value = feature[exposure_class]
        if not exposure_value or isinstance(exposure_value, QPyNullVariant):
            exposure_value = 'NULL'

        flat_table.add_value(
            value,
            aggregation_id=aggregation_value,
            hazard_id=hazard_value,
            exposure_class=exposure_value
        )

        # We summarize every absolute values.
        for field, field_definition in absolute_values.iteritems():
            value = feature[field]
            if not value or isinstance(value, QPyNullVariant):
                value = 0
            field_definition[0].add_value(
                value,
                aggregation_id=aggregation_value,
                hazard_id=hazard_value
            )

    hazard_keywords = aggregate_hazard.keywords['hazard_keywords']
    classification = hazard_keywords['classification']

    for area in aggregate_hazard.getFeatures(request):
        aggregation_value = area[aggregation_id]
        feature_hazard_id = area[hazard_id]
        if not feature_hazard_id or isinstance(
                feature_hazard_id, QPyNullVariant):
            feature_hazard_id = not_exposed_class['key']
        feature_hazard_value = area[hazard_class]
        total = 0
        for i, val in enumerate(unique_exposure):
            sum = flat_table.get_value(
                aggregation_id=aggregation_value,
                hazard_id=feature_hazard_id,
                exposure_class=val
            )
            total += sum
            aggregate_hazard.changeAttributeValue(area.id(), shift + i, sum)

        affected = post_processor_affected_function(
            classification=classification, hazard_class=feature_hazard_value)
        affected = tr(unicode(affected))
        aggregate_hazard.changeAttributeValue(
            area.id(), shift + len(unique_exposure), affected)

        aggregate_hazard.changeAttributeValue(
            area.id(), shift + len(unique_exposure) + 1, total)

        for i, field in enumerate(absolute_values.itervalues()):
            value = field[0].get_value(
                aggregation_id=aggregation_value,
                hazard_id=feature_hazard_id
            )
            aggregate_hazard.changeAttributeValue(
                area.id(), shift + len(unique_exposure) + 2 + i, value)

    aggregate_hazard.commitChanges()

    aggregate_hazard.keywords['title'] = (
        layer_purpose_aggregate_hazard_impacted['name'])
    if qgis_version() >= 21800:
        aggregate_hazard.setName(aggregate_hazard.keywords['title'])
    else:
        aggregate_hazard.setLayerName(aggregate_hazard.keywords['title'])
    aggregate_hazard.keywords['layer_purpose'] = (
        layer_purpose_aggregate_hazard_impacted['key'])

    check_layer(aggregate_hazard)
    return aggregate_hazard