KẾT QUẢ NGHIÊN CỨU BAN ĐẦU VỀ TỐC ĐỘ CHUYỂN DỊCH KIẾN TẠO HIỆN ĐẠI TRÊN BIỂN ĐÔNG

PHAN TRỌNG TRỊNH1, NGÔ VĂN LIÊM1, NGUYỄN TUẤN ANH2, VY QUỐC HẢI1, TRẦN ĐÌNH TÔ1, NGUYỄN VĂN HƯỚNG1, HOÀNG QUANG VINH1, BÙI VĂN THƠM1, NGUYỄN ĐĂNG TÚC1, NGUYỄN QUANG XUYÊN1, VŨ TUẤN HÙNG1, NGUYỄN HUY THỊNH1, TRẦN QUỐC HÙNG1, LÊ MINH TÙNG1, BÙI THỊ THẢO1, NGUYỄN VIỆT TIẾN1, ĐINH VĂN THẾ1               

1Viện Địa chất, Viện Khoa học và Công nghệ Việt Nam;
2Trung tâm Viễn thám Quốc gia, Bộ Tài nguyên và Môi trường

Tóm tắt: Bài báo trình bày kết quả nghiên cứu mới về chuyển dịch kiến tạo hiện đại trên biển Đông nhờ phân tích kết quả đo lặp GPS trong khoảng thời gian 2007-2008. Phương pháp đã được kiểm nghiệm bằng việc xử lý số đo GPS tại một số trạm GPS tại châu Âu bằng phiên bản 5.0 của phần mềm Bernese. Kết quả tính của chúng tôi trùng khít với kết quả tính toán của Viện Thiên văn thuộc trường Đại học Bern, Thuỵ Sĩ, với cùng số liệu đo. Các trạm đo GPS Láng, Bạch Long Vĩ, Song Tử Tây, Côn Đảo đã được liên kết với các trạm đo GPS trong hệ thống đo GPS quốc tế là BAKO, PIMO và KUNM. Sau khi tính chuyển dịch tương đối với giả định điểm Bạch Long Vĩ cố định, chúng tôi đã xác định chuyển dịch và tốc độ chuyển dịch tuyệt đối của các điểm trong hệ tọa độ toàn cầu IGS05. Kết quả nhận được cho thấy trạm Láng chuyển dịch về phía đông với tốc độ 43,4 mm/năm, chuyển dịch về phía nam với tốc độ 8,3 mm/năm; trạm Bạch Long Vĩ - về phía đông là 32,2 mm/năm và về phía nam là 8,3 mm/năm; trạm Song Tử Tây - về phía đông là 17,4 mm/năm và về phía nam là 14,2 mm/năm; trạm Côn Đảo - về phía đông là 16,1 mm/năm và về phía nam là 12,1 mm/năm. Sai số theo mỗi hướng chuyển dịch dao động từ 0,5 đến 0,6 mm/năm.


I. MỞ ĐẦU

