Luận văn Nghiên cứu lập trình tính biến dạng xoay trong một tam giác địa động lực

1 TỔNG QUAN VỀ XÂY DỰNG CHƯƠNG TRÌNH Để xây dựng chương trình tính biến dạng xoay trong một tam giác địa động lực này trước tiên ta cần xác định các dữ liệu đầu vào và dữ liệu đầu ra. Ở đây dữ liệu đầu vào của chương trình là file số liệu là một file Excel có chứa các dữ liệu về các điểm đo tốc độ dịch chuyển trong đó tại mỗi điểm có chứa các thông số tọa độ, vận tốc dịch chuyển theo phương Bắc và Đông, sai số vận tốc dịch chuyển (Hình 3.3). Và thêm một file đầu vào nữa là file Excel chứa dữ liệu chọn tam giác (Hình 3.4), file chọn tam giác này do người sử dụng tự chọn để tính biến dạng cho những tam giác trong lưới đo một cách trực quan nhất. Trong hình 3.3 là một ví dụ về file số liệu đầu vào để tính biến dạng, file này lấy số liệu từ bảng “ Vận tốc chuyển dịch tuyệt đối của các điểm đo GPS trong khu vực nghiên cứu tính biến dạng trong tài liệu luận văn tiến sĩ “BIẾN DẠNG KIẾN TẠO HIỆN ĐẠI KHU VỰC BIỂN ĐÔNG VIỆT NAM THEO SỐ LIỆU DỊCH CHUYỂN GPS” của nghiên cứu sinh Nguyễn Văn Hướng tại Phòng Địa Động Lực – Viện Địa Chất. Trong hình 3.4 là file chọn tam giác địa động lực cho tính biến dạng, File này gồm các trường đó là tên tam giác (tên này sẽ được thể hiện tại trọng tâm của tam giác) và 3 trường là số ký hiệu của 3 đỉnh tam giác, mỗi trường là một ký hiệu cho một điểm của tam giác được chọn. Ký hiệu này là số thứ tự trong file Excel của file đầu vào cho tính toán biến dạng. Các con số ký hiệu ở 3 trường này là không được trùng nhau.Ví dụ như tại file chọn tam giác, tam giác có tên TG_01 sẽ là 3 điểm BLV1, DOHO,KUNM trong file số liệu đầu vào tạo thành. 2,4,7 là ký hiệu trong 3 trường của file chọn tam giác là số thứ tự của 3 điểm BLV1, DOHO, KUNM trong file số liệu đầu vào. Tên điểm đo Kinh độ Vỹ độ Vận tốc chuyển dịch theo hướng Đông Vận tốc chuyển dịch theo hướng Bắc Sai số Vận tốc chuyển dịch theo hướng Đông,Bắc Hình 2.3: File đầu vào theo số liệu GPS Biển Đông Tên tam giác chọn 3 đỉnh của tam giác được ký hiệu là số thứ tự của điểm đo trong file đầu vào Hình 2.4: File chọn tam giác địa động lực Để tính biến dạng xoay trong một tam giác địa động lực theo công thức (1.14) thì ta phải làm tuần tự các bước như trong mục 1.1.2.2(Chương 1). Trong xây dựng chương trình thì đầu tiên chương trình cần phải đọc hiểu được dữ liệu trong file đầu vào (Hình 2.3) và kết quả đầu ra là một file kết quả cho ta biết được giá trị biến dạng xoay và sai số của nó trong từng tam giác. Kết quả này được ghi vào một file có định dạng “.txt “. Các thuật toán sẽ được chỉ rõ ở mục 2.2.2 (Chương 2) và cách sử dụng và giao diện của chương trình được chỉ rõ tại chương 3 của đồ án. 2.2.2 THUẬT TOÁN TRONG XÂY DỰNG CHƯƠNG TRÌNH Chương trình được xây dựng sử dụng cho tính biến dạng xoay trong một tam giác địa động lực được thực hiện bằng thuật toán tổng hợp như sau (Hình 2.5) : Yes Bắt đầu Nhập dữ liệu từ file Excel Yes/No Hiển thị lưới tốc độ chuyển dịch No Hiển thị lưới điểm đo Nhập dữ liệu chọn tam giác Tính ω=12∂vx∂y-∂vy∂x Tính sai số Xuất file dữ liệu kết quả Hiển thị lưới tam giác Hiển thị ký hiệu biến dạng xoay Kết thúc Hình 2.5: Minh họa thuật toán chương trình Dưới đây là những thuật toán chi tiết trong chương trình : 2.2.2.1 Thuật toán lấy dữ liệu và hiển thị điểm lưới đo Thuật toán này có thể mô tả như sau : Trước tiên ta có một khai báo một biến có cấu trúc DI() có dạng cấu trúc tên là Diem có cách khai báo như sau: Type Diem tendiem as string x as float y as float vx as float vy as float svx as float svy as float End Type Global DI() as Diem Trong này có chứa các biến x,y,vx,vy,svx,svy có nghĩa lần lượt là các biến lưu kinh độ, vĩ độ, vận tốc dịch chuyển theo hướng Đông, vận tốc chuyển dịch theo hướng Bắc bản đồ, sai số vận tốc dịch chuyển theo hướng Đông, Bắc bản đồ.Biến DI() này sẽ lưu dữ liệu từ file Excel số liệu có chứa nhưng thông tin trên (Hình 2.3), dữ liệu được lấy lần lượt theo hàng và DI(i) với i tương ứng là số hàng lưu trên DI() cho tới khi kết thúc file số liệu thì thôi: Sub HienDiem Filename=FileOpenDlg("","","XLS","Hay nhap file so lieu cho tinh bien dang") Register Table Filename Type "XLS" Into "C:\BienDang.Tab" Dim awin as integer Dim tenlop as String awin=frontwindow() tenlop=LayerInfo(awin,1,LAYER_INFO_NAME) if tenlop"Cosmetic1" then Set Coordsys table tenlop Set Map Layer 1 Editable on else tenlop=LayerInfo(1,2,LAYER_INFO_NAME) Set Coordsys table tenlop Set Map Layer 2 Editable on End if Dim i as integer i=1 Redim DI(i) Open table "C:\BienDang.Tab" Fetch First From BienDang DI(i).tendiem = BienDang.col1 DI(i).x= BienDang.col2 DI(i).y= BienDang.col3 DI(i).vx= BienDang.col4 DI(i).vy= BienDang.col5 DI(i).svx= BienDang.col6 DI(i).svy= BienDang.col7 Do while not EOT(BienDang) DI(i).tendiem = BienDang.col1 DI(i).x= BienDang.col2 DI(i).y= BienDang.col3 DI(i).vx= BienDang.col4 DI(i).vy= BienDang.col5 DI(i).svx= BienDang.col6 DI(i).svy= BienDang.col7 i=i + 1 Redim DI(i) Fetch Next From BienDang Loop Call Thehiendiem(DI()) Close Table BienDang Alter Menu Item BienDang Enable Alter Menu Item HienDiem Disable End Sub Phía trên là dòng code khai báo của hàm có tên là Hienthi trong hàm này có dòng lệnh Call có tác dụng gọi một hàm con là Thehiendiem(D()). Hàm này có tác dụng lấy dữ liệu từ biến D() vừa đã được chứa dữ liệu từ file Excel (Hình 2.3) để hiển thị ra bản đồ của Mapinfo những thông số của dữ liệu được lấy là vị trí của điểm, tên điểm,và khi chạy hết hàm Hienthi sẽ cho ra bản đồ Mapinfo thông tin của lưới các điểm như vậy : Sub Thehiendiem(DI() as Diem) Dim i as integer Dim tl as float Dim hoi as logical tl=5/50 Dim ten as string For i=1 to (ubound(DI)-1) Create Point (DI(i).x,DI(i).y) Symbol (34,Black,5) ten = DI(i).tendiem Create text ten (DI(i).x, DI(i).y) (DI(i).x + 0.2,DI(i).y + 0.2) Font MakeFont ("Helvetica",1,12,BLUE,WHITE) Next Hoi = Ask(" Ban co muon the hien toc do dich chuyen","Co","Khong") If Hoi = True then For i=1 to (ubound(DI)-1) Call Tocdodichchuyen(DI(i).x, DI(i).y,DI(i).vx, DI(i).vy,DI(i).svx, DI(i).svy,tl) Next End If End Sub 2.2.2.2 Thuật toán hiển thị tốc độ dịch chuyển Trong phần trên tại hàm Thehiendiem(DI() as Diem) có lệnh gọi tới hàm Tocdodichchuyen(DI(i).x, DI(i).y,DI(i).vx, DI(i).vy,DI(i).svx, DI(i).svy,tl) thì đây là một hàm thể hiện tốc độ chuyển dịch của tại một điểm đo có chứa các thông số về thông tin kinh độ x, vỹ độ y, vận tốc chuyển dịch theo hướng Đông, Bắc là vx,vy và sai số vận tốc chuyển dịch theo hướng Đông, Bắc là svx,svy, tỷ lệ bản đồ tl : Sub Tocdodichchuyen(x as float,y as float,vx as float, vy as float, svx as float, svy as float,tle as float) Dim x1,y1,x2,y2,x3,y3 as float x1= x + tle*vx y1= y + tle*vy x2= x1 + svx*tle*3 y2= y1 - svy*tle*3 x3= x1 - svx*tle*3 y3= y1 + svy*tle*3 Create Ellipse (x2,y2) (x3,y3) Pen MakePen(1,2,GREEN) Brush MakeBrush(01,Black,Black) Create Line (x,y) (x1,y1) Pen MakePen(1,59,Black) End Sub 2.2.2.3 Thuật toán lấy dữ liệu chọn tam giác và hiển thị lưới tam giác Sau khi ta có phần hiển thị và lấy dữ liệu của lưới điểm đo thì công việc tiếp theo sẽ là việc chọn tam giác địa động lực. Qua phần hiển thị điểm ta có cách nhìn trực quan từ đó chọn các tam giác vào file Excel như hình 2.4. Dữ liệu tam giác này sẽ được lưu vào một biến có cấu trúc tên là Tamgiac: Type Tamgiac Tendiem as string d1 as integer d2 as integer d3 as integer End Type Trong đó biến Tendiem là kiểu chuỗi lưu tên của tam giác được chọn, d1,d2,d3 là 3 biến lưu tương ứng 3 điểm trên một tam giác được chọn. Dưới đây là code của hàm có tác dụng lấy dữ liệu tam giác trong file dữ liệu chọn tam giác: Sub Laydulieutamgiac(TGD() as Tamgiac) Filename=FileOpenDlg("","","XLS","Hay nhap file so lieu cho chon tam giac") Register Table Filename Type "XLS" Into "C:\Chontamgiac.Tab" Dim i as integer i=1 Redim TGD(i) Open table "C:\Chontamgiac.Tab" Fetch First From Chontamgiac TGD(i).Tendiem = Chontamgiac.col1 TGD(i).d1 = Chontamgiac.col2 TGD(i).d2 = Chontamgiac.col3 TGD(i).d3 = Chontamgiac.col4 Do while not EOT(Chontamgiac) TGD(i).Tendiem = Chontamgiac.col1 TGD(i).d1 = Chontamgiac.col2 TGD(i).d2 = Chontamgiac.col3 TGD(i).d3 = Chontamgiac.col4 i=i + 1 Redim TGD(i) Fetch Next From Chontamgiac Loop Close Table Chontamgiac End Sub 2.2.2.4 Thuật toán nhân ma trận bậc 4 Do Mapbasic không hỗ trợ mảng hai chiều do đó hạn chể trong việc tính toán trong ma trận.Trong chương trình này ta chỉ cần tính toán đối với các ma trận bậc 4.Do đặc điểm này mà ta sẽ khai báo một kiểu cấu trúc ma trận bậc 4 như sau : Type Matranb4 h1 as float h2 as float h3 as float h4 as float End Type Trong đó h1,h2,h3,h4 được ví như các hàng của ma trận, khai báo một biến chuỗi kiểu Matranb4 ta sẽ có được một biến dữ liệu ma trận 4 hàng và nhiều cột. Hàm có tên NhanMT(MT1() as Matranb4, MT2() as Matranb4, MTKQ() as Matranb4) có tác dụng nhân 2 ma trận bậc 4 đầu vào là MT1(), MT2() và cho ra một kết quả là MTKQ(). Code của hàm này như sau : Sub NhanMT(MT11() as Matranb4, MT22() as Matranb4, MTN() as Matranb4) Dim k1,k2,k3,k4 as float Dim i,s as integer s = ubound(MT22) Redim MTN(s) For i = 1 to ubound(MT22) k1 = MT11(1).h1 * MT22(i).h1 k2 = MT11(2).h1 * MT22(i).h2 k3 = MT11(3).h1 * MT22(i).h3 k4 = MT11(4).h1 * MT22(i).h4 MTN(i).h1 = k1+k2+k3+k4 Next For i = 1 to ubound(MT22) k1 = MT11(1).h2 * MT22(i).h1 k2 = MT11(2).h2 * MT22(i).h2 k3 = MT11(3).h2 * MT22(i).h3 k4 = MT11(4).h2 * MT22(i).h4 MTN(i).h2 = k1+k2+k3+k4 Next For i = 1 to ubound(MT22) k1 = MT11(1).h3 * MT22(i).h1 k2 = MT11(2).h3 * MT22(i).h2 k3 = MT11(3).h3 * MT22(i).h3 k4 = MT11(4).h3 * MT22(i).h4 MTN(i).h3 = k1+k2+k3+k4 Next For i = 1 to ubound(MT22) k1 = MT11(1).h4 * MT22(i).h1 k2 = MT11(2).h4 * MT22(i).h2 k3 = MT11(3).h4 * MT22(i).h3 k4 = MT11(4).h4 * MT22(i).h4 MTN(i).h4 = k1+k2+k3+k4 Next End Sub 2.2.2.5 Thuật toán chuyển vị ma trận bậc 4 Đây là một hàm có tác dụng với đầu vào là một ma trận bậc 4 sẽ cho ra kết quả là ma trận chuyển vị của nó có dạng ma trận bậc 4 theo phép chuyển đổi ma trận chuyển vị ở chương 2 phần ma trận chuyển vị. Hàm này được khai báo có tên Chuyenvi(MTBD() as Matranb4, MTCV() as Matranb4) Trong đó MTBD() là biến có kiểu cấu trúc Matranb4 là ma trận dữ liệu đầu vào, và kết quả của chuyển vị này là ma trận dữ liệu được lưu vào biến MTCV() có cấu trúc Matranb4. Hàm này có code như sau : Sub Chuyenvi(MTBD() as Matranb4, MTCV() as Matranb4) Dim a1,a2,a3,a4 as float Dim i as integer Redim MTCV(4) For i = 1 to 4 a1=MTBD(i).h1 a2=MTBD(i).h2 a3=MTBD(i).h3 a4=MTBD(i).h4 If i=1 then MTCV(1).h1 = a1 MTCV(2).h1 = a2 MTCV(3).h1 = a3 MTCV(4).h1 = a4 ElseIf i=2 then MTCV(1).h2 = a1 MTCV(2).h2 = a2 MTCV(3).h2 = a3 MTCV(4).h2 = a4 ElseIf i=3 then MTCV(1).h3 = a1 MTCV(2).h3 = a2 MTCV(3).h3 = a3 MTCV(4).h3 = a4 ElseIf i=4 then MTCV(1).h4 = a1 MTCV(2).h4 = a2 MTCV(3).h4 = a3 MTCV(4).h4 = a4 End If Next End Sub 2.2.2.6 Thuật toán tính định thức ma trận bậc 4 Để tính định thức của ma trận bậc 4 theo công thức (1.22) trước tiên ta tạo một hàm tính định thức bậc 3 rồi tính định thức bậc 4 theo như phần trình bày tại mục 1.2.3. Hàm tính định thức bậc 4 được khai báo có tên như sau: Det44(MT4() as Matranb4,det as float) và hàm con của nó tính định thức cho ma trận bậc 3 có tên được khai báo là Det33(MT3() as Matranb4,det as float). Sub Det44(MT4() as Matranb4,det as float) Dim MT33() as Matranb4 Dim dt,b1,b2,b3,k1,k2,k3,k4 as float Dim i as integer Redim MT33(3) For i=1 to 3 b1 = MT4(i+1).h2 b2 = MT4(i+1).h3 b3 = MT4(i+1).h4 MT33(i).h1 = b1 MT33(i).h2 = b2 MT33(i).h3 = b3 Next Call Det33(MT33(),dt) k1 = dt* MT4(1).h1 For i=1 to 4 If i = 1 then b1 = MT4(i).h2 b2 = MT4(i).h3 b3 = MT4(i).h4 MT33(i).h1 = b1 MT33(i).h2 = b2 MT33(i).h3 = b3 ElseIf i>2 then b1 = MT4(i).h2 b2 = MT4(i).h3 b3 = MT4(i).h4 MT33(i-1).h1 = b1 MT33(i-1).h2 = b2 MT33(i-1).h3 = b3 End If Next Call Det33(MT33(),dt) k2 = dt* MT4(2).h1 For i=1 to 4 If i=1 or i=2 then b1 = MT4(i).h2 b2 = MT4(i).h3 b3 = MT4(i).h4 MT33(i).h1 = b1 MT33(i).h2 = b2 MT33(i).h3 = b3 ElseIf i>3 then b1 = MT4(i).h2 b2 = MT4(i).h3 b3 = MT4(i).h4 MT33(i-1).h1 = b1 MT33(i-1).h2 = b2 MT33(i-1).h3 = b3 End If Next Call Det33(MT33(),dt) k3 = dt* MT4(3).h1 For i=1 to 3 b1 = MT4(i).h2 b2 = MT4(i).h3 b3 = MT4(i).h4 MT33(i).h1 = b1 MT33(i).h2 = b2 MT33(i).h3 = b3 Next Call Det33(MT33(),dt) k4 = dt* MT4(4).h1 det = k1 - k2 + k3 -k4 End Sub '....................................................................... Sub Det33(MT3() as Matranb4,det as float) Dim a1,a2,a3,a4,a5,a6 as float a1= MT3(1).h1 * MT3(2).h2 * MT3(3).h3 a2= MT3(1).h1 * MT3(3).h2 * MT3(2).h3 a3= MT3(2).h1 * MT3(1).h2 * MT3(3).h3 a4= MT3(2).h1 * MT3(3).h2 * MT3(1).h3 a5= MT3(3).h1 * MT3(1).h2 * MT3(2).h3 a6= MT3(3).h1 * MT3(2).h2 * MT3(1).h3 det = a1 - a2 - a3 + a4 + a5 -a6 End Sub Trong đó hàm Det33 xử lý tính định thức ma trận đầu vào là MT3() cho kết quả định thức của ma trận này là det. Hàm Det44 sẽ tính định thức cho ma trận đầu vào MT4() và kết quả là det. 2.2.2.7 Thuật toán tính ma trận nghịch đảo bậc 4 Do việc Mapbasic không hỗ trợ mảng 2 chiều do vậy việc tính ma trận nghịch đảo tương đối dài dòng. Hàm có tên được khai báo là Ngichdao(MT4() as Matranb4, MTND() as Matranb4) có tác dụng tính ma trận nghịch đảo của ma trận đầu vào được lưu tại biến tên MT4() và kết quả của nó được lưu tại biến MTND(), cả hai biến này đều có dạng biến cấu trúc của Matranb4. Cách thức tính toán theo phần lý thuyết tính ma trận nghịch đảo tại chương 1. Dưới đây là code của hàm tính ma trận nghịch đảo bậc 4 : Sub Ngichdao(MT4() as Matranb4, MTND() as Matranb4) Dim MT33(),ND() as Matranb4 Dim k,dt,b1,b2,b3 as float Dim i as integer Redim MT33(3) Redim MT4(4) Redim ND(4) Redim MTND(4) '1111111111111111111111111111111111111111111111 For i=1 to 3 b1 = MT4(i+1).h2 b2 = MT4(i+1).h3 b3 = MT4(i+1).h4 MT33(i).h1 = b1 MT33(i).h2 = b2 MT33(i).h3 = b3 Next Call Det33(MT33(),dt) k = dt ND(1).h1 = k For i=1 to 4 If i = 1 then b1 = MT4(i).h2 b2 = MT4(i).h3 b3 = MT4(i).h4 MT33(i).h1 = b1 MT33(i).h2 = b2 MT33(i).h3 = b3 ElseIf i>2 then b1 = MT4(i).h2 b2 = MT4(i).h3 b3 = MT4(i).h4 MT33(i-1).h1 = b1 MT33(i-1).h2 = b2 MT33(i-1).h3 = b3 End If Next Call Det33(MT33(),dt) k = dt ND(2).h1 = -k For i=1 to 4 If i=1 or i=2 then b1 = MT4(i).h2 b2 = MT4(i).h3 b3 = MT4(i).h4 MT33(i).h1 = b1 MT33(i).h2 = b2 MT33(i).h3 = b3 ElseIf i>3 then b1 = MT4(i).h2 b2 = MT4(i).h3 b3 = MT4(i).h4 MT33(i-1).h1 = b1 MT33(i-1).h2 = b2 MT33(i-1).h3 = b3 End If Next Call Det33(MT33(),dt) k = dt ND(3).h1 = k For i=1 to 3 b1 = MT4(i).h2 b2 = MT4(i).h3 b3 = MT4(i).h4 MT33(i).h1 = b1 MT33(i).h2 = b2 MT33(i).h3 = b3 Next Call Det33(MT33(),dt) k = dt ND(4).h1= -k '222222222222222222222222222222222222222222 For i=1 to 3 b1 = MT4(i+1).h1 b2 = MT4(i+1).h3 b3 = MT4(i+1).h4 MT33(i).h1 = b1 MT33(i).h2 = b2 MT33(i).h3 = b3 Next Call Det33(MT33(),dt) k = dt ND(1).h2 = -k For i=1 to 4 If i = 1 then b1 = MT4(i).h1 b2 = MT4(i).h3 b3 = MT4(i).h4 MT33(i).h1 = b1 MT33(i).h2 = b2 MT33(i).h3 = b3 ElseIf i>2 then b1 = MT4(i).h1 b2 = MT4(i).h3 b3 = MT4(i).h4 MT33(i-1).h1 = b1 MT33(i-1).h2 = b2 MT33(i-1).h3 = b3 End If Next Call Det33(MT33(),dt) k = dt ND(2).h2 = k For i=1 to 4 If i=1 or i=2 then b1 = MT4(i).h1 b2 = MT4(i).h3 b3 = MT4(i).h4 MT33(i).h1 = b1 MT33(i).h2 = b2 MT33(i).h3 = b3 ElseIf i>3 then b1 = MT4(i).h1 b2 = MT4(i).h3 b3 = MT4(i).h4 MT33(i-1).h1 = b1 MT33(i-1).h2 = b2 MT33(i-1).h3 = b3 End If Next Call Det33(MT33(),dt) k = dt ND(3).h2 = -k For i=1 to 3 b1 = MT4(i).h1 b2 = MT4(i).h3 b3 = MT4(i).h4 MT33(i).h1 = b1 MT33(i).h2 = b2 MT33(i).h3 = b3 Next Call Det33(MT33(),dt) k = dt ND(4).h2= k '33333333333333333333333333333333333333 For i=1 to 3 b1 = MT4(i+1).h1 b2 = MT4(i+1).h2 b3 = MT4(i+1).h4 MT33(i).h1 = b1 MT33(i).h2 = b2 MT33(i).h3 = b3 Next Call Det33(MT33(),dt) k = dt ND(1).h3 = k For i=1 to 4 If i = 1 then b1 = MT4(i).h1 b2 = MT4(i).h2 b3 = MT4(i).h4 MT33(i).h1 = b1 MT33(i).h2 = b2 MT33(i).h3 = b3 ElseIf i>2 then b1 = MT4(i).h1 b2 = MT4(i).h2 b3 = MT4(i).h4 MT33(i-1).h1 = b1 MT33(i-1).h2 = b2 MT33(i-1).h3 = b3 End If Next Call Det33(MT33(),dt) k = dt ND(2).h3 = -k For i=1 to 4 If i=1 or i=2 then b1 = MT4(i).h1 b2 = MT4(i).h2 b3 = MT4(i).h4 MT33(i).h1 = b1 MT33(i).h2 = b2 MT33(i).h3 = b3 ElseIf i>3 then b1 = MT4(i).h1 b2 = MT4(i).h2 b3 = MT4(i).h4 MT33(i-1).h1 = b1 MT33(i-1).h2 = b2 MT33(i-1).h3 = b3 End If Next Call Det33(MT33(),dt) k = dt ND(3).h3 = k For i=1 to 3 b1 = MT4(i).h1 b2 = MT4(i).h2 b3 = MT4(i).h4 MT33(i).h1 = b1 MT33(i).h2 = b2 MT33(i).h3 = b3 Next Call Det33(MT33(),dt) k = dt ND(4).h3= -k '4444444444444444444444444444444444444444444 For i=1 to 3 b1 = MT4(i+1).h1 b2 = MT4(i+1).h2 b3 = MT4(i+1).h3 MT33(i).h1 = b1 MT33(i).h2 = b2 MT33(i).h3 = b3 Next Call Det33(MT33(),dt) k = dt ND(1).h4 = -k For i=1 to 4 If i = 1 then b1 = MT4(i).h1 b2 = MT4(i).h2 b3 = MT4(i).h3 MT33(i).h1 = b1 MT33(i).h2 = b2 MT33(i).h3 = b3 ElseIf i>2 then b1 = MT4(i).h1 b2 = MT4(i).h2 b3 = MT4(i).h3 MT33(i-1).h1 = b1 MT33(i-1).h2 = b2 MT33(i-1).h3 = b3 End If Next Call Det33(MT33(),dt) k = dt ND(2).h4 = k For i=1 to 4 If i=1 or i=2 then b1 = MT4(i).h1 b2 = MT4(i).h2 b3 = MT4(i).h3 MT33(i).h1 = b1 MT33(i).h2 = b2 MT33(i).h3 = b3 ElseIf i>3 then b1 = MT4(i).h1 b2 = MT4(i).h2 b3 = MT4(i).h3 MT33(i-1).h1 = b1 MT33(i-1).h2 = b2 MT33(i-1).h3 = b3 End If Next Call Det33(MT33(),dt) k = dt ND(3).h4 = -k For i=1 to 3 b1 = MT4(i).h1 b2 = MT4(i).h2 b3 = MT4(i).h3 MT33(i).h1 = b1 MT33(i).h2 = b2 MT33(i).h3 = b3 Next Call Det33(MT33(),dt) k = dt ND(4).h4= k Dim ND2() as Matranb4 Dim det as float Call Det44(MT4(),det) Call Chuyenvi(ND(),ND2()) For i= 1 to 4 MTND(i).h1 = ND2(i).h1/det MTND(i).h2 = ND2(i).h2/det MTND(i).h3 = ND2(i).h3/det MTND(i).h4 = ND2(i).h4/det Next End Sub 2.2.2.8 Thuật toán tính và hiển thị biến dạng xoay và sai số Đây là thuật toán đòi hỏi sự tổng hợp của các thuật toán trình bày ở trên trong việc tính toán để đưa ra được kết quả là giá trị của biến dạng xoay trong một tam giác địa động lực, từ đó tính cho tất cả các tam giác trong lưới tam giác đã được chọn.Đó là để tính được giá trị biến dạng xoay ω trong một tam giác thì đầu tiên phải sử dụng hàm lấy dữ liệu từ file số liệu đầu vào, số liệu này được lưu vào biến DI() như một ma trận lưu dữ liệu. Đây là số liệu thô ban đầu gồm tất cả các dữ liệu liên quan đến số liệu đầu vào cho việc tính toán biến dạng và sai số cho biến dạng xoay. Từ số liệu thô này ta phải chuyển dữ liệu để có thể tính được ma trận F như theo công thức (1.15) là F=A-1Vij, ma trận A và Vij là 2 ma trận cần tìm, ma trận A là ma trận khoảng cách, và Vij là ma trận giá trị vận tốc đều có thể tính được theo công thức tính biến dạng trình bày ở mục 1.1.5.2, ma trận A các giá trị khoảng cách sẽ được tính bằng hàm tính khoảng cách trong Mapbasic. Và việc tính giá trị của 2 ma trận này sẽ được thực hiện theo hàm có tên ChuyenDuLieuMatran(MTDL() as Diem,VIJ() as Matranb4,A() as Matranb4,MS() as Matranb4) Trong đó có biến đầu vào là MTDL() sẽ chứa dữ liệu đầu vào ban đầu, và các biến đầu ra là :VIJ() là biến lưu giá trị của ma trận Vij, A() là biến lưu giá trị của ma trận A, MS() là ma trận lưu giá trị sai sổ được tính cho sai số vận tốc dịch chuyển.Hàm này được định nghĩa trong code của chương trình như sau : Sub ChuyenDuLieuMatran(MTDL() as Diem,VIJ() as Matranb4,A() as Matranb4,MS() as Matranb4) Redim MS(1) Redim VIJ(1) Redim A(4) VIJ(1).h1=(MTDL(1).vx - MTDL(2).vx)/1000 VIJ(1).h2=(MTDL(1).vy - MTDL(2).vy)/1000 VIJ(1).h3=(MTDL(1).vx - MTDL(3).vx)/1000 VIJ(1).h4=(MTDL(1).vy - MTDL(3).vy)/1000 MS(1).h1=abs((MTDL(1).svx - MTDL(2).svx)/1000) MS(1).h2=abs((MTDL(1).svy - MTDL(2).svy)/1000) MS(1).h3=abs((MTDL(1).svx - MTDL(3).svx)/1000) MS(1).h4=abs((MTDL(1).svy - MTDL(3).svy)/1000) Dim kx12,ky12,kx13,ky13 as float kx12 = distance(MTDL(1).x,MTDL(1).y,MTDL(2).x,MTDL(1).y,"m") ky12 = distance(MTDL(1).x,MTDL(1).y,MTDL(1).x,MTDL(2).y,"m") kx13 = distance(MTDL(1).x,MTDL(1).y,MTDL(3).x,MTDL(1).y,"m") ky13 = distance(MTDL(1).x,MTDL(1).y,MTDL(1).x,MTDL(3).y,"m") If MTDL(2).x > MTDL(1).x then kx12 = -kx12 End If If MTDL(2).y > MTDL(1).y then ky12 = -ky12 End If If MTDL(3).x > MTDL(1).x then kx13 = -kx13 End If If MTDL(3).y > MTDL(1).y then ky13 = -ky13 End If A(1).h1 = kx12 A(2).h1 = ky12 A(3).h2 = kx12 A(4).h2 = ky12 A(1).h3 = kx13 A(2).h3 = ky13 A(3).h4 = kx13 A(4).h4 = ky13 End Sub Sau khi tìm được ma trận A và Vij phục vụ cho việc tìm ma trận F thì ta sẽ áp dụng thuật toán tính nghịch đảo của ma trận A và thuật toán nhân 2 ma trận để tính F theo công thức (1.15). Tìm được ma trận F từ đó ta tính được biến dạng xoay theo công thức (1.14) là ω=12∂vx∂y-∂vy∂x, Sử dụng hàm lấy dữ liệu chọn tam giác ta sẽ tính được các giá trị ω trong các tam giác trong lưới được chọn. Việc này được thực hiện tuần tự trong từng tam giác một kèm theo cách tính sai số của trong từng tam giác, Tính các giá trị biến dạng xoay rồi ra cũng sẽ thể hiện lưới tam giác ra bản đồ của Mapinfo bằng hàm có tên Biendangxoay(x as float, y as float, w as float) trong code của chương trình Hàm này có tác dụng nhận các biến đầu vào là kinh độ x, vỹ độ y, giá trị biến dạng xoay w để thể hiện biến dạng xoay ra bản đồ mapinfo giống trong hình 1.2((c),(d)). Hàm này sẽ xét các điều kiện giá trị biến dạng xoay w là lớn hay bé, ngược chiều hay cùng chiều kim đồng hồ mà có những góc xoay hay kiểu xoay khác nhau được thể hiện trên bản đồ Mapinfo theo tỉ lệ của bản đồ.Hàm này được định nghĩa trong code của chương trình như sau : Sub Biendangxoay(x as float, y as float, w as float) Dim tle,sw as float tle = 500 sw = w*(10^9) Dim y1,x2,y2 as float y1=y + 1 Dim o,ob,ob1,ob2,ob3,oo as object o = CreatePoint(x,y) ob1 = CreateLine(x,y,x,y1) ob2 = RotateAtPoint(ob1,-5,Centroid(o)) x2 = ObjectGeography(ob2, OBJ_GEO_LINEENDX) y2 = ObjectGeography(ob2, OBJ_GEO_LINEENDY) ob3 = CreateLine(x,y1,x2,y2) ob = combine(ob1,ob2) ob = combine(ob,ob3) oo= RotateAtPoint(ob,-5,Centroid(o)) If sw>0 then If sw <= 1 then Insert Into Lonlat(obj) values (ob) ElseIf sw >1 and sw <=3 then oo= combine(ob,oo) Insert Into Lonlat(obj) values (oo) ElseIf sw>3 and sw <=5 then oo = combine (ob,oo) ob1 = RotateAtPoint(ob,-10,Centroid(o)) oo = combine (ob,oo) Insert Into Lonlat(obj) values (oo) Else oo = combine (ob,oo) ob1 = RotateAtPoint(ob,-30,Centroid(o)) oo = combine (ob,oo) Insert Into Lonlat(obj) values (oo) End If Else oo= RotateAtPoint(ob,5,Centroid(o)) If sw >= -1 then Insert Into Lonlat(obj) values (oo) ElseIf sw >1 and sw <=3 then oo= combine(ob,oo) Insert Into Lonlat(obj) values (oo) ElseIf sw=-5 then oo = combine (ob,oo) ob1 = RotateAtPoint(ob,10,Centroid(o)) oo = combine (ob,oo) Insert Into Lonlat(obj) values (oo) Else oo = combine (ob,oo) ob1 = RotateAtPoint(ob,30,Centroid(o)) oo = combine (ob1,oo) Insert Into Lonlat(obj) values (oo) End If End if End Sub Tiếp theo sẽ là phần ghi dữ liệu đó là hàm ghi file dữ liệu Print File trong Mapbasic và việc này se được ghi tuần tự khi mỗi giá trị biến dạng xoay và sai số được tính trong một tam giác.File này sẽ được lưu vào đường dẫn chứa file chọn tam giác. Công việc tính toán giá trị biến dạng xoay của lưới tam giác trong vùng nghiên cứu, việc in file kết quả, xuất hiển thị ra bản đồ sẽ được xử lý bằng hàm có tên BienDang trong code của chương trình. Hàm này được định nghĩa như sau: Sub BienDang Dim x,y,vx,vy,svx,svy,ex,ey,exy,e1,e2,w,g as float Dim i as integer Call Laydulieutamgiac(TG()) Dim DDI() as Diem ReDim DDI(3) Dim d as integer Dim tle as float tle = tile/10 Dim se1,se2,sw,sg,sex,sey,sexy,xg,yg as float Dim V(),A(),NA(),MTN(),MS(),MTN1() as Matranb4 Dim st2,st3,st4,st6,st7,st8 as string*15 Dim st,st1,st5,st9 as string*8 Open File "ketqua.txt" for Output as #1 Print #1, "TenDagiac"+" "+" LON "+" " + " LAT "+" " + "BienDangXoay_W"+" "+ " SaiSo_W " Print #1, " Name "+" "+ " (Độ) "+" " + " (Độ) "+" " + " (Radian/Năm) "+" "+ " (Radian/Năm) " For i=1 to (ubound(TG)-1) d = TG(i).d1 DDI(1).tendiem = DI(d).tendiem DDI(1).x=DI(d).x DDI(1).y=DI(d).y DDI(1).vx=DI(d).vx DDI(1).vy=DI(d).vy DDI(1).svx=DI(d).svx DDI(1).svy=DI(d).svy d = TG(i).d2 DDI(2).tendiem = DI(d).tendiem DDI(2).x=DI(d).x DDI(2).y=DI(d).y DDI(2).vx=DI(d).vx DDI(2).vy=DI(d).vy DDI(2).svx=DI(d).svx DDI(2).svy=DI(d).svy d = TG(i).d3 DDI(3).tendiem = DI(d).tendiem DDI(3).x=DI(d).x DDI(3).y=DI(d).y DDI(3).vx=DI(d).vx DDI(3).vy=DI(d).vy DDI(3).svx=DI(d).svx DDI(3).svy=DI(d).svy Call Trongtam(DDI(),xg,yg) Create text TG(i).tendiem (xg, yg) (xg + 0.2,yg + 0.2) Font MakeFont ("Helvetica",1,12,Black,Green) Call ChuyenDuLieuMatran(DDI(), V(),A(),MS()) Call Ngichdao(A(),NA()) Call NhanMT(NA(),V(), MTN()) Call NhanMT(NA(),MS(),MTN1()) ex=MTN(1).h1 ey=MTN(1).h4 exy = (MTN(1).h2 + MTN(1).h3)/2 e1= (ex+ey)/2 + sqr(exy^2 + ((ex-ey)/2)^2) e2= (ex+ey)/2 - sqr(exy^2 + ((ex-ey)/2)^2) w= (MTN(1).h2 - MTN(1).h3)/2 g=atn(2*exy/(ey-ex))/2 g=g*(180/3.141) sex=abs(MTN1(1).h1) sey=abs(MTN1(1).h4) sexy = (abs(MTN1(1).h2) + abs(MTN1(1).h3))/2 se1= (sex+sey)/2 + sqr(sexy^2 + ((sex+sey)/2)^2) se2= abs((sex+sey)/2 - sqr(sexy^2 + (abs(sex-sey)/2)^2)) sw= abs(abs(MTN1(1).h2) - abs(MTN1(1).h3))/2 sg=atn(2*sexy/abs(sey-sex))/2 sg=sg*(180/3.141) st= Str$(xg) st1= Str$(yg) st2= Str$(e1) st3= Str$(e2) st4= Str$(w) st5= Str$(g) st6= Str$(se1) st7= Str$(se2) st8= Str$(sw) st9= Str$(sg) Print #1, TG(i).Tendiem + " "+ st + " " + st1 + " " +" " + st4 +" " + st8 g = -g ' Call Vemuiten(xg,yg,e1,e2,g,tle) Call Biendangxoay(xg,yg,w) Next Close File #1 Alter Menu Item BienDang Disable Alter Menu Item HienDiem Enable End Sub Ngoài ra còn một số hàm nhỏ khác phục vụ cho các hàm chính được trình bày ở trên. Các hàm này được chỉ rõ trong code của chương trình. 2.2.2 ĐÁNH GIÁ VỀ THUẬT TOÁN TRONG CHƯƠNG TRÌNH Cũng do Mapbasic không hỗ trợ mảng 2 chiều do đó thuật toán còn khá dài và hạn chế tính đa dạng của nó trong tính toán (như trong ngôn ngữ lập trình C có hỗ trợ mảng 2 chiều nên việc tính ma trận trong thuật toán rất ngắn gọn mà sử dụng cho việc tính cho các ma trận lớn nhiều bậc khác nhau rất linh hoạt mà ở đây chỉ hạn chế trong tính ma trận bậc 4). Từ đó mà những thuật toán trong chương trình tương đối khó hiểu và dài dòng. Tuy nhiên nó cũng đáp ứng được những yêu cầu đề ra để tính toán thể hiện các số liệu biến dạng xoay trong lưới tam giác địa động lực của khu vực nghiên cứu. Việc phân tách các hàm xử lý cho từng vấn đề giúp người lập trình dễ dàng trong việc sửa chữa và nâng cấp. CHƯƠNG 3 : GIAO DIỆN VÀ SỬ DỤNG CHƯƠNG TRÌNH Để minh họa cho cách sử dụng chương trình em sẽ minh họa chương trình sẽ tính biến dạng xoay trong một tam giác địa động lực theo số liệu trong luận văn tiến sĩ “BIẾN DẠNG KIẾN TẠO HIỆN ĐẠI KHU VỰC BIỂN ĐÔNG VIỆT NAM THEO SỐ LIỆU DỊCH CHUYỂN GPS” của nghiên cứu sinh Nguyễn Văn Hướng tại Phòng Địa Động Lực – Viện Địa Chất mà trong luận văn này sử dụng để tính toán biến dạng. 3.1 GIAO DIỆN CHƯƠNG TRÌNH Chương trình có giao diện trong khung đỏ ở hình 3.1 Hình 3.1: Giao diện của chương trình Chương trình khi chạy lên Mapinfo nằm ở ngoài cùng của menu chính, có tên là Tinh Bien Dang (Hình 3.1), khi click chuột vào tool này nó sẽ sổ xuống là các lựa chọn : Tinh thực hiện việc chạy chương trình cho tính biến dạng, Ket Thuc C.trinh thực hiện thoát tool khỏi Mapinfo và Thoat thực hiện việc thoát tất cả cửa sổ Mapinfo. Khi lựa chọn menu Tinh (Hình 3.1) thì ở đây có 2 lựa chọn nhỏ đó là Hien Thi Diem có nhiệm vụ hiển thị các điểm hay tốc độ chuyển dịch ra cửa sổ bản đồ Mapinfo và Tinh Bien Dang thực hiện việc tính biến dạng xoay, in file dữ liệu kết quả, hiển thị lưới tam giác và biến dạng xoay trong tam giác ra cửa sổ bản đồ của Mapinfo. 3.2 SỬ DỤNG CHƯƠNG TRÌNH Để minh họa cho việc chạy chương trình ta kiểm định bằng việc chạy thử một ví dụ. Đó là tính biến dạng xoay cho khu vực biển Đông theo ssố liệu GPS trong đề tài “BIẾN DẠNG KIẾN TẠO HIỆN ĐẠI KHU VỰC BIỂN ĐÔNG VIỆT NAM THEO SỐ LIỆU DỊCH CHUYỂN GPS” của nghiên cứu sinh Nguyễn Văn Hướng tại Phòng Địa Động Lực – Viện Địa Chất. Bảng dữ liệu đầu vào gồm các số liệu như sau : STT Tên điểm Kinh độ Vỹ độ Ve Vn SVe SVn 1 STT1 114,331 11,429 23,40 -10,57 1,44 1,34 2 BLV1 107,723 20,128 30,12 -12,46 1,36 1,29 3 CDA1 106,652 8,692 21,76 -9,84 1,39 1,32 4 DOHO 106,616 17,507 26,63 -9,48 1,37 1,31 5 HOCM 106,560 10,849 22,00 -13,75 1,47 1,36 6 NTUS 103,680 1,346 27,41 -7,05 0,87 0,85 7 KUNM 102,797 25,030 30,63 -17,77 0,80 0,79 8 PIMO 121,078 14,636 -29,00 4,70 0,90 0,90 9 TNML 120,987 24,798 28,66 -9,62 0,80 0,80 10 TNSM 116,725 20,702 29,70 -12,90 0,20 0,10 11 BRG1 120,601 18,520 -48,50 9,80 0,40 0,40 12 BTS3 121,963 20,438 -38,50 27,90 0,90 0,60 13 S102 121,558 22,037 -38,80 39,30 0,30 0,30 14 S23R 120,606 22,645 -22,90 -10,60 0,10 0,10 15 LNI1 106,539 11,822 25,83 -13,47 2,50 2,19 16 JB21 110,306 25,186 33,01 -13,49 0,07 0,15 17 XIAM 118,082 24,449 32,10 -14,40 0,10 0,10 18 YONG 112,335 16,834 30,10 -11,60 0,10 0,10 19 TABA 108,891 0,863 29,83 -9,20 0,95 1,26 20 BRUN 115,031 4,966 28,66 -10,65 2,10 3,30 21 PUER 118,851 10,086 32,00 -13,60 1,00 0,70 Bảng 3.1: Số liệu tốc độ chuyển dịch theo số liệu GPS khu vực biển Đông Trong đó Ve,Vn lần lượt là giá trị tốc độ chuyển dịch theo hướng Đông và Bắc, đơn vị là mm/năm. SVe và SVn lần lượt các sai số của Ve và Vn,đơn vị là mm/năm. Ngoài bảng số liệu ta còn có dữ liệu đầu vào nữa đó là một ảnh mô hình số độ cao DEM (Hình 3.2) cho khu vực Biển Đông để làm nền cho hiển thị cửa sổ bản đồ trong Mapinfo. Hình 3.2: Mô hình số độ cao DEM khu vực Biển Đông Trước khi sử dụng chương trình đầu tiên ta tạo một file số liệu Excel (Hình 2.3) theo số liệu ở bảng 3.1 (Lưu ý: file Excel lưu định dạng của Excel 97 -2003 có đuôi.XLS) Sau khi tạo file dữ liệu đầu vào xong ta khởi động Mapinfo. Chọn Tool\Run MapBasic Program (Hình 3.3). Sau đó sẽ xuất hiện hộp thoại chọn file chạy chương trình (Hình 3.4). Ta chọn đường dẫn tới chương trình được lưu có tên biendang.MBX rồi chọn Open.Ta được giao diện như trong hình 3.5. Hình 3.3: Giao diện chạy chương trình MapBasic Hình 3.4: Hộp thoại chọn chương trình chạy Hình 3.5: Giao diện chương trình khi khởi động Trước khi chạy chương trình ta tạo một bản đồ mới có hệ tọa tọa địa lý (Lon/Lat) được minh họa bằng các hình dưới đây: Hình 3.6: Giao diện lựa chọn tạo một bản đồ mới Ta chọn File\New Table... để tạo một bản đồ mới.Sẽ xuất hiện hộp thoại lựa chọn tạo một Table mới (Hình 3.7): Hình 3.7: Hộp thoại tạo một Table mới Hộp thoại xuất hiện ta click chọn Create để tạo một Table mới.Sẽ xuất hiện hộp thoại lựa chọn (Hình 3.8). Ta chọn Add Field để tạo một trường bản đồ, Click chọn Projection để chọn hệ qui chiếu và hệ tọa độ chỏ Table. Hình 3.8: Hộp thoại lựa chọn tạo Table mới Hình 3.9: Hộp thoại lựa chọn hệ tọa độ cho Table Ở đây ở mục Category ta chọn hệ qui chiếu là hệ tọa độ địa lý Longitude /Latitude Và mục Category Member là lựa chọn hệ tọa độ ta chọn Longitude /Latitude.Sau đó ta click Ok sẽ xuất hiện hộp thoại chọn đường dẫn và lưu tên Table (Hình 3.10). Hình 3.10: Chọn đường dẫn và lưu tên Table mới Ta click chọn Save lúc này ở trong thư mục được chọn sẽ được tạo các file Lonlat.DAT,Lonlat.ID,Lonlat.Map, và Lonlat.TAB là những tệp liên quan đến table, dữ liệu, bản đồ, liên kết của Table có tên Lonlat.Lúc này cửa sổ bản đồ trong Mapinfo sẽ xuất hiện như hình 3.11. Ta khởi chạy chương trình Tinh Bien Dang bằng chọn Tinh Bien Dang\Tinh\Hien Thi Diem. Lúc này chương trình sẽ xuất hiện một hộp thoại chọn đường dẫn tới file số liệu đầu vào tính biến dạng, file này được chọn có định dạng.xls. Khi chọn được file dữ liệu đầu vào cho tính toán biến dạng ta chọn Open (Hình 3.12), Khi đó sẽ xuất hiện một hộp thoại hỏi “ Ban co muon the hien toc do dich chuyen” Chọn “ Co” để hiện thị tốc độ chuyển dịch, Chọn “ Khong” để chỉ hiển thị lưới điểm (Hình 3.13). Hình 3.11: Giao diện sau khi tạo Table mới tên Lonlat Hình 3.12: Hộp thoại chọn file số liệu đầu vào Hình 3.13: Hộp thoại lựa chọn có thể hiện tốc độ dịch chuyển hay không Nếu ta chọn “Khong” thì chương trình chạy sẽ hiển thị toàn bộ các điểm đo ra cửa sổ bản đồ của Mapinfo (Hình 3.14) gồm ký hiệu điểm, tên điểm đo và vị trí của nó trên cửa sổ bản đồ, khi phóng to hình ta nhìn rõ hơn (Hình 3.15). Hình 3.14: Lưới điểm đo khi chạy chương trình Khi phóng to một góc trên cửa sổ bản đồ này ta thấy rõ hơn tên điểm và ký hiệu của điểm đo trong cửa sổ bản đồ này (Hình 3.15). Hình 3.15: Phóng to một góc hình 3.14 Trong hình 3.15 ta quan sát được 4 điểm đo đó là LNI1, HOCM, CDA1,STT1 là 4 điểm đo có trong dữ liệu file đầu vào, được ký hiệu điểm bởi một chấm đen và chữ có màu xanh,nền trắng.Trên hộp thoại hỏi (Hình 3.13) nếu ta chọn “Co” thì trên cửa sổ bàn đồ sẽ xuất hiện các đường dịch chuyển trong lưới đo (Hình 3.16). Hình 3.16: Lưới tốc độ chuyển dịch Khi phóng to một góc của bản đồ thể hiện lưới dịch chuyển ta sẽ quan sát được rõ hơn hình dáng mũi tên thể hiện tốc độ chuyển dịch (Hình 3.17). Hình 3.17: Phóng to một góc hình 3.16 Trên hình 3.17 ta quan sát thấy tại vị trí các điểm đo có một đường mũi tên màu đen và hình elip máu xanh lục xuất hiện ở đầu mũi tên. Ở đây ta có thể thấy hình mũi tên chính là đường thể hiện tốc độ chuyển dịch tại vị trí điểm đo được thể hiện trên bản đồ. Hướng của mũi tên chính là hướng dịch chuyển của vị trí điểm đo, còn hình elip màu xanh lục thể hiện sai số tốc độ dịch chuyển tại vị trí điểm đấy. Elip và đường mũi tên đều được phóng đại theo cùng một tỷ lệ nhất định (vì thông số dịch chuyển rất bé nên khi hiển thị trên bản đồ cần phóng đại để có thể quan sát được). Để minh họa lưới dịch chuyển được quan sát thực tế ta cho thể hiện thêm bản đồ số độ cao DEM của khu vực Biển Đông lên bản đồ sẽ cho ta nhìn một cách trực quan hơn (Hình 3.18). Cách thức thêm bản đồ mô hình số độ cao vào làm nền minh họa sẽ được trình bày ở phía dưới. Hình 3.18: Lưới tốc độ dịch chuyển trên nền bản đồ độ cao khu vực biển Đông Sau khi chạy chương trình trong lựa chọn menu Hien Thi Diem, lúc này menu Tinh Bien Dang sẽ được kích hoạt đậm màu và menu Hien Thi Diem bị vô hiệu hóa và mờ đi (Hình 3.19).Để tiếp tục chạy chức năng Tinh Bien Dang ta tạo file dữ liệu chọn tam giác (Quan sát lưới điểm đo hình 3.14) để có tiếp tục (Hình 2.4).Sau khi tạo file dữ liệu chọn tam giác ta chọn menu Tinh Bien Dang để tiếp tục chạy chương trình (Hình 3.18). Hình 3.19: Menu của chương trình Khi chạy chọn Tinh Bien Dang trong Menu Tinh (Hình 2.19) Chương trình sẽ xuất hiện một hộp thoại chọn đường dẫn tới file dữ liệu chọn tam giác vừa mới được tạo (Hình 3.20). Hình 3.20: Hộp thoại chọn mở file dữ liệu chọn tam giác Ta sẽ chọn file có tên “File chon tam giac.xls” Rồi click chọn Open (Hình 3.20). Sau khi khởi chạy chương trình cửa sổ bản đồ lúc này hiển thị kết quả là một lưới tam giác có các cạnh màu đỏ được liên kết với nhau (Lưới tam giác này được chọn theo người dùng đã lưu vào file số liệu chọn tam giác(Hình 2.4)) (Hình 3.21). Tại mỗi tam giác tại vị trí trọng tâm là ký hiệu của biến dạng xoay, ký hiệu này tùy thuộc vào giá trị của biến dạng xoay mà được thể hiện khác nhau. Khi phóng to một góc bản đồ ta sẽ quan sát rõ hơn cách thể hiện biến dạng xoay trong một tam giác (Hình 3.23). Hình 3.21: Lưới tam giác thể hiện biến dạng xoay Trong đó các ký hiệu của biến dạng xoay trong hình 3.21 có ý nghĩa sau (Hình 3.22): -1 -3 -5 <-5 +1 +3 +5 >+5 Ngược chiều kim đồng hồ Đơn vị (nano radian/năm) Cùng chiều kim đồng hồ Đơn vị (nano radian/năm) Hình 3.22: Ký hiệu biến dạng xoay Hình 3.23: Phóng to một góc bản đồ của hình 3.21 Hình 3.21 thể hiện kết quả hiển thị trên cửa sổ bản đồ của Mapinfo bằng kết quả của hiển thị lưới tam giác và thể hiện tốc độ xoay.Cùng đó chương trình cũng xuất một file dữ liệu kết quả của tính toán tốc độ biến dạng xoay và sai số của nó, file này có tên mặc định là “ketqua.txt” được lưu vào đường dẫn có chứa file dữ liệu đầu vào. File dữ liệu kết quả này cho ta biết các thông số về tên tam giác được chọn, kinh độ,vĩ độ trọng tâm của tam giác đó, giá trị biến dạng xoay và sai số của nó trong tam giác được tính. File dữ liệu kết quả được minh họa ở hình 3.24 : Hình 3.24: File kết quả tính biến dạng xoay và sai số Ở hình 3.21 thể hiện lưới tam giác cùng ký hiệu thể hiện tốc độ biến dạng xoay của khu vực tính toán, vấn đề đặt ra là chèn một nền bản đồ thực tế để minh họa cho khu vực nghiên cứu.Trong phần đầu chương 3 đã có mô tả dữ liệu đầu vào gồm có một ảnh số mô hình độ cao (DEM) của khu vực Biển Đông và khu vực lân cận. Nhiệm vụ bây giờ là đưa file ảnh lên cửa sổ Mapinfo làm nền cho kết quả sau việc chạy chương trình.Dưới đây là trình tự các bước thực hiện: Sau khi chạy chương trình đã có kết quả được hiển thị trên cửa sổ bản đồ trong Mapinfo, tiếp theo ta sẽ mở một Table có chứa ảnh số mô hình độ cao đã được Register (đăng ký). Đầu tiên từ thanh menu chính của Mapinfo ta chọn File\Open (Hình 3.25) một hộp thoại sẽ xuất hiện chọn đường dẫn tới file Table có chứa dữ liệu bản đồ mô hình số độ cao, ta chọn Open (Hình 3.26) lúc này cửa sổ bản đồ hiện thị bản đồ mô hình số độ cao cùng lưới tam giác trong đó (Hình 3.27). Hình 3.25: Menu mở một Table trong Mapinfo Hình 3.26: Hộp thoại chọn đường dẫn tới Table cần mở Hình 3.27: Lưới tam giác trên bản đồ mô hình số độ cao Khi đã có một kết quả cuối cùng trên hình 3.27 ta sẽ in bản đồ với những tỷ lệ khác nhau và tùy trên những bản giấy có khổ in khác nhau tùy thuộc vào lựa chọn của người dùng, những công cụ này có hỗ trợ sẵn trên Mapinfo. CHƯƠNG 4: ĐÁNH GIÁ VÀ PHÂN TÍCH KẾT QUẢ CHƯƠNG TRÌNH ĐẠT ĐƯỢC 4.1 ƯU ĐIỂM VÀ HẠN CHẾ CỦA CHƯƠNG TRÌNH 4.1.1 ƯU ĐIỂM Chương trình chạy nhanh và ổn định trên nền của phần mềm Mapinfo Giao diện dễ nhìn, phần mềm tương đối dễ sử dụng đối với người dùng Kết quả tính toán biến dạng xoay có độ chính xác cao (phù hợp với mục tiêu tính toán đề ra) Các số liệu đầu ra mang tính thực tế cao (vừa xuất luôn ra cửa sổ dữ liệu kiểu bản đồ và file dữ liệu kết quả) 4.1.2 HẠN CHẾ Chương trình vẫn chạy dựa trên là một tool của Mapinfo, không chạy trực tiếp trên nền Window (do vậy cần phải cài đặt Mapinfo trước khi sử dụng) Vẫn chưa bắt được hết các lỗi về dữ liệu cũng như lỗi hiển thị Chương trình vẫn còn dạng sơ khai, chưa có những công cụ để lựa chọn thuộc tính trong việc nhập cũng như xuất dữ liệu Kết quả tính sai số vẫn còn chưa chuẩn xác. Trên là những hạn chế trong chương trình, tuy nhiên những hạn chế này là do hạn hẹp về thời gian cũng như kiến thức, nếu có thời gian tiếp tục nghiên cứu đề tài này em tin tưởng rằng trong tương lai gần sẽ khắc phục được những hạn chế nói trên. 4.2 PHÂN TÍCH ĐÁNH GIÁ KẾT QUẢ ĐẠT ĐƯỢC Để đánh giá kết quả đạt được của chương trình, ta so sánh kết quả của chương trình tính được với kết quả tính của bộ phần mềm QOCA, kết quả so sánh này được thể hiện trong bảng 4.1 : STT Tên tam giác Kinh độ trọng tâm (degere) Vỹ độ trọng tâm (degere) Biến dạng xoay tính theo chương trình (radian/năm) Biến dạng xoay tính theo QOCA (radian/năm) 1 TG_01 105.712 20.8883 4.53437e-009 4.53e-09 2 TG_02 106.942 23.448 -9.18381e-010 -1.01e-09 3 TG_03 110.121 20.716 1.7459e-009 1.77 e-09 4 TG_04 113.122 20.9073 1.60446e-009 1.61 e-09 5 TG_05 115.038 23.4457 3.71338e-009 3.71 e-09 6 TG_06 108.891 18.1563 6.90789e-009 6.88 e-09 7 TG_07 111.094 15.2567 8.76603e-009 8.97 e-09 8 TG_08 114.464 16.3217 4.53523e-009 5.55 e-09 9 TG_09 109.162 13.586 -1.21725e-009 -1.18 e-09 10 TG_10 109.174 10.6477 4.36526e-009 1.23 e-09 11 TG_11 108.221 7.15567 -4.71115e-009 -4.81 e-09 12 TG_12 108.967 4.546 -1.77868e-009 -1.81 e-09 13 TG_13 112.751 5.75267 -2.45351e-009 -2.49 e-09 14 TG_14 116.071 8.827 3.363e-010 2.83 e-10 15 TG_15 118.087 12.0503 -5.83877e-008 -5.85 e-08 16 TG_16 117.378 15.589 1.29965e-009 1.53 e-09 17 TG_17 119.468 17.9527 -6.91522e-008 -7.18 e-08 18 TG_18 119.628 20.4197 -1.19375e-008 - 19 TG_19 119.757 22.5123 2.00909e-008 - 20 TG_20 118.598 23.3163 -3.48413e-009 - Bảng 4.1: Bảng so sánh kết quả biến dạng xoay với phần mềm QOCA Nhìn vào bảng so sánh kết quả tính biến dạng xoay trong chương trình ta thấy các tam giác có số thứ tự từ 1-17 là có sự lựa chọn về tam giác địa động lực là giống nhau nên ta chỉ so sánh kết quả tính được trong những tam giác này. Kết quả tính theo phần mềm QOCA được thể hiện trong chuyên đề tiến sĩ : “BIẾN DẠNG KIẾN TẠO HIỆN ĐẠI KHU VỰC BIỂN ĐÔNG VIỆT NAM THEO SỐ LIỆU DỊCH CHUYỂN GPS” của nghiên cứu sinh Nguyễn Văn Hướng.Theo PGS.TS Phan Trọng Trịnh-Viện Địa Chất thì các kết quả này nằm trong sự chênh lệch e-10 nên có thể chấp nhận đối với chương trình. Và kết quả mang độ chính xác tương đối cao. Trừ điểm tam giác TG_10 thì độ chênh lệch là tương đối lớn (3.12e-09) tuy vậy kết quả này vẫn có thể chập nhận được bởi vẫn nằm trong phạm vi sai số. Về đánh giá sai số, do kết quả tính toán của chương trình còn khác biệt khá lớn với kết quả tính trên phần mềm QOCA do vậy em không đưa vào so sánh, tuy vậy phần tính toán sai số là phần tương đối phức tạp nên hiện tại em chưa thể hoàn thành phần tính toán sai số này. Đánh giá chung thì chương trình đã đạt được mục tiêu đề ra là tính biến dạng xoay trong một tam giác địa động lực, và kết quả tính toán của nó có thể đưa vào sử dụng trong thực tế. KẾT LUẬN VÀ KIẾN NGHỊ Chương trình theo nhận định chung đã đáp ứng được yêu cầu tính biến dạng xoay và sai số của nó trong một tam giác địa động lực. Chương trình tương đối dễ sử dụng và dễ hiểu đối với người dùng. Kết quả của nó mang độ chính xác và tin cậy cao. Tuy nhiên hạn chế về thời gian cũng như kiến thức còn hạn hẹp nên kết quả này vẫn chưa được hiệu chỉnh tuy nhiên nó cũng đáp ứng được mục đích của tính biến dạng xoay trong một tam giác địa động lực. Chương trình vẫn chưa đa dạng trong phần lựa chọn các công cụ (như việc chọn màu chữ, cỡ chữ, đường, font... vẫn để chế độ mặc định khách quan của người lập trình viên). Tuy vậy chương trình vẫn dừng lại ở tính biến dạng xoay cho một tam giác, chưa tính được đối với các đa giác lớn hơn (tứ giác, ngũ giác...) và việc chọn tam giác phải mang tính chất thủ công, chưa tự động chọn được lưới tam giác.Và tính biến dạng còn thiếu phần nội suy, đó là tại bất kỳ đâu trong lưới tính biến dạng có thể xác định được giá trị của biến dạng xoay chứ không phải chỉ ở mỗi trọng tâm tam giác tính được. Tất cả những yếu tố trên đều là những hạn chế của chương trình, tuy nhiên để khắc phục được nó cần những kiến thức chuyên sâu về ngôn ngữ lập trình Mapbasic, sử dụng phần mềm Mapinfo, kiến thức về trắc địa, địa thống kê,cơ học, toán học...và cần một thời gian tương đối dài để thục hiện. Do sự hiểu biết có hạn của mình và thời gian hạn hẹp nên em chưa thể khắc phục được những điểm hạn chế này. Nếu còn thời gian tiếp tục nghiên cứu về đề tài này em có những hướng phát triển chương trình như sau: Xây dựng chương trình trên một nền tảng riêng biệt có thể chạy trực tiếp trên Window mà không phụ thuộc vào các phần mềm khác (Mapinfo). Điều này đòi hỏi một kiến thức trắc địa thâm sâu vì nếu tính trên nền tảng riêng như vậy đòi hỏi phải tự xây dựng một hàm tính khoảng cách giữa 2 điểm bất kỳ trên trái đất và trên những hệ tọa độ chiếu khác nhau (Điều này thực hiện trên Mapinfo rất dễ dàng, đây chính là ưu điểm của nó). Chương trình xây dựng sẽ nghiên cứu tính nhiều giá trị biến dạng kiến tạo trong vỏ Trái Đất chứ không phải biến dạng xoay, như biến dạng chính, biến dạng trượt,... và nghiên cứu tính biến dạng trong một đa giác bất kỳ chứ không phải trong một tam giác. Nghiên cứu thuật toán Delaunay để có thể lựa chọn các đa giác một cách tự động mang độ chính xác cao trong lựa chọn tính toán biến dạng. Nghiên cứu phương pháp nội suy Kringing để có thể tính biến dạng tại bất kỳ vị trí nào trong lưới tính biến dạng. Và nghiên cứu các phương pháp nội suy khác để áp dụng cho từng khu vực có điều kiện khác nhau một cách hợp lý. Một vấn đề nữa cần nghiên cứu đó là hiệu chỉnh biến dạng theo tỷ lệ của tam giác địa động lực,.... Để hoàn thành được đề tài tốt nghiệp “ Nghiên cứu tính biến dạng xoay trong một tam giác địa động lực” này em xin chân thành cám ơn PGS Phan Trọng Trịnh, trưởng phòng Địa Động Lực – Viện Địa Chất, người đã hướng dẫn chỉ bảo tận tình giúp em hoàn thành đề tài tốt nghiệp của mình một cách tốt nhất. Cùng đó em cũng xin cám ơn các thầy các cô trong bộ môn Tin Học Địa Chất đã tạo điều kiện cho em hoàn thành tốt đồ án tốt nghiệp của mình, một lần nữa em xin chân thành cám ơn!