โฟโตแกรมเมตรีและผลลัพธ์

ทำไมต้องรู้

โปรแกรมทำแผนที่สมัยใหม่กดปุ่มเดียวก็ได้ผลลัพธ์สวยงาม แต่เมื่อแผนที่แยกเป็นสองก้อน กองวัสดุมีรูโหว่ หรือผู้ใช้ขอ “ความสูงพื้นดิน” ใต้ต้นไม้ ผู้ที่ไม่เข้าใจว่าภาพกลายเป็นแผนที่ได้อย่างไรจะแก้ปัญหาไม่ได้ และอาจส่งไฟล์ที่ชื่อถูกแต่ความหมายผิด โมดูลนี้อธิบายกลไกเท่าที่จำเป็นต่อการตัดสินใจ

ความลึกมาจาก parallax

ภาพถ่ายหนึ่งใบบอกทิศทางจากกล้องไปยังวัตถุ แต่ไม่บอกระยะ เมื่อกล้องย้ายตำแหน่ง จุดเดียวกันจะปรากฏในตำแหน่งภาพที่ต่างกัน เรียกว่า parallax ลองยกนิ้วขึ้นหน้าตาแล้วหลับตาทีละข้าง นิ้วที่อยู่ใกล้จะกระโดดมากกว่าวัตถุที่อยู่ไกล นี่คือหลักเดียวกับที่โฟโตแกรมเมตรีใช้

สำหรับภาพคู่ที่ถ่ายจากความสูงเดียวกัน ห่างกัน $B$ (ระยะฐาน) ตำราโฟโตแกรมเมตรีของ Wolf และคณะให้ความสัมพันธ์ $Z=Bf/p$ เมื่อ $p$ คือ parallax บนระนาบภาพ

BASE_M, FOCAL_MM = 16.0, 8.8
PIXEL_MM = 13.2 / 5472

def depth_m(parallax_px):
    return BASE_M * FOCAL_MM / (parallax_px * PIXEL_MM)

ground_px, roof_px = 729.6, 810.7
z_ground, z_roof = depth_m(ground_px), depth_m(roof_px)
print(f"ground {z_ground:.1f} m below camera, roof {z_roof:.1f} m, roof height {z_ground - z_roof:.1f} m")
print(f"1 px parallax error at 80 m changes depth by about {80 ** 2 / (BASE_M * FOCAL_MM) * PIXEL_MM:.2f} m")

ground 80.0 m below camera, roof 72.0 m, roof height 8.0 m
1 px parallax error at 80 m changes depth by about 0.11 m

Structure from Motion

ในงานโดรน เราไม่ได้ใช้ภาพคู่เดียว แต่ใช้ภาพหลายร้อยใบ Structure from Motion (SfM) ประมาณตำแหน่งและทิศของกล้องทุกตัวพร้อมรูปทรงของฉากไปพร้อมกัน วิธีนี้พัฒนาจากงานคอมพิวเตอร์วิทัศน์ เช่น Snavely และคณะ (2006) และถูกนำมาใช้ในงานภูมิศาสตร์อย่างแพร่หลายตามที่ Westoby และคณะ (2012) สรุปไว้

1. ภาพและแบบจำลองกล้อง อ่านขนาดภาพ ทางยาวโฟกัส และความบิดเบี้ยวของเลนส์
2. หาจุดเด่นและจับคู่ เช่น มุมของรอยบนพื้น จุดที่ตามได้ข้ามหลายภาพเรียก tie point คู่ที่ขัดกับเรขาคณิตจะถูกคัดออก
3. SfM และ bundle adjustment OpenSfM เริ่มจากภาพคู่แรก แล้วเพิ่มภาพทีละใบ คำนวณจุดสามมิติ และปรับค่าทั้งหมดร่วมกันเป็นระยะด้วย bundle adjustment เพื่อให้จุดที่ฉายกลับลงภาพตรงกับตำแหน่งที่เห็นจริงมากที่สุด
4. จับคู่หนาแน่น ประมาณความลึกของพิกเซลจำนวนมาก ได้ point cloud หนาแน่น
5. ผูกพิกัดและสร้างผลลัพธ์ ใช้ตำแหน่งกล้องหรือ GCP แล้วสร้าง orthomosaic, DSM หรือ mesh