Chuyển động kiến tạo hiện đại là nguồn gốc sâu xa gây ra các  thảm hoạ như động đất, sóng thần, ... có ảnh hưởng to lớn tới sinh mạng con người cũng như sự phát triển kinh tế của một vùng hay lãnh thổ. Về mặt khoa học, cơ chế của biến dạng thạch quyển luôn là chủ đề của các cuộc thảo luận quốc tế. Để xác định được tốc độ chuyển dịch kiến tạo hiện đại, các phương pháp trắc địa truyền thống từng được sử dụng, như phương pháp đo thuỷ chuẩn và phương pháp tam giác đạc. Trong quy mô nhỏ, các phương pháp trên có độ chính xác cao, nhưng tỏ ra hạn chế trên một quy mô rộng lớn. Để liên kết trên diện rộng các phương pháp trắc địa không gian như DOPPLER, VLBI (Very Long Baseline Interference), phương pháp định vị toàn cầu GPS đã được áp dụng. Công nghệ GPS, với việc sử dụng các máy thu 2 tần số và các phần mềm xử lý chuyên dụng có kết hợp với các dữ liệu bổ trợ được mô hình hoá, như tầng điện ly, tầng đối lưu, mô hình khí quyển, mô hình thuỷ triều, đã giúp các tính toán về vị trí và vận tốc chuyển dịch của vỏ Trái đất đạt tới sai số cỡ milimet trong phạm vi rộng tới hàng ngàn km. Để tính toán với độ chính xác cao các chuyển dịch kiến tạo hiện đại bằng công nghệ GPS, các trạm đo GPS quốc tế (International GPS Station - IGS) đã xây dựng, thu và xử lý số liệu liên tục của mạng lưới các trạm IGS trên phạm vi toàn cầu và đã tính được hướng cũng như vận tốc chuyển dịch của các mảng kiến tạo của Trái đất ở nhiều quy mô và mức độ khác nhau với độ chính xác ngày càng cao. Nhờ những số liệu đo liên tục tại các trạm IGS, cũng như các tham số bổ trợ về tốc độ chuyển dịch của các trạm IGS, ta có thể xác định được tốc độ chuyển dịch tuyệt đối của các trạm GPS. Nghiên cứu về khu vực Nam và Đông Nam Châu Á, đề án GEODYSSEA (Geodynamics of South and Southeast Asia), thông qua 3 chu kỳ đo 1994, 1996 và 1998, đã xác định được  tốc độ và hướng chuyển dịch tuyệt đối của vỏ Trái đất trong khu vực này với sai số 4-7 mm theo mặt bằng và 10 mm theo theo chiều cao [3, 7]. Trung Quốc đã thiết lập mạng lưới quan trắc chuyển dịch của vỏ Trái đất từ 1997 với 27 trạm đo liên tục và 1100 điểm đo không liên tục. Các đợt đo 1999, 2001 và 2004 cho thấy biến dạng trên cao nguyên Tây Tạng, các rìa của nó, đới Himalaya và Alty Tagh đã hấp thụ 90% chuyển dịch tương đối giữa mảng Ấn-Úc và mảng Châu Á [9, 15, 16]. Ở rìa đông của cao nguyên Tây Tạng, chuyển dịch hướng về phía đông và suy giảm nhanh ở rìa tây của Tứ Xuyên, trong khi Bắc Vân Nam chuyển dịch về đông nam và nam Vân Nam chuyển dịch chuyển thành nam - đông nam. Tại mảng Nam Trung Hoa, tốc độ chuyển dịch về phía đông trong khoảng 6-10 mm/năm [15, 16]. Trận động đất xảy ra ở Tứ Xuyên ngày 12/5/2008 với chấn cấp 7,9 là kết quả hấp thụ của chuyển dịch về phía đông qua ranh giới đứt gãy chờm nghịch ở rìa tây Tứ Xuyên.

 Hiện nay, các kết quả tính toán số liệu GPS đo liên tục 24/24 giờ đã đạt tới độ chính xác 0,5-1,2 mm theo phương ngang và 1-3 mm theo phương thẳng đứng [11]. Trong điều kiện thu tín hiệu tốt, mạng lưới trạm GPS hợp lý, máy thu được đặt cố định, chúng ta còn có thể xử lý với kết quả có sai số chỉ cỡ 0,3-0,4 mm theo phương ngang và 1,2 mm theo phương thẳng đứng [1].

II. KIỂM NGHIỆM PHƯƠNG PHÁP VÀ QUY TRÌNH TÍNH TOÁN

Để kiểm nghiệm khả năng làm chủ công nghệ tính toán GPS với độ chính xác cao, chúng tôi đã sử dụng số liệu đo tại một số trạm GPS khu vực châu Âu để tính trên phần mềm Bernese 5.0. Đây là phiên bản mới nhất có khả năng xử lý các trị số đo GPS với độ chính xác cao được xây dựng bởi Viện Thiên văn thuộc Trường Đại học Tổng hợp Bern, Thụy Sĩ. 

