ค้นหาข้อมูลในบล็อก (Search in this blog)

Friday, September 25, 2009

การแก้ปัญหาระบบสมการเชิงเส้นด้วย Octave และตัวอย่างสำหรับการแปลงแบบสัมพรรค (ของฟรีและดีก็มีในโลก)

วันนี้ขอแนะนำวิธีการแก้ปัญหาระบบสมการแบบใช้เมทริกซ์เข้ามาช่วย (ลืมเรื่องสมการเชิงเส้นและเมทริกซ์กันรึยังครับ แฮะๆ) ตัวอย่างผมมีระบบสมการเชิงเส้นแบบสองตัวแปร (ตัวไม่รู้ค่า 2 ตัว) ดังนี้
2x + y = 19 ----- (1)
x + y = 14 ----- (2)
ถ้าเป็นสมัยมัธยมต้นเราก็จะใช้วิธีย้ายข้างสมการ ดังนี้ครับ
เริ่มต้นจากย้ายข้าง สมการที่ (1) จะได้ y = 19 + 2x แล้วแทนค่า y ในสมการที่ (2) จะได้ค่า x เท่ากับ 5
หลังจากนั้นให้แทนค่า x ในสมการที่ (1) หรือ (2) ก็ได้แล้วจะได้ค่า y = 9
ดังนั้นจากระบบสมการดังกล่าวเราจะได้ค่า x และ y เท่ากับ 5 และ 9 ตามลำดับ

จากตัวอย่างดังกล่าวจะเห็นได้ว่าขั้นตอนไม่ค่อยซับซ้อนเท่าไหร่นะครับ คราวนี้เรามาดูตัวอย่างที่มีการใช้งานจริงในด้านภูมิสารสนเทศนะครับ
สมมุติผมต้องการแปลงระบบพิกัดของข้อมูลภาพ (Image Georeferencing) โดยอาศัยการแปลงแบบสัมพรรคแบบสองมิติ (2D Affine Transformation) ซึ่งก็คือสมการโพลีโนเมียลกำลังหนึ่ง (1st order polynomial) นั่นเอง ก่อนอื่นผมต้องคำนวณก่อนว่าหากใช้สมการโพลีโนเมียลกำลังหนึ่งนั้นต้องรังวัดจุดควบคุมภาคพื้นดินและรังวัดบนภาพทั้งหมดกี่จุด (วิธีการคำนวณดูได้จากกระทู้ก่อนหน้านี้ ที่่นี่ ) เพื่อความสะดวกและง่ายต่อการอธิบายผมขออนุญาตใช้จำนวนสมการเท่ากับจำนวนตัวไม่รู้ค่า (Unknown) หรือแค่ให้พอหาคำตอบได้เท่านั้นแต่ไม่เหมาะกับการทำงานจริง (การทำงานจริงควรรังวัดให้เกินพอนะครับ จำนวนจุดยิ่งเยอะและค่า residual น้อยๆยิ่งดีครับ)
สรุป นะครับ หากต้องการแปลงระบบพิกัดด้วยการแปลงสัมพรรคแบบสองมิติ นั้นต้องใช้ จุดควบคุม 3 จุด (ขั้นต่ำ)
สำหรับครั้งนี้ผมสมมุติว่าผมรังวัดค่าพิกัดภาพ (Image Coordinate) และพิกัดพื้นดิน (Ground Coornidate) ได้ดังนี้ครับ

โดย (ix,iy) คือพิกัดที่วัดได้บนภาพหรือพิกัดภาพนั่นเอง ส่วน (gX,gY) คือค่าพิกัดพื้นดิน (ระบบใดๆก็ได้ ตัวอย่างนี้ผมขอใช้ระบบพิกัดภาพเทียบกับระบบพิกัดกริด UTM 47N) และจากระบบสมการโพลีโนเมียลกำลังหนึ่งหรือสมการการแปลงแบบสัมพรรคซึ่งมีรูปสมการดังนี้

Xi = a (xi) + b (yi) + c
Yi = d (xi) + e (yi) + f

จากสมข้างต้นจะเห็นได้ว่าถ้าเราต้องการทราบพิกัดของจุดภาพใดๆ(จุดใดก็ได้บนภาพ) ว่ามีค่าพิกัดพื้นดิน (UTM 47N) เป็นเท่าไหร่ เราก็สามารถทำได้เพียงใส่ค่าพิกัดภาพ (xi,yi) แทนลงไปในสมการก็จะได้ค่าพิกัดพื้นดิน (Xi,Yi) ทันที
ง่ายไหมครับ ดูผ่านๆหมือนง่ายนะครับ แต่จากสมการข้างต้นสิ่งที่เรายังไม่รู้คือค่าของสัมประสิทธิ์ (Coeficient) ต่างๆ ซึ่งก็คือค่า a,b,c,d,e,f นั่นเอง แล้วเราจะหาได้อย่างไรละครับ (ค่าสัมประสิทธิ์เหล่านี้จะได้จากการหาความสัมพันธ์ระหว่างระบบพิกัดทั้งของโดยใช้ค่าพิกัดของจุดร่วมหรือ Common points ที่เรารังวัดได้จากทั้งระบบพิกัดภาพและระบบพิกัดพื้นดิน)

วิธีการก็คือเราสามารถหาได้โดยอาศัยความสัมพันธ์ของพิกัดจุดภาพและพิกัดพื้นดิน (จากรังวัดจุด) จากตารางตัวอย่างข้างบนครับ ขั้นตอนต่อไป ผมจะทำการแทนค่าพิกัดของจุดร่วมที่รังวัดระบบบนภาพและระบบพิกัพื้นดิน ซึ่งจะได้ระบบสมการใหม่ดังนี้ (วัด 3 จุดจะได้ 6 สมการ)
499271.67=a(8445.75)+b(6611.75)+c ----(1)
2074132.01=d(8445.75)+e(6611.75)+f ----(2)
487786.57=a(3436.19)+b(8011.69)+c ----(3)
2073436.38=d(3436.19)+e(8011.69)+f ----(4)
497860.20=a(8433.13)+b(9457.38)+c ----(5)
2068443.25=d(8433.13)+e(9457.38)+f ----(6)

ระบบสมการข้างต้นจะเห็นว่าจำนวนสมการ (6 สมการ) เท่ากับจำนวนตัวไม่รู้ค่า (6 ตัวคือ a,b,c,d,e,f)
ดังนั้นสามารถหาคำตอบ!!!
แต่............................. ถ้าจะให้มานั่งย้ายข้างสมการไปมาหรือใช้ความรู้ชั้นมัธยมต้นคงจะไม่ไหวแล้วละครับ เพราะว่าตัวเลขมากขึ้นจำนวนสมการก็มากขึ้นตามไปด้วย แล้วเราจะทำยังไงดี ???

วันนี้ผมมีตัวอย่างวิธีการแก้ปัญหาระบบสมการโดยอาศัยเมทริกซ์เข้ามาช่วยมานำเสนอครับ
จากระบบสมการก่อหน้านี้ผมขอจัดระบบสมการใหม่อีกทีและแทนในเมทริกซ์ต่างๆในรูป AX = B ดังนี้ครับ
เมทริกซ์ A บรรจุสัมประสิทธิ์ของตัวไม่รู้ค่า มีขนาด 6x6
เมทริกซ์ X บรรจุตัวไม่รู้ค่า มีขนาด 6x1
เมทริกซ์ B บรรจุค่าพิกัดพื้นดิน มีขนาด 6x1
หลังจากนั้นเราจะใช้คุณสมบัติของเมทริกซ์ ซึ่งไม่สามารถย้ายข้างสมการแบบสมการทั่วๆไปได้ เช่น ย้ายจากฝั่งซ้าย(การคูณ)ไปฝั่งขวา(กลายเป็นการหาร) เป็นต้น
คุณสมบัติของเมทริกซ์นั้นจะย้ายข้างสมการได้แต่เปลี่ยนจากคูณฝั่งซ้ายไปเป็นการคูณด้วยอินเวอร์สที่ฝั่งขวา เช่นจากตัวอย่างข้างต้น AX=B จะเป็น X= Inverse (A)*B ซึ่งความยากคือการหาอินเวอร์สนี่แหละครับ (ซึ่งมีหลายวิธี เช่น Gauss-Jordan เป็นต้น) ครั้งนี้ผมขออนุญาตให้ทุกท่านใช้เครื่องมือฟรี (โปรแกรม) ที่นิยมอย่างมากในแวดวง Open Source ซึ่งคือเจ้า Octave นั่นเอง ซึ่งก็โหลดได้จาก ที่นี่ ครับ (โหลดแล้วติดตั้งกันเลยนะครับ จะได้ทดลองทำตามกันเลย)
หลังจากติดตั้งแล้วให้ทำตามขั้นตอนดังนี้
1) ป้อนค่าเมทริกซ์ A ดังนี้

>>a=[8445.75 6611.75 1 0 0 0;0 0 0 8445.75 6611.75 1;3436.19 8011.69 1 0 0 0;0 0 0 3436.19 8011.69 1;8433.13 9457.38 1 0 0 0;0 0 0 8433.13 9457.38 1]

### สำหรับการป้อนค่าสมาชิกของเมทริกซ์ทนั้นให้ป้อนค่าตัวเลขที่ละแถวซึ่งจะแบ่งค่าของแต่ละคอลัมน์ด้วยการเว้นว่าง (Space) และเมื่อป้อนค่าสมาชิกในแต่ละแถวแล้วการขึ้บรรทัดใหม่นั้นให้ใช้เครื่องหมาย ;

2) ป้อนค่าของเมทริกซ์ B ดังนี้

>>b=[499271.67;2074132.01;487786.57;2073436.38;497860.20;2068443.25]

3) ทำการคำนวณค่าสัมประสิทธ์จากสมการ X=Inverse (A) * (B) โดยพิมพ์คำสั่งในโปรแกรมดังนี้
>>inv (a) * b

4)ผลลัพธ์ที่ได้หมายถึงค่าของสัมประสิทธิ์ของแต่ละตัว (a,b,c,d,e,f) ตามลำดับ

จากตัวอย่างการคำนวณทำให้เราได้ความสัมพันธ์ระหว่างพิกัดภาพและพิกัดพื้นดินดังนี้

Xi = 2.22330667 (xi) - 0.71668202 (yi) + 485888.70006051
Yi = -0.42032228 (xi) - 2.00098553 (yi) + 2090911.96301988

นั่นหมายความว่า จุดภาพที่ (1,1) นั้นค่าพิกัดพื้นดินเท่ากับ (485890.21,2090909.54) และจุดภาพที่ (1,2) มีค่าพิกัดพื้นดินเท่ากับ (485889.49,2090907.54) ... จนถึงจุดภาพสุดท้าย

No comments:

Post a Comment