📎 🤹🏿 🔊 نمذجة GPR 🎅🏾 ✍🏿 🔀

إن Georadar (جهاز تقني إذاعي للسبر تحت السطحي ، GPR ، رادار اختراق الأرض) ، والذي يستخدم حاليًا على نطاق واسع - بدءًا من تعيين ثقوب الأرانب ودراسة السحالي إلى البحث عن الألغام ، لا يزال متعة باهظة الثمن إلى حد ما.

شاشة GPR (إطار من برنامج التلفاز البريطاني "Command of the Time")

لكن من الممكن تقييم قدراتها ودراسة جوانب مختلفة من تفاعل المجال الكهرومغناطيسي لجيورادار مع البيئة دون الحصول على / استئجار جهاز "حديدي". سوف تساعدنا حزمة gprMax ( gpr - من اختصار GPR ، Max - الأحرف الأولى من اسم James Clerk Maxwell ، الذي وضع أسس الديناميكا الكهربائية) ، الموزعة تحت رخصة GNU GPL v3.
مؤلفو هذا المشروع ، الذي بدأ في عام 1996 ، هم كريج وارن من جامعة نورثمبريا وأنتونيس جيانوبولوس من جامعة أدنبرة. تم تطوير الحزمة في الأصل في C ثم إعادة كتابتها بالكامل في مجموعات Python 3 / Cython.

يتطلب تثبيت الحزمة برامج مجمعة مثبتة تدعم OpenMP (أدوات إنشاء Microsoft Visual C ++ 2015 (يوصى بهذا الإصدار!) لنظام Windows / gcc لنظام التشغيل Linux) ، ومكتبة NumPy ، ومترجم Cython. بعد التنزيل من المستودع على GitHub وفك الشفرة المصدرية للمشروع ، انتقل إلى المجلد الجذر وقم بتنفيذ الأوامر:

python setup.py build python setup.py install

"كبداية سريعة" ، نفكر في العمل مع الحزمة باستخدام مثال بسيط ثنائي الأبعاد - هوائي الإرسال T لرادار نابض (نبضة GPR) يصدر نبضة كهرمغنطيسية ، يصل جزء من طاقتها مباشرة إلى هوائي الاستقبال R في شكل موجة مباشرة (DW - موجة مباشرة) ، وبعضها يخترق عبر الرمل ، ينعكس من سطح أسطوانة التوصيل ويصل هوائي الاستقبال في شكل موجة عاكسة (الموجة المنعكسة (RW)):

تنسيق ملف الإدخال
قم بإنشاء مجلد الطُرز في المجلد الجذر للمشروع ، والذي نضع فيه hello.in الملف النصي الذي يحتوي على الأوامر لتنفيذ المحاكاة (الأوامر التالية تتوافق مع الإصدار الحالي (الثالث) من المشروع).

كل فريق لديه النموذج:

 #: _1 _2 _3 ...

يمكن كتابة أمر واحد فقط على سطر واحد ، ويجب أن يكون الحرف الأول من السطر الذي يحتوي على الأمر #.

يمكن أن تكون الأوامر مصحوبة بتعليقات:

##

ترتيب الأوامر مهم للأوامر الخاصة ببناء الكائنات - يتم تنفيذ هذه الأوامر بالترتيب الذي تظهر به في ملف الإدخال.

شكل نبض

النبض الكهرومغناطيسي المنبعث من georadar يستمر بضع كسور من النانو ثانية ، وعلاوة على ذلك ، غالبا ما تستخدم ثلاثة أشكال النبض:

فترة موجة جيبية واحدة (جيب)
زخم غاوسي (غاوسي)
تتناسب القبعة المكسيكية (ricker) مع المشتق الثاني للدالة الغوسية ، يشبه شكل الموجة الدافعة هذه سمبريرو (استخدم هذا النموذج النبضي بواسطة عالم الجيوفيزيائي الأمريكي نورمان ريكير في عام 1953 لدراسة الإشارات السيزمية)

على سبيل المثال ، نختار نبضة غوسية (نوع النبض - غاوسي) بتردد مركزي

f_{c} = 1 G H z = 1 c d o t 10^{9} H z

$f_c = 1 \\ GHz = 1 \ cdot 10 ^ 9 Hz$ يحددها الأمر:

 #waveform: gaussian 1 1e9 pulse

(1 - سعة النبض الشرطي ، تسمية النبض النبضي)

في هذه الحالة ، يتم وصف الزخم المستخدم في المحاكاة بواسطة التعبير:

W (t) = e ^ {- 2 \ cdot {\ pi} ^ 2 \ cdot {f_c} ^ 2 {(t- {1 \ over {f_c}}) ^ 2}

$W (t) = e ^ {- 2 \ cdot {\ pi} ^ 2 \ cdot {f_c} ^ 2 {(t- {1 \ over {f_c}}) ^ 2}$

نموذج البيئة ونظام الإحداثيات

في النمذجة ثنائية الأبعاد ، تنقسم المساحة المدروسة إلى خلايا ذات حجم معين ، ويظهر نظام إحداثيات النموذج على هذا النحو - محاور X و Y تشكل المستوى المحسوب (مع عرض

w

$w$ و طويل القامة

ح

$ح$ ) ، طول النموذج على طول المحور Z له قيمة مساوية لخطوة أخذ العينات

د ل ت ا

$\ دلتا$ .

عند اختيار خطوة أخذ العينات ، يمكنك اتباع قاعدة الإبهام - يجب ألا يتجاوز حجم الخطوة عُشر أصغر طول موجي تمت دراسته في النموذج ("10 خلايا لكل طول موجي"):

D e l t a l e 0 ، 1 c d o t {l a m b d a}_{m i n}

$\ Delta \ le 0،1 \ cdot {\ lambda} _ {min}$

لتحديد طول الموجة ، يحتاج المرء إلى معرفة الحد الأقصى للتردد الذي يؤخذ في الاعتبار في طيف الإشارة المنبعثة وسرعة الموجة في الوسط قيد النظر.

سرعة الانتشار للموجة الكهرومغناطيسية في وسط مع ثابت عازل نسبي

{e p s i l o n}_{r}

${\ epsilon} _r$ ، معبراً عنه بالسنتيمتر لكل نانو ثانية - تحدد وحدة السرعة المعتمدة في الرادار بالتعبير:

v = 30 o v e r s q r t {e p s i l o n}_{r}

$v = {30 \ over {\ sqrt {{\ epsilon} _r}}}$

يتم تحديد الطول الموجي بالسنتيمتر من خلال التعبير:

l a m b d a = v o v e r f

$\ lambda = {v \ over f}$

(

و

$و$ - التردد في غيغاهرتز).
لعرض شكل النبض الغوسي وطيفه ، يمكنك استخدام الأمر:

 python -m tools.plot_source_wave gaussian 1 1e9 5e-9 1e-12 -fft

(gaussian - نوع النبضة ، 1 - سعة النبضة ، 1e9 - التردد المركزي (1 GHz) ، 5e-9 - مدة عرض النبضة (5 ns) ، 1e-12 - الخطوة الزمنية (1 ps))

وفقا لطيف نبض غاوسي مع

f_{s} = 1 G H z

$f_s = 1 \ GHz$ تحديد أن في -40 ديسيبل التردد

f a p p r o x 3 G H z

$f \ approx 3 \ GHz$ .
في هذا المثال ، الوسيط الذي يقع فيه الكائن المجس ، نختار الرمال الجافة ذات السماحية النسبية

{e p s i l o n}_{r} = 3

${\ epsilon} _r = 3$ .
العثور على سرعة انتشار الموجة الكهرومغناطيسية في الرمال:

v = 30 o v e r s q r t 3 = 17.3 c m / n s

$v = {30 \ over {\ sqrt {3}}} = 17.3 \ cm / ns$

حدد الطول الموجي في الرمال:

l a m b d a = v o v e r f = 17.3 o v e r 3 = 5.8 س م = 58 م م د و ل ا

$\ lambda = {v \ over f} = {17.3 \ over {3}} = 5.8 سم = 58 مم دولا$

بناءً على ذلك ، نختار الخطوة التي هي نفسها لجميع المحاور (

D e l t a = {D e l t a}_{x} = {D e l t a}_{y} = {D e l t a}_{z}

$\ Delta = {\ Delta} _x = {\ Delta} _y = {\ Delta} _z$ ) ويساوي 2 مم = 0.002 م لأسباب الراحة (عدد صحيح من الخطوات يناسب 1 سم):

 #dx_dy_dz: 0.002 0.002 0.002

قصر المساحة المحاكاة على مستطيل العرض

w

$w$ يساوي 80 سم = 0.8 متر وارتفاع

ح

$ح$ يساوي 60 سم = 0.6 م:

 #domain: 0.60 0.60 0.002

(بالنسبة للنموذج ثنائي الأبعاد ، يجب الإشارة إلى سمك يساوي خطوة واحدة (0.002))

يحدد حجم منطقة المحاكاة وحجم الخطوة المكانية عدد خلايا النموذج ، وبالتالي ، متطلبات ذاكرة الكمبيوتر.

وصفنا الرمال مع التوصيل محددة

س ي ج م ا = 0.01 c m / m

$\ سيجما = 0.01 \ \ cm / m$ وثابت عازل ثابت

{e p s i l o n}_{r} = 3

${\ epsilon} _r = 3$ الأوامر:

 #material: 3 0.01 1 0 sand

(1 يتوافق مع النفاذية المغناطيسية النسبية

{m u}_{r}

${\ mu} _r$ مساوية للوحدة (بدون خواص مغناطيسية) ، 0 - لا خسارة مغناطيسية ، وتسمية الرمل لهذه المادة).

املأ الرمل بمعظم المساحة المحاكاة (من y = 0 إلى y = 38 cm = 0.38 m):

 #box: 0 0 0 0.80 0.38 0.002 sand

(0 0 0 - إحداثيات الزاوية اليسرى السفلى ، 0.80 0.38 0.002 - إحداثيات الزاوية اليمنى العليا (0.002 - خطوة أخذ العينات))
الباقي هو مساحة حرة (تسمية free_space) ، مكافئة تقريبًا في خصائص الهواء (

{e p s i l o n}_{r} = 1

${\ epsilon} _r = 1$ .

{m u}_{r} = 1

${\ mu} _r = 1$ .

s i g m a = 0

${\ sigma} = 0$ ).
يتم عرض حدود منطقة المحاكاة على أنها شرط حدود الامتصاص (ABC).

كهدف ، نختار اسطوانة من موصل مثالي (يعكس الإشعاع الكهرومغناطيسي بالكامل) بنصف قطر 6 سم = 0.06 م مع مركز يقع عند نقطة مع إحداثيات س = 25 سم = 0.25 م و ص = 10 سم = 0.1 م:

 #cylinder: 0.250 0.1 0 0.250 0.1 0.002 0.06 pec

(بيك هو مادة موصلة تماما)

هوائيات

تم تجهيز GPR المحاكاة بهوائيين - الإرسال والاستقبال.
في دراسة الحالة الخاصة بنا ، تخيل هوائي الإرسال مع ثنائي القطب Hertz طوله يساوي خطوة أخذ العينات (في حالة ثلاثية الأبعاد ، يمكنك تحديد هوائي من مكتبة واسعة) ملقاة في الرمال (العمل في اتصال مع المتوسط سبر) على مسافة 5 سم إلى يسار وسط المنطقة (س = 35 سم = 0.35 م ، ص = 38 سم = 0.38 م):

 #hertzian_dipole: z 0.35 0.380 0.0 pulse

(z هو محور الاستقطاب ثنائي القطب (للحالة ثنائية الأبعاد (وضع TMD 2D) z فقط صالح) ، والنبض هو تسمية شكل النبضة التي يشعها الهوائي)

يقع هوائي الاستقبال عادةً على مسافة ثابتة صغيرة من جهاز الاستقبال ، والذي يُسمى قاعدة وحدة الهوائي (يُسمى هذا الخيار الخاص بالموضع النسبي للهوائيات "الإزاحة المشتركة"). كقاعدة ، حدد مسافة 10 سم ، وبالتالي فإن الإحداثي الأفقي 35 + 10 = 45 سم = 0.45 م (5 سم إلى يمين منتصف المنطقة):

 #rx: 0.45 0.38 0.0

فاصل المحاكاة

يتم تحديد اختيار إطار زمني للنمذجة بحيث يكون للإشارة المنعكسة من الهدف وقت للوصول إلى هوائي الاستقبال.

نحدد الوقت التقريبي المطلوب للإشارة في الحالة قيد النظر ، مع أخذ المسافة من الهوائيات إلى الهدف

ح = 18 س م د و ل ا

$ح = 18 \ \ سم دولا$ :

t a p p r o x 2 c d o t h o v e r v = 2 c d o t 18 o v e r 17.3 = 2.1 n s

$t \ approx {2 \ cdot h \ over {v}} = {2 \ cdot 18 \ over 17.3} = 2.1 \ \ ns$
نظرًا لأن الجزء العلوي من نبضة غاوس الرادار بتردد مركزي قدره 1 جيجاهرتز قد تم تغييره بالنسبة لبداية المحور الزمني بمقدار 1 نانوثانية ، فإننا نختار نافذة زمنية تبلغ 5 نانوثانية:

 #time_window: 5e-9

تصميم

وبالتالي ، فإن محتويات ملف الإدخال هي كما يلي:

hello.in

 #domain: 0.80 0.60 0.002 #dx_dy_dz: 0.002 0.002 0.002 #time_window: 5e-9 #material: 3 0.01 1 0 sand #waveform: gaussian 1 1e9 pulse #hertzian_dipole: z 0.35 0.380 0.0 pulse #rx: 0.45 0.38 0.0 #box: 0 0 0 0.80 0.38 0.002 sand #cylinder: 0.40 0.1 0 0.40 0.1 0.002 0.06 pec

نبدأ عملية النمذجة:

 python -m gprMax models\hello.in

لإجراء المحاكاة ، يتم استخدام طريقة المجال الزمني للفرق المحدود (FDTD ، المجال الزمني للفرق المحدود) (تم اقتراح الخوارزمية الأساسية للطريقة بواسطة Kane Yee) ، وبعدها يتم تسمية الخلايا التي قسمت النموذج إليها بخلايا Yee . يتم الحصول على حل رقمي في المجال الزمني عن طريق حل معادلات ماكسويل لكل خلية.

في الحالة ثنائية الأبعاد (وضع TMZ 2D) ، يتم حساب المكون فقط

E_{z}

$E_z$ المجال الكهربائي ومكوناته

H_{x}

$H_x$ و

H_{y}

$H_y$ المجال المغناطيسي.

في حالة تجاوز مقدار الذاكرة المتوفرة ، يتم إصدار رسالة حول نقص الذاكرة لإجراء المحاكاة:

إذا كان أي من الكائنات خارج نطاق المحاكاة ، يتم عرض رسالة خطأ:

لإجراء المحاكاة باستخدام النموذج الموصوف ، استغرق الأمر حوالي 56 ميغابايت فقط من ذاكرة الوصول العشوائي (إذا قمت بتقليل الخطوة بمقدار النصف - حتى 1 مم - تزداد متطلبات الذاكرة إلى 99 ميجابايت).

بعد اكتمال المحاكاة ، يظهر ملف hello.out في مجلد الطرز ، والذي يحتوي على نتائج المحاكاة بتنسيق HDF5 ، المصمم خصيصًا لتخزين البيانات الرقمية.

بناء المسار

لتصور النتائج ، نقوم بتصميم الآثار:

 python -m tools.plot_Ascan models\hello.out

يعرض كل مسار (A-scan) يتم تقديمه في النافذة التي يتم فتحها رسمًا بيانيًا لأحد مكونات الحقل الكهرومغناطيسي في موقع هوائي الاستقبال:

إن الموجة المباشرة التي تأتي مباشرة من هوائي الإرسال (DW) والموجة المنعكسة من الهدف (RW) تكون مرئية على المسارات.

يرتبط المحور الأفقي للوقت بعمق الهدف الذي يعكس الإشارة من خلال سرعة الموجة الكهرومغناطيسية في المادة.

ولكن ماذا يحدث إذا وضعت أمر الأسطوانة أمام أمر الصندوق في ملف الإدخال؟

اختفت الإشارة المنعكسة من الأسطوانة - امتص الرمل الأسطوانة (هذا مثال على أهمية الترتيب الذي تم به بناء الكائنات).

بناء الشخصية

لكن الأكثر إفادة هو الراداروغرام (radargram) - المظهر الجانبي (B-scan) ، وهو مزيج من العديد من المسارات التي بنيت عند تحريك georadar على طول اتجاه معين - وهذا هو الإجراء ذاته لتحريك عربة مع رادار على طول المنطقة المدروسة.

تغيير وصف الهوائيات عن طريق تحريكها إلى بداية المحور الأفقي:

 #hertzian_dipole: z 0.10 0.38 0.0 pulse #rx: 0.20 0.38 0.0

وضعنا الخطوة لتحريك الهوائيات تساوي 1 سم = 0.01 م:

 #src_steps: 0.01 0 0 #rx_steps: 0.01 0 0

وبالتالي ، فإن محتويات ملف الإدخال هي كما يلي:

hello.in

 #domain: 0.80 0.60 0.002 #dx_dy_dz: 0.002 0.002 0.002 #time_window: 5e-9 #material: 3 0.01 1 0 sand #waveform: gaussian 1 1e9 pulse #hertzian_dipole: z 0.10 0.38 0.0 pulse #rx: 0.20 0.38 0.0 #src_steps: 0.01 0 0 #rx_steps: 0.01 0 0 #box: 0 0 0 0.80 0.38 0.002 sand #cylinder: 0.40 0.1 0 0.40 0.1 0.002 0.06 pec

قم بتشغيل المحاكاة في وضع الدُفعات:

 python -m gprMax models\hello.in -n 50

(50 هو عدد الخطوات التي يتحرك بها الرادار).

بعد البدء ، يتم إجراء المحاكاة بالتسلسل لمواقع GPR 50:

بعد نهاية المحاكاة ، يوجد 50 ملفًا hello1.out ... hello50.out في مجلد الطراز.
ادمج هذه الملفات في ملف hello_merged.out مع الأمر:

 python -m tools.outputfiles_merge models/hello

بناء ملف تعريف:

 python -m tools.plot_Bscan models\hello_merged.out Ez

(Ez هو مكون المجال الكهرومغناطيسي الذي نبني به الملف الشخصي - المكون المحول مباشرة إلى الجهد)

كما ترون ، يتم عرض الشريط الأفقي الناجم عن الموجة المباشرة في الملف الشخصي أعلاه ، والإشارة التشعبية المميزة الناتجة عن الموجة المنعكسة وتظهر التغيير في المسافة من الهدف إلى GPR عند نقلها أدناه.

توضح وسيلة الإيضاح الموجودة على اليمين الألوان المطابقة وشدة المجال.

E_{z}

$E_z$ (كما هو موضح على المسار):
أحمر -

E_{z} > 0

$E_z> 0$
أبيض -

E_{z} = 0

$E_z = 0$
الأزرق -

E_{z} < 0

$E_z <0$
من خلال تحليل هذه الملفات الشخصية ، يمكننا استخلاص استنتاجات حول عمق وحجم وحتى شكل الأهداف ، كما تستخدم الشبكات العصبية لهذا الغرض.

طرق وملف تعريف على شاشة GPR (إطار من برنامج التلفاز البريطاني "Command of the Time")

لكن الانعكاس يحدث ليس فقط من الأجسام الموصلة ، ولكن أيضًا عند حدود طبقتين مع ثوابت عازلة مختلفة.

قم بإنشاء طبقة رملية ثانية مع عازل ثابت في الجزء السفلي من النموذج

{e p s i l o n}_{r} = 9

${\ epsilon} _r = 9$ :

hello.in

 #domain: 0.8 0.6 0.002 #dx_dy_dz: 0.002 0.002 0.002 #time_window: 5e-9 #material: 3 0.01 1 0 sand_1 #material: 9 0.01 1 0 sand_2 #waveform: gaussian 1 1e9 pulse #hertzian_dipole: z 0.10 0.38 0.0 pulse #rx: 0.20 0.38 0.0 #src_steps: 0.01 0 0 #rx_steps: 0.01 0 0 #box: 0 0 0 0.80 0.38 0.002 sand_1 #box: 0 0 0 0.80 0.10 0.002 sand_2

كما يمكن أن يرى ، أسفل "التتبع" للموجة المباشرة (DW) ، ظهر قطاع خطي من الموجة (RW) ينعكس من واجهة الطبقتين الرمليتين.

خارج نطاق هذا المثال ، تتميز ميزات gprMax مثل النمذجة ثلاثية الأبعاد والأشكال السطحية المعقدة ونماذج الهوائي المفصلة والمحاسبة عن تشتت الموجات الكهرومغناطيسية ... علاوة على ذلك ، يمكن استخدام الحزمة ليس فقط لمحاكاة تشغيل GPR ، ولكن أيضًا لدراسة انتشار الموجات الكهرومغناطيسية في بيئات مختلفة ، و تسريع العمليات الحسابية باستخدام تقنية NVIDIA CUDA ، والتي توفر زيادة عشرة أضعاف في السرعة مقارنة بالحوسبة المتوازية على وحدة المعالجة المركزية باستخدام OpenMP. لزيادة مرونة النمذجة ، يمكنك وضع كتل الكود في Python في ملف الإدخال.

بعض الأمثلة على استخدام حزمة gprMax:

روابط مفيدة:

الموقع الرسمي gprMax
دليل مستخدم GprMax
gprMax على يوتيوب

نمذجة GPR

More articles: