opencv auto georeferencing scanned map

Question

I have the following sample image:

where I am trying to locate the coords of the four corners of the inner scanned map image like:

i have tried with something like this:

import cv2
import numpy as np

# Load the image
image = cv2.imread('sample.jpg')

# Convert to grayscale
gray = cv2.cvtColor(image, cv2.COLOR_BGR2GRAY)

# Apply Gaussian blur
blurred = cv2.GaussianBlur(gray, (5, 5), 0)

# Perform edge detection
edges = cv2.Canny(blurred, 50, 150)

# Find contours in the edge-detected image
contours, _ = cv2.findContours(edges.copy(), cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)

# Iterate through the contours and filter for squares or rectangles for contour in contours:
perimeter = cv2.arcLength(contour, True)
approx = cv2.approxPolyDP(contour, 0.04 * perimeter, True)

if len(approx) == 4:
    x, y, w, h = cv2.boundingRect(approx)
    aspect_ratio = float(w) / h
    
    # Adjust this threshold as needed
    if aspect_ratio >= 0.9 and aspect_ratio <= 1.1:
        cv2.drawContours(image, [approx], 0, (0, 255, 0), 2)

# Display the image with detected squares/rectangles
cv2.imwrite('detectedy.png', image)

but all i get is something like this:

UPDATE:

so i have found this code that should do what i require:

import cv2
import numpy as np
from subprocess import call

file = "samplebad.jpg"
img = cv2.imread(file)
orig = img.copy()

# sharpen the image (weighted subtract gaussian blur from original)
'''
https://stackoverflow.com/questions/4993082/how-to-sharpen-an-image-in-opencv
larger smoothing kernel = more smoothing
'''
blur = cv2.GaussianBlur(img, (9,9), 0)
sharp = cv2.addWeighted(img, 1.5, blur, -0.5, 0)

# convert the image to grayscale
#       gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
gray = cv2.cvtColor(sharp, cv2.COLOR_BGR2GRAY)

# smooth  whilst keeping edges sharp
'''
(11) Filter size: Large filters (d > 5) are very slow, so it is recommended to use d=5 for real-time applications, and perhaps d=9 for offline applications that need heavy noise filtering.
(17, 17) Sigma values: For simplicity, you can set the 2 sigma values to be the same. If they are small (< 10), the filter will not have much effect, whereas if they are large (> 150), they will have a very strong effect, making the image look "cartoonish".
These values give the best results based upon the sample images
'''
gray = cv2.bilateralFilter(gray, 11, 17, 17)

# detect edges
'''
(100, 200) Any edges with intensity gradient more than maxVal are sure to be edges and those below minVal are sure to be non-edges, so discarded. Those who lie between these two thresholds are classified edges or non-edges based on their connectivity. If they are connected to "sure-edge" pixels, they are considered to be part of edges.
'''
edged = cv2.Canny(gray, 100, 200, apertureSize=3, L2gradient=True)
cv2.imwrite('./edges.jpg', edged)

# dilate edges to make them more prominent
kernel = np.ones((3,3),np.uint8)
edged = cv2.dilate(edged, kernel, iterations=1)
cv2.imwrite('edges2.jpg', edged)

# find contours in the edged image, keep only the largest ones, and initialize our screen contour
cnts, hierarchy = cv2.findContours(edged.copy(), cv2.RETR_TREE, cv2.CHAIN_APPROX_SIMPLE)
cnts = sorted(cnts, key=cv2.contourArea, reverse=True)[:10]
screenCnt = None

# loop over our contours
for c in cnts:

    # approximate the contour
    peri = cv2.arcLength(c, True)
    approx = cv2.approxPolyDP(c, 0.02 * peri, True)

    # if our approximated contour has four points, then we can assume that we have found our screen
    if len(approx) > 0:
        screenCnt = approx
        print(screenCnt)
        cv2.drawContours(img, [screenCnt], -1, (0, 255, 0), 10)
        cv2.imwrite('contours.jpg', img)
        break

# reshaping contour and initialise output rectangle in top-left, top-right, bottom-right and bottom-left order
pts = screenCnt.reshape(4, 2)
rect = np.zeros((4, 2), dtype = "float32")

# the top-left point has the smallest sum whereas the bottom-right has the largest sum
s = pts.sum(axis = 1)
rect[0] = pts[np.argmin(s)]
rect[2] = pts[np.argmax(s)]

# the top-right will have the minumum difference and the bottom-left will have the maximum difference
diff = np.diff(pts, axis = 1)
rect[1] = pts[np.argmin(diff)]
rect[3] = pts[np.argmax(diff)]

# compute the width and height  of our new image
(tl, tr, br, bl) = rect
widthA = np.sqrt(((br[0] - bl[0]) ** 2) + ((br[1] - bl[1]) ** 2))
widthB = np.sqrt(((tr[0] - tl[0]) ** 2) + ((tr[1] - tl[1]) ** 2))
heightA = np.sqrt(((tr[0] - br[0]) ** 2) + ((tr[1] - br[1]) ** 2))
heightB = np.sqrt(((tl[0] - bl[0]) ** 2) + ((tl[1] - bl[1]) ** 2))

# take the maximum of the width and height values to reach our final dimensions
maxWidth = max(int(widthA), int(widthB))
maxHeight = max(int(heightA), int(heightB))