Reprojection error วัดว่าจุดที่ฉายกลับลงภาพห่างจากตำแหน่งที่เห็นกี่พิกเซล เป็นความสอดคล้องภายใน ไม่ใช่ความคลาดบนพื้นโลก แบบจำลองที่ reprojection error ต่ำยังเลื่อนทั้งผืนได้ถ้าการผูกพิกัดผิด

ผลลัพธ์ห้าชนิด

OpenDroneMap ส่งออก orthophoto และ DEM เป็น GeoTIFF, point cloud เป็น LAZ และ mesh เป็น OBJ ส่วน DSM และ DTM ต้องเปิดตัวเลือกให้สร้าง

กล้อง RGB เห็นเฉพาะสิ่งที่แสงสะท้อนมาถึง ใต้พุ่มไม้ทึบจึงอาจไม่มีภาพพื้นดินเลย ส่วน LiDAR วัดระยะจากเวลาที่แสงเลเซอร์เดินทางไปกลับ และอาจได้จุดพื้นดินผ่านช่องว่างของพุ่มไม้ แต่ต้องมีการสอบเทียบและระบบนำทางที่ดี USGS กำหนดข้อกำหนดของข้อมูล LiDAR ไว้ใน Lidar Base Specification ไฟล์ point cloud ที่หน้าตาเหมือนกันจากสองวิธีจึงมีหลักฐานและข้อจำกัดต่างกัน

ปริมาตรกองวัสดุ

ปริมาตรคำนวณจากกริดความสูง: แต่ละช่องมีพื้นที่ $A_i$ และความสูงผิวลบความสูงฐาน ผลลัพธ์ขึ้นกับ ระดับฐาน ที่เลือกมาก และช่องที่ไม่มีข้อมูลต้องรายงาน ไม่ใช่นับเป็นศูนย์

CELL_AREA = 2.0 * 2.0          # m² ต่อช่อง
dsm = [  # ความสูงผิว (m) None = ไม่มีข้อมูล
    [100.0, 101.0, 102.0, 100.0],
    [101.0, 103.0, 104.0, 101.0],
    [100.0, 102.0, None, 100.0],
]

def volume(base_level):
    cells = [z for row in dsm for z in row if z is not None]
    return sum(CELL_AREA * max(z - base_level, 0.0) for z in cells)

missing = sum(z is None for row in dsm for z in row)
for base in (100.0, 100.5):
    print(f"base {base} m: volume {volume(base):.0f} m3")
print(f"cells without data: {missing} of {sum(len(r) for r in dsm)} -> report, do not treat as zero")

base 100.0 m: volume 56 m3
base 100.5 m: volume 42 m3
cells without data: 1 of 12 -> report, do not treat as zero

ปฏิบัติการประจำโมดูล

ข้อผิดพลาดที่พบบ่อย

สรุป

• ความลึกมาจาก parallax เมื่อกล้องย้ายตำแหน่ง ความสูงไวต่อความคลาดของ parallax
• SfM ประมาณกล้องและรูปทรงร่วมกัน bundle adjustment ลดความต่างเมื่อฉายกลับ reprojection error เป็นความสอดคล้องภายใน
• เลือกผลลัพธ์ให้ตรงคำถาม DSM คือผิวบนสุด DTM คือพื้นดิน
• ปริมาตรขึ้นกับระดับฐาน และต้องเปิดเผยพื้นที่ที่ไม่มีข้อมูล

แบบฝึกตรวจความเข้าใจ

1. $B$ = 20 m, $f$ = 8.8 mm, parallax 2.2 mm ความลึกเท่าใด
2. reprojection error 0.8 พิกเซล แปลว่าแผนที่คลาดไม่เกิน 0.8 พิกเซลบนพื้นใช่หรือไม่
3. ผู้ใช้ต้องการออกแบบทางระบายน้ำในพื้นที่มีต้นไม้ ควรใช้ DSM หรือ DTM และต้องระวังอะไร
4. กริด 5 ช่อง ช่องละ 4 m² สูงเหนือฐาน 1, 1, 2, 2, 3 m ปริมาตรเท่าใด
5. ทำไมช่องที่ไม่มีข้อมูลจึงห้ามนับเป็นศูนย์