How to use the rasterio.mask.mask function in rasterio

def test_mask_with_extra_padding(basic_image_2x2, basic_image_file, basic_geometry):
    """Output should have 2 extra pixels compared to the standard padded mask"""

    geometries = [basic_geometry]
    with rasterio.open(basic_image_file) as src:
        masked, transform = mask(src, geometries, crop=True, pad=True, pad_width=2)

    assert masked.shape == (1, 7, 7)
    assert np.array_equal(masked[0], basic_image_2x2[0:7, 0:7])

def test_mask_with_maximum_padding(basic_image_2x2, basic_image_file, basic_geometry):
    """Output should not break if too much padding is requested"""
    geometries = [basic_geometry]
    with rasterio.open(basic_image_file) as src:
        masked, transform = mask(src, geometries, crop=True, pad=True, pad_width=10)

    assert masked.shape == (1, 10, 10)
    assert np.array_equal(masked[0], basic_image_2x2[0:10, 0:10])

transform=transform,
            driver='GTiff'
        )
        with rasterio.open(file_path, 'w', **profile) as output:
            output.write(new_data.astype(profile['dtype']))

        bbox = self.bbox

        if bbox is not None:
            bbox = box(*bbox)
            geo = gpd.GeoDataFrame({'geometry': bbox}, index=[0], crs=from_epsg(4326))
            geo = geo.to_crs(crs=profile['crs'])
            bbox = geo.geometry

            with rasterio.open(file_path, 'r') as merged:
                new_data, transform = rasterio.mask.mask(dataset=merged, shapes=bbox, all_touched=True, crop=True)

            # profile.pop('tiled', None)
            profile.update(
                height=new_data.shape[1],
                width=new_data.shape[2],
                transform=transform,
            )
            with rasterio.open(file_path, 'w', **profile) as clipped:
                clipped.write(new_data)

        with rasterio.open(file_path) as f:
            geometry = util.bbox2poly(f.bounds.left, f.bounds.bottom, f.bounds.right, f.bounds.top, as_shapely=True)
        update_metadata(catalog_entry, quest_metadata={'geometry': geometry.to_wkt()})

        return {'datasets': new_dset, 'catalog_entries': catalog_entry}

# test, save to disk
    kwargs = src.meta
    kwargs.update(
        dtype=rasterio.uint8,
        count=1,
        width=burn_boxes_width,
        height = burn_boxes_height,
        transform=new_transform)
    with rasterio.open('test_6_albers.tif', 'w', **kwargs) as dst:
        dst.write_band(1, out_label.astype(rasterio.uint8))

    # mask, get pixels cover by polygons, set all_touched as True
    polygons_json = [mapping(item) for item in adj_polygons]
    out_image, out_transform = mask(src, polygons_json, nodata=0, all_touched=True, crop=True)

    #test: output infomation
    print('out_transform', out_transform)
    print('out_image',out_image.shape)

    # test: save it to disk
    out_meta = src.meta.copy()
    out_meta.update({"driver": "GTiff",
                     "height": out_image.shape[1],
                     "width": out_image.shape[2],
                     "transform": out_transform})   # note that, the saved image have a small offset compared to the original ones (~0.5 pixel)
    save_path = "masked_of_polygon_%d.tif"%(idx+1)
    with rasterio.open(save_path, "w", **out_meta) as dest:
        dest.write(out_image)

>>> from earthpy.io import path_to_example
        >>> # Clip an RGB image to the extent of Rocky Mountain National Park
        >>> rmnp = gpd.read_file(path_to_example("rmnp.shp"))
        >>> with rio.open(path_to_example("rmnp-rgb.tif")) as raster:
        ...     src_image = raster.read()
        ...     out_image, out_meta = es.crop_image(raster, rmnp)
        >>> out_image.shape
        (3, 265, 281)
        >>> src_image.shape
        (3, 373, 485)
    """
    if isinstance(geoms, gpd.geodataframe.GeoDataFrame):
        clip_extent = [extent_to_json(geoms)]
    else:
        clip_extent = geoms
    out_image, out_transform = mask(
        raster, clip_extent, crop=True, all_touched=all_touched
    )
    out_meta = raster.meta.copy()
    out_meta.update(
        {
            "driver": "GTiff",
            "height": out_image.shape[1],
            "width": out_image.shape[2],
            "transform": out_transform,
        }
    )
    return out_image, out_meta

progress.setPercentage(int(i * total))
                    lotId = feat['properties'][lotIdFieldName]

                    # Skip features with invalid geometries
                    if not feat['geometry']:
                        continue

                    # Skip features that are not inside rasters bounds
                    poly = shape(feat['geometry'])
                    if not bbox.contains(poly):
                        continue

                    # Calculate change percentage
                    try:
                        cdImg, _ = rasterio.mask.mask(cdDs, [feat['geometry']], crop=True)
                        img, _ = rasterio.mask.mask(imgDs, [feat['geometry']], crop=True)
                    except ValueError as err:
                        progress.setText(
                            self.tr('Error on lot id {}: {}. Skipping').format(lotId, err))
                        continue

                    # Skip features with no pixels in raster (too low resolution?)
                    totalPixels = np.sum(img[0] > 0)
                    if totalPixels == 0:
                        progress.setText(
                            self.tr('Lot {} has no pixels? Skipping...').format(lotId))
                        continue

                    count = np.sum(cdImg[0] > 0)
                    perc = count / float(totalPixels)
                    changeDetected = perc >= selectionThreshold

bottom = max(prediction_dataset.bounds.bottom,test_labels_dataset.bounds.bottom)
            right = min(prediction_dataset.bounds.right,test_labels_dataset.bounds.right)
            top = min(prediction_dataset.bounds.top,test_labels_dataset.bounds.top)
            
            common_bbox = {
                        "type": "Polygon",
                        "coordinates": [[
                            [left, bottom],
                            [left, top],
                            [right, top],
                            [right, bottom],
                            [left, bottom]]]}
                        
            # Read data from only the overlapping area
            y_pred, transform = rasterio.mask.mask(prediction_dataset, common_bbox, crop=True)
            y_true, transform = rasterio.mask.mask(test_labels_dataset, common_bbox, crop=True)

            ### Let's print the shapes of both images. They need to be the same
            print("Shape of the Prediction image clip: ",  y_pred.shape)
            print("Shape of the Validation image clip: ",  y_true.shape)

            ### Calculate the f1, precision and recall with Solaris
            ### prop_threshold determines what value is the dividing number on the predicted image. Below that will get 0, over it 1 so change it to your liking
            f1, precision, recall = sol.eval.pixel.f1(y_true,y_pred,prop_threshold=0.001,verbose=True)
            #print("F1 score: {}, Precision: {}, Recall: {}".format(f1,precision,recall))

left = max(prediction_dataset.bounds.left,test_labels_dataset.bounds.left)
            bottom = max(prediction_dataset.bounds.bottom,test_labels_dataset.bounds.bottom)
            right = min(prediction_dataset.bounds.right,test_labels_dataset.bounds.right)
            top = min(prediction_dataset.bounds.top,test_labels_dataset.bounds.top)
            
            common_bbox = {
                        "type": "Polygon",
                        "coordinates": [[
                            [left, bottom],
                            [left, top],
                            [right, top],
                            [right, bottom],
                            [left, bottom]]]}
                        
            # Read data from only the overlapping area
            y_pred, transform = rasterio.mask.mask(prediction_dataset, common_bbox, crop=True)
            y_true, transform = rasterio.mask.mask(test_labels_dataset, common_bbox, crop=True)
            
            # Reshape data for scikit-learn
            y_pred2 = y_pred.reshape(-1)
            y_true2 = y_true.reshape(-1)
            
            # If results of binary classification, reclassify the data based on the treshold.
            if no_of_classes == 2: 
                y_pred2[(y_pred2 >= prediction_treshold)] = 1
                y_pred2[(y_pred2 < prediction_treshold)] = 0
                y_pred2 = y_pred2.astype('int')
                        
            print('Confusion Matrix')    
            print(confusion_matrix(y_true2, y_pred2))
            print('Classification Report')
            print(classification_report(y_true2, y_pred2, zero_division=0))

with rasterio.open(save_path, "w", **out_meta) as dest:
                dest.write(out_image)
        pass
    else:
        # for the case it overlap more than one raster, need to produce a mosaic
        tmp_saved_files = []

        for k_img,image_path in enumerate(image_list):
            with rasterio.open(image_path) as src:
                polygon_json = mapping(selected_polygon)
                if brectangle:
                    # polygon_box = selected_polygon.bounds
                    polygon_json = mapping(selected_polygon.envelope)  # shapely.geometry.Polygon([polygon_box])

                # crop image and saved to disk
                out_image, out_transform = mask(src, [polygon_json], nodata=dstnodata, all_touched=True, crop=True)

                tmp_saved = os.path.splitext(save_path)[0] +'_%d'%k_img + os.path.splitext(save_path)[1]
                # test: save it to disk
                out_meta = src.meta.copy()
                out_meta.update({"driver": "GTiff",
                                 "height": out_image.shape[1],
                                 "width": out_image.shape[2],
                                 "transform": out_transform})  # note that, the saved image have a small offset compared to the original ones (~0.5 pixel)
                with rasterio.open(tmp_saved, "w", **out_meta) as dest:
                    dest.write(out_image)
                tmp_saved_files.append(tmp_saved)

        # mosaic files in tmp_saved_files
        mosaic_args_list = ['gdal_merge.py', '-o', save_path,'-n',str(dstnodata),'-a_nodata',str(dstnodata)]
        mosaic_args_list.extend(tmp_saved_files)
        if basic.exec_command_args_list_one_file(mosaic_args_list,save_path) is False:

# Create a bounding box with Shapely
    bbox = box(minx, miny, maxx, maxy)

    # Insert the bbox into a GeoDataFrame
    geo = gpd.GeoDataFrame({'geometry': bbox}, index=[0]) #, crs=gpkg_crs['init'])

    # Re-project into the same coordinate system as the raster data
    # geo = geo.to_crs(crs=raster.crs.data)

    # Get the geometry coordinates by using the function.
    coords = getFeatures(geo)

    # clip the raster with the polygon
    try:
        out_img, out_transform = mask(dataset=raster, shapes=coords, crop=True)
    except ValueError as e:  # if gpkg's extent outside raster: "ValueError: Input shapes do not overlap raster."
        # TODO: warning or exception? if warning, except must be set in images_to_samples
        raise

    out_meta = raster.meta.copy()
    out_meta.update({"driver": "GTiff",
                     "height": out_img.shape[1],
                     "width": out_img.shape[2],
                     "transform": out_transform})
    out_tif = Path(raster.name).parent / f"{Path(raster.name).stem}_clipped{Path(raster.name).suffix}"
    dest = rasterio.open(out_tif, "w", **out_meta)
    if debug:
        print(f"DEBUG: writing clipped raster to {out_tif}")
        dest.write(out_img)

    return out_img, dest