# construct our destination points which will be used to map the screen to a top-down, "birds eye" view
dst = np.array([
    [0, 0],
    [maxWidth - 1, 0],
    [maxWidth - 1, maxHeight - 1],
    [0, maxHeight - 1]], dtype = "float32")

# calculate the perspective transform matrix and warp the perspective to grab the screen
M = cv2.getPerspectiveTransform(rect, dst)
warp = cv2.warpPerspective(orig, M, (maxWidth, maxHeight))

#       cv2.imwrite('./cvCropped/frame/' + file, warp)

# crop border off (85px is empirical)
#       cropBuffer = 85     # this is for the old (phone) images
cropBuffer = 105    # this is for those taken by Nick
height, width = warp.shape[:2]
cropped = warp[cropBuffer:height-cropBuffer, cropBuffer:width-cropBuffer]

# output the result
cv2.imwrite('cropped.jpg', cropped)

but because these old scanned maps have a fold in them it fails and only detects one side like:

is there a way to somehow get opencv to ignore the center region?

Answer 1

Regarding the specific question: I don't think that there is an open cv function to ignore the center region of an image. But you could always change the color of a part of the image. E.g. you could cover the center region with black right after the image preparation and before detecting the contours:

# set color in an area of size -/+ 5% of the image width left and right of the horizontal center of the image
img_width = gray.shape[1]
colored_width = 5 # in percent
greyscale_color = 0  # black
edged[:, int(img_width/2 - img_width*colored_width/100):int(img_width/2 + img_width*colored_width/100)] = greyscale_color

Running contour detection would than result in:

That said, I do not think the fold is the actual issue here. Using your code, I do detect contours on both pages, left and right. I guess, it is more an issue of selecting the right contours to extract the points.

To solve the issue, I would like to suggest another approach that does not require contour approximation, and most of all does not depend on finding one contour that contains all relevant points (the corners).

(*** Updated for high resolution image***)

# setup
bw_threshold = 180  # threshold for creating a binary (black and white image)
padding_size = 10  # used to color around the image margin, in percent of image length and width
colored_center_width = 10 # used to hide fold
k_val = 25   # value to define size of the structuring elements
contour_retrieval_mode = cv2.RETR_LIST  
contour_approx_mode = cv2.CHAIN_APPROX_SIMPLE

# image loading and preparation
orig = cv2.imread(img_path, cv2.IMREAD_UNCHANGED)
img_grey = cv2.cvtColor(orig.copy(), cv2.COLOR_BGR2GRAY)
# converting the image to black and white
_, binary_img = cv2.threshold(img_grey, bw_threshold, 255, cv2.THRESH_BINARY_INV)
# create a black padding along the image border to hide lines where the paper ends
binary_img = create_padding(binary_img, padding_size, greyscale_color=0)
# cover center part of image to hide page fold
binary_img = color_image_center(binary_img, colored_center_width)

# detect vertical contours
k_val = 60
v_kernel = cv2.getStructuringElement(cv2.MORPH_RECT, (1,k_val))
v_img = cv2.morphologyEx(binary_img, cv2.MORPH_OPEN, v_kernel, iterations=2)
contours, _ = cv2.findContours(v_img, contour_retrieval_mode, contour_approx_mode)

# find corners
v_min = np.inf
v_max = - np.inf
h_min = np.inf
h_max = -np.inf
for c in contours:
    if v_min > c[:, 0, 0].min():
        v_min = c[:, 0, 0].min()
    if v_max < c[:, 0, 0].max():
        v_max = c[:, 0, 0].max()
    if h_min > c[:, 0, 1].min():
        h_min = c[:, 0, 1].min()
    if h_max < c[:, 0, 1].max():
        h_max = c[:, 0, 1].max()

# original image with corners
corners = np.array([[[v_min, h_min]], [[v_min, h_max]], [[v_max, h_max]],[[v_max, h_min]]])
corner_img = cv2.drawContours(orig.copy(), corners, -1, (0,255,0), 5)

# helper functions
 
def create_padding(img, padding_size = 10, greyscale_color=0): 
    # margin size in percentage of image width and length
    img_height = img.shape[0]
    img_width = img.shape[1]
    # top
    img[0:int(img_height*padding_size/100)] = greyscale_color
    # bottom
    img[img_height-int(img_height*padding_size/100): img_height] = greyscale_color
    # left
    img[:, 0:int(img_width*padding_size/100)] = greyscale_color
    # right
    img[:, img_width - int(img_width*padding_size/100): img_width] = greyscale_color
    return img

def color_image_center(img, colored_width = 10, greyscale_color=0): 
    # colored with in percentage of image width 
    img_width = img.shape[1]
    img[:, int(img_width/2 - img_width*colored_width/2/100):int(img_width/2 + img_width*colored_width/2/100)] = greyscale_color
    return img

Compared to the original post:

• I covered the center part to avoid detecting the fold in the page
• I increased the kernel size to avoid detecting (the relatively) small vertical lines in the legend (especially on the right side of the image)
• Horizontal contours are not used. Problem here were the long horizontal lines in the lower left legend which would have needed to be identified as irrelevant.
• Additional step to select contours not longer required.