Dữ liệu được xử lý, tính toán trong ví dụ đối sánh này bao gồm số liệu thu của mạng lưới 8 trạm IGS khu vực Châu Âu (Hình 1), trong đó 3 trạm MATE, ONSA và VILL được sử dụng làm trạm tham chiếu trong hệ tọa độ quốc tế (International Terrestrial Reference Frame ITRF2000/IGS00). Trong 5 trạm còn lại có 2 trạm FFMJ và ZIMJ được trang bị máy thu có tích hợp cả tư liệu của các vệ tinh GPS và GLONASS, 2 trạm  BRUS và PTBB được sử dụng loại máy thu ASHTECH Z-XII3T. Khoảng cách giữa các trạm thu liền kề nằm trong khoảng 300-1200 km ngoại trừ 2 trạm thu (ZIMM và ZIMJ) ở Zimmerwald, Thụy Sĩ có khoảng cách là 14 m. Để phù hợp với tính toán của Viện Thiên văn trường Đại học Tổng hợp Bern, chúng tôi xử lý số liệu của 4 ngày đo, trong đó có 2 ngày của năm 2002 (ngày GPS 143 và 144) và 2 ngày của năm 2003 (ngày GPS 138 và 139). Để xử lý được dữ liệu trong tính toán này, chúng tôi đã sử dụng các sản phẩm cuối cùng của Trung tâm IGS cung cấp trên Internet như: lịch vệ tinh chính xác, mô hình tầng điện ly, file hiệu chỉnh giữa P1-C1, P1-P2 đối với vệ tinh và máy thu, tọa độ và vận tốc chuyển dịch của các trạm IGS trong hệ tọa độ toàn cầu IGS00, ... Cần lưu ý là trong khoảng thời gian từ năm 2000 đến nay, hệ tọa độ toàn cầu thay đổi từ hệ tọa độ toàn cầu IGS97 sang hệ tọa độ toàn cầu IGS00 vào ngày 2/12/2001, chuyển từ hệ tọa độ toàn cầu IGS00 sang hệ tọa độ toàn cầu IGB00 từ ngày 11/1/2004 và chuyển IGB00 sang IGS05 từ ngày 21/12/2006. Từ trước ngày 4/11/2000, thời gian chỉnh đồng hồ vệ tinh trong khoảng 15 phút, còn từ sau ngày đó đến nay thời gian chỉnh đồng hồ vệ tinh trong khoảng 5 phút [6, 12].


Hình 1. Vị trí các trạm GPS khu vực Châu Âu được sử dụng trong tính toán.

Sau lần bình sai thứ nhất, chúng tôi đã ước tính được giá trị tọa độ của các trạm đo và các tham số tầng đối lưu. Tiếp đến, chúng tôi đã thực hiện các bước trung gian như: xử lý các phân sai bậc hai, giải quyết các số nguyên đơn trị, hoàn thiện việc bình sai lưới cho một ngày đo, tạo ra file hệ phương trình chuẩn NQ0 phục vụ cho việc xử lý kết hợp cả đợt đo, kiểm tra tính nhất quán của tọa độ các điểm IGS chuẩn, kiểm tra sự liên tục của tọa độ cùng một điểm qua các ngày đo. Kết thúc công việc tính toán của một đợt đo, chúng tôi tính được tọa độ chính xác của các trạm đo và xác định được tốc độ chuyển dịch của các trạm đo đó. Kết quả tính cuối cùng, chúng tôi thực hiện đối sánh với kết quả tính của Viện Thiên văn, Đại hoc Bern, thể hiện qua các Bảng 1a-b và Bảng 2a-b dưới đây. Đối sánh các Bảng 1a và 1b, 2a và 2b, ta thấy có sự trùng khít về các kết quả tính toán của chúng tôi và Viện Thiên văn, Đại học Bern cả về vị trí các điểm cũng như tốc độ chuyển dịch. Sai số trong các tính toán này mà phần mềm thông báo chỉ rất nhỏ, khoảng 0,3-0,4 mm đối với thành phần nằm ngang và 1,2 mm đối với thành phần thẳng đứng [1]. Điều đó thể hiện tính chính xác cũng như việc làm chủ quy trình xử lý, tính toán của phần mềm Bernese phiên bản 5.0 được sử dụng trong tính toán này.

Bảng 1a. Kết quả tính tốc độ chuyển dịch chu kỳ 2002-2003 dạng đầy đủ

(xử lý tại Viện Địa chất)

Bảng 1b. Kết quả tính tốc độ chuyển dịch chu kỳ 2002-2003 dạng đầy đủ
(xử lý tại Viện Thiên văn, Đại hoc Bern)

Bảng 2a. Kết quả tính vận tốc dạng rút gọn (xử lý tại Viện Địa chất)

Bảng 2b. Kết quả tính vận tốc dạng rút gọn (xử lý tại Viện Thiên văn, Đại hoc Bern)

III. TỐC ĐỘ CHUYỂN DỊCH KIẾN TẠO HIỆN ĐẠI TRÊN BIỂN ĐÔNG

Chúng tôi đã tiến hành đo 2 đợt tại các trạm Láng, Bạch Long Vĩ, Song Tử Tây và Côn Đảo trong 2 năm 2007 và 2008. Trong mỗi đợt đo, chúng tôi đã tiến hành đo liên tục 7 ca, mỗi ca 23 giờ 40 phút. Cơ sở dữ liệu được sử dụng trong tính toán này, ngoài dữ liệu của trạm GPS Láng (LANG), Bạch Long Vĩ (BLV1), Song Tử Tây (STT1), Côn Đảo (CDA1), chúng tôi sử dụng dữ liệu đo liên tục của 3 trạm IGS (PIMO, BACO, WUHN) và trạm TCMS, riêng trạm TCMS được sử dụng dữ liệu thuần tuý như trạm đo Láng, Bạch Long Vĩ, Song Tử Tây và Côn Đảo. Các dữ liệu của tổ chức GPS quốc tế phục vụ địa động lực (IGS) như lịch vệ tinh chính xác, mô hình tầng điện ly, file hiệu chỉnh giữa P1-C1, P1-P2 đối với vệ tinh và máy thu, tọa độ cũng như vận tốc chuyển dịch của các trạm IGS trong hệ quy chiếu toàn cầu  IGS05 đã được sử dụng trong tính toán. Khoảng cách gần nhất giữa hai trạm là 223 km (Láng - Bạch Long Vĩ).

1. Chuyển dịch tương đối

Giả sử cố định chuyển dịch tại đảo Bạch Long Vĩ, chúng tôi đã tính được chuyển dịch của Côn Đảo về phía tây với tốc độ 17,9 mm/năm và về phía nam với tốc độ 2,6 mm/năm. Sai số theo phương đông-tây và bắc-nam đều trong khoảng 0,6 mm/năm. Trạm Láng chuyển dịch  về phía đông với tốc độ 11,1 mm/năm, sai số 0,5 mm/năm, chuyển dịch về phía nam với tốc độ 1,2 mm/năm, sai số 0,5 mm/năm. Tốc độ chuyển dịch của Song Tử Tây về phía tây là 15,3 mm/năm với sai số 0,5 mm/năm về phía nam là 5,1 mm/năm với sai số 0,6 mm/năm. Chi tiết tính toán được trình bày trên Bảng 3.

2. Chuyển dịch tuyệt đối

Sử dụng hệ tọa độ toàn cầu ITRF2005/IGS05, với tốc độ đã biết của các trạm IGS  PIMO, BACO, WUHN, chúng ta có thể tính được chuyển dịch tuyệt đối của các trạm đo. Trạm Láng chuyển dịch về phía đông với tốc độ 43,4 mm/năm, sai số 0,5 mm/năm, chuyển dịch về phía nam với tốc độ 8,3 mm/năm, sai số 0,6 mm/năm. Tốc độ chuyển dịch của trạm Bạch Long Vĩ về phía đông là 32,2 mm/năm, sai số 0,5 mm/năm, đồng thời chuyển dịch về phía nam với tốc độ 8,3 mm/năm với sai số 0,6 mm/năm. Trạm Song Tử Tây chuyển dịch về phía đông với tốc độ 17,4 mm/năm, sai số 0,6 mm/năm và chuyển dịch về phía nam với tốc độ 14,2 mm/năm với sai số 0,5 mm/năm. Tốc độ chuyển dịch về phía đông của trạm Côn Đảo là 16,1 mm/năm với sai số 0,6 mm/năm, chuyển dịch về phía nam 12,1 mm/năm với sai số 0,5 mm/năm. Kết quả tính toán chi tiết được trình bày trên Bảng 4.

Bảng 3. Kết quả tính tốc độ chuyển dịch tương đối so với điểm Bạch Long Vĩ cố định

Bảng 4. Kết quả tính tốc độ chuyển dịch tuyệt đối trong hệ toàn cầu IGS05

6
 
Cần lưu ý sai số tính toán theo phần mềm Bernese 5.0 có thể nhỏ hơn sai số thực tế, điều này cần làm sáng tỏ trong những tính toán chi tiết hơn tiếp theo.

IV. ĐẶC ĐIỂM BIẾN DẠNG KIẾN TẠO HIỆN ĐẠI TRÊN BIỂN ĐÔNG

Mặc dù chỉ là kết quả đo lặp của 2 kỳ đo trong khoảng thời gian 2007-2008, chúng ta đã có thể rút ra một số nhận xét về đặc điểm biến dạng của biển Đông Việt Nam.

Tiếp tục với xu thế chuyển dịch về phía đông - đông nam đã quan sát thấy trên đất liền của Việt Nam [7-9], chúng ta quan sát thấy toàn bộ các trạm đo GPS đều chuyển dịch về phía đông - đông nam. Kết quả trên cũng phù hợp với quan sát ở đảo Hải Nam, Quảng Tây, Quảng Đông cũng như toàn rìa đông nam Trung Quốc [1, 4, 10, 15, 16]. Điều này cho thấy biến dạng trên biển Đông chịu sự chi phối chủ yếu của va chạm giữa các mảng Ấn-Úc và Âu-Á.

 

Hình 2. Vị trí, hướng và tốc độ chuyển dịch của các trạm GPS trên biển Đông.

7
 
1. Sự suy giảm tốc độ chuyển dịch theo hướng từ tây sang đông của các trạm đo GPS phía bắc (Láng, Bạch Long Vĩ, Hải Nam) cho thấy, hiện nay vịnh Bắc Bộ bị biến dạng nén và chịu xiết ép theo phương á vĩ tuyến hoặc lệch một chút theo phương TTB-ĐĐN. Trường lực này không thuận lợi cho việc phát triển hệ thống đứt gãy đang hoạt động tách giãn phương á kinh tuyến. Các hướng chính và giá trị chính của trục ứng suất - biến dạng sẽ được chúng tôi chính xác hoá ở các chu kỳ đo sau. Theo tính toán sơ bộ của chúng tôi, tốc độ biến dạng nén tính từ trạm Láng tới trạm Bạch Long Vĩ đạt giá trị xấp xỉ 5.10-8 mm/năm.

2. Phía bắc biển Đông đang đóng lại theo phương TTB-ĐĐN với tốc độ cỡ 77 mm/năm.  Hướng của véc tơ chuyển dịch tại Láng, Bạch Long Vĩ, Hải Nam hầu như ngược với hướng véc tơ chuyển dịch ở PIMO, phản ánh hướng chuyển dịch của mảng Bắc biển Đông cắm dưới Philippines tại trũng Manila về phía Đ-ĐN.  

3. Các trạm đo GPS phía nam (Song Tử Tây, Côn Đảo) có hướng chuyển dịch về phía ĐN cho thấy chế độ địa động lực ở phía nam biển Đông đã thay đổi so với phần phía bắc biển Đông. Tốc độ chuyển dịch ngang nhỏ hơn ở phía bắc. Biển Đông ở phần phía nam không bị đóng lại. Tốc độ biến dạng nhỏ hơn phía bắc Biển Đông.  

4. Mặc dù sai số theo phần mềm Bernese 5.0 thông báo rất nhỏ, nhưng đây chỉ là tính toán sơ bộ của 2 kỳ đo đầu nên véc tơ chuyển dịch sẽ còn được hiệu chỉnh khi bổ sung thêm các trạm IGS khác và các trạm đo ven bờ biển Việt Nam, cũng như thay đổi phương án tính tối ưu hơn.

V. KẾT LUẬN

Mặc dù mới chỉ qua 2 kỳ đo năm 2007-2008, những nét cơ bản về chuyển dịch kiến tạo hiện đại trên biển Đông đã được xác định với tốc độ đóng ở phía bắc biển Đông trong khoảng 77 mm/năm. Hoạt động xiết ép theo phương á vĩ tuyến sẽ cản trở chuyển dịch của các đứt gãy thuận có phương á kinh tuyến. Vai trò của va chạm giữa mảng Ấn-Úc đối với mảng Âu-Á đóng vai trò chủ đạo đối với biến dạng của biển Đông. Hướng chuyển dịch thay đổi từ ĐĐN ở phần phía bắc biển Đông chuyển sang ĐN ở phần phía nam biển Đông. Biến dạng xiết ép giảm ở phần phía nam biển Đông.

Kết quả của chúng tôi còn có ý nghĩa đặc biệt về chủ quyền lãnh thổ của Việt Nam với các tên Bạch Long Vĩ (BLV1), Côn Đảo (CDA1), Song Tử Tây (STT1) từ nay sẽ xuất hiện trên bản đồ địa động lực thế giới.

Bài viết này là kết quả ban đầu của Đề tài trọng điểm cấp nhà nước KC09.11/06-10. Chúng tôi chân thành cám ơn sự giúp đỡ của Bộ Tư lệnh Hải quân, Phòng Công binh, Hải quân Vùng 4, Ban Chỉ huy đảo Song Tử Tây, lãnh đạo các huyện đảo Côn Đảo và Bạch Long Vĩ.

VĂN LIỆU

1. Astronomical Institute, 2004. Tutorial of Bernese GPS Software Version 5.0. University of Bern, Switzerland.

2. Bakiiz H., J. Wang, Y. Chen, 1999. A regional GPS network solution fop monitoring deformations of the Southeastern Eurasian Plate. GPS Solutions, 2/4 : 44-55.

3. Becker M., E. Reinhart, S.B. Nordin, D. Angermann, G. Michel, C. Reigber, 2000. Improving the velocity field in South and Southeast Asia: The third round of GEODYSSEA. Earth Planet Space, 52 : 721-726.

4. Brendan J. Meade, 2007. Present-day kinematics at the India-Asia collision zone. Geol. Soc. of America, 35/1 : 81-84.

5. Habrich H., 2005. New developments in EPN analysis. EUREF Symp. in Vienna, Austria.

6. Kenneth L. Senior, Æ. Jim, R. Ray, Æ. Ronald, L. Beard, 2008. Characterization of periodic variations in the GPS satellite clocks. GPS Solut., DOI 10.1007/s10291-008-0089-9.

7. Michel G.W. et al., 2001. Crustal motion and block behaviour in SE-Asia from GPS measurements. Earth and Planetary Sci. Lett., 187 : 239-244.

8. Nguyễn Tuấn Anh, Phan Trọng Trịnh và nnk., 2007. Báo cáo “Xây dựng hệ thống các điểm trắc địa sử dụng công nghệ GPS độ chính xác cao trong việc quan trắc biến dạng lớp vỏ Trái đất và cảnh báo thiên tai tại khu vực Việt Nam”. Lưu trữ Trung tâm TTKHQG, Hà Nội.

9. Ngô Văn Liêm, Phan Trọng Trịnh, Nguyễn Tuấn Anh, Hoàng Quang Vinh, 2008. Ứng dụng công nghệ GPS trong việc xác định chuyển dịch kiến tạo hiện đại, biến dạng mặt đất và công trình. Địa kỹ thuật, 2 (in press).

10. Niu Z., Wang M. et al., 2005. Contemporary velocity field of crustal movement of Chinese mainland from GPS measurments. Chinese Sci. Bull., 50 : 1-3.

11. Schaer S., D. Ineichen, E. Brockmann, 2005. EUREF LAC analysis at Swisstopo/Code using Bernese Software V5.0. In: Torres J.A. and H. Hornik (Eds). Subcommission for the European Reference Frame (EUREF) Publication in preparation. Vienne.

12. Teferle F. N., Æ.E.J. Orliac, Æ.R.M. Bingley, 2007. An assessment of Bernese GPS software precise point positioning using IGS final products for global site velocities. GPS Solut., 11 : 205-213.

13. Trần Đình Tô, Nguyễn Trọng Yêm, 2004. Chuyển động hiện đại vỏ Trái đất lãnh thổ Việt Nam theo số liệu đo GPS. TC Các khoa học về TĐ, 26/4 : 579-586. Hà Nội.

14. Vigny C., A. Socquet, C. Rangin, N. Chamot-Rooke, M. Pubellier, M.N. Bouin, G. Bertrand, M. Becker, 2003. Present-day crustal deformation around Sagaing fault, Myanmar. J. of Geoph. Res., 108/B11, 2533 doi:10.1029/2002JB001999.

15. Vy Quốc Hải, 2004. So sánh kết quả xử lý số liệu GPS của lưới địa động lực bằng phần mềm GPS 2.35 và Bernese 4.2. TC Các khoa học về TĐ, 26/4 : 418-425. Hà Nội.

16. Wang Qi, Pei-Zhen Zhang, J.T. Freymueller et al., 2001. Present-day crustal deformation in China constrained by GPS measurements. Science, 294.

17. Zhang P.Zh., Zh.K. Shen, M. Wang, W.J. Gan et al., 2004. Continuous deformation of the Tibetan Plateau from GPS data. Geol. Soc. of America,  32/9 : 809-812.