هشدار: رویه‌های شرح داده شده در این آموزش برای هماهنگ‌سازی داده‌های ETM+ و OLI لندست قدیمی هستند و هنگام کار با داده‌های بازتاب سطحی مجموعه 2 لندست توصیه یا ضروری نیستند. برای اطلاعات بیشتر، به بخش سوالات متداول در مورد هماهنگ‌سازی لندست مراجعه کنید. داده‌های مجموعه 1 منسوخ شده‌اند و کد موجود در این آموزش دیگر کار نخواهد کرد یا به‌روزرسانی نخواهد شد.
این آموزش مربوط به هماهنگ‌سازی بازتاب سطحی ماهواره ETM+ لندست با بازتاب سطحی ماهواره OLI لندست است. این آموزش موارد زیر را ارائه می‌دهد:
- یک تابع تبدیل طیفی
- توابعی برای ایجاد داده‌های آماده برای تحلیل
- یک مثال مصورسازی سری‌های زمانی
این راهنما به عنوان یک راهنمای جامع برای هماهنگ‌سازی و مصورسازی سری‌های زمانی منطقه‌ای بیش از ۳۵ سال از داده‌های لندست در نظر گرفته شده است که می‌تواند بلافاصله در منطقه (یا مناطق) مورد نظر شما اعمال شود.
توجه داشته باشید که ضرایب تبدیل ETM+ به OLI برای TM نیز اعمال می‌شود. به همین دلیل، در این آموزش، اشاره به ETM+ مترادف با TM است.
درباره لندست
لندست یک برنامه تصویربرداری ماهواره‌ای است که از سال ۱۹۷۲ تصاویر زمین با وضوح متوسط را جمع‌آوری می‌کند. به عنوان طولانی‌ترین برنامه رصد زمین مبتنی بر فضا، این برنامه یک رکورد زمانی ارزشمند برای شناسایی روندهای مکانی-زمانی در تغییرات چشم‌انداز ارائه می‌دهد. داده‌های ابزار نقشه‌بردار موضوعی (TM)، نقشه‌بردار موضوعی پیشرفته پلاس (ETM+) و تصویرگر عملیاتی زمین (OLI) در این آموزش استفاده می‌شوند. آنها ارتباط نزدیکی با هم دارند و ادغام آنها در یک سری زمانی ثابت که یک رکورد پیوسته از سال ۱۹۸۴ تا به امروز، با ریتم ۱۶ روز برای هر حسگر با وضوح مکانی ۳۰ متر، تولید می‌کند، نسبتاً آسان است. ابزار اسکنر چندطیفی، رکورد لندست را تا سال ۱۹۷۲ گسترش می‌دهد، اما در این آموزش استفاده نمی‌شود. داده‌های آن کاملاً متفاوت است و ادغام با حسگرهای بعدی را چالش برانگیز می‌کند. برای مثال‌هایی از هماهنگی در بین همه حسگرها، به Savage و همکاران (۲۰۱۸) و Vogeler و همکاران (۲۰۱۸) مراجعه کنید.
چرا هماهنگ‌سازی
روی و همکاران (۲۰۱۶) نشان می‌دهند که بسته به کاربرد، تفاوت‌های کوچک اما بالقوه معناداری بین ویژگی‌های طیفی لندست ETM+ و OLI وجود دارد. دلایلی که ممکن است بخواهید مجموعه داده‌ها را هماهنگ کنید عبارتند از: تولید یک سری زمانی طولانی که Landsat TM، ETM+ و OLI را در بر می‌گیرد، تولید کامپوزیت‌های درون سالانه نزدیک به تاریخ برای کاهش اثرات مشاهدات از دست رفته از شکاف‌های SLC-off ETM+ و پوشش ابر/سایه، یا افزایش فراوانی مشاهدات در یک سری زمانی. لطفاً برای اطلاعات بیشتر به نسخه خطی پیوند داده شده در بالا مراجعه کنید.
دستورالعمل‌ها
توابع
در ادامه مجموعه‌ای از توابع مورد نیاز برای هماهنگ‌سازی ETM+ با OLI و تولید داده‌های آماده برای تحلیل ارائه شده است که در بخش کاربرد این آموزش برای مصورسازی تاریخچه طیفی-زمانی یک پیکسل استفاده خواهند شد.
هماهنگ‌سازی
هماهنگ‌سازی از طریق تبدیل خطی فضای طیفی ETM+ به فضای طیفی OLI طبق ضرایب ارائه شده در جدول 2 ضرایب رگرسیون OLS روی و همکاران (2016) حاصل می‌شود. ضرایب مربوط به باند در فرهنگ لغت زیر با ثابت‌های تصویر شیب ( slopes ) و عرض از مبدا ( itcps ) تعریف شده‌اند. توجه داشته باشید که مقادیر عرض از مبدا y در 10000 ضرب می‌شوند تا با مقیاس‌بندی داده‌های بازتاب سطحی USGS Landsat مطابقت داشته باشند.
var coefficients = {
itcps: ee.Image.constant([0.0003, 0.0088, 0.0061, 0.0412, 0.0254, 0.0172])
.multiply(10000),
slopes: ee.Image.constant([0.8474, 0.8483, 0.9047, 0.8462, 0.8937, 0.9071])
};
نام باندهای ETM+ و OLI برای یک پنجره پاسخ طیفی یکسان، یکسان نیستند و نیاز به استانداردسازی دارند. توابع زیر فقط باندهای بازتاب و باند pixel_qa را از هر مجموعه داده انتخاب می‌کنند و آنها را بر اساس محدوده طول موجی که نشان می‌دهند، تغییر نام می‌دهند.
// Function to get and rename bands of interest from OLI.
function renameOli(img) {
return img.select(
['B2', 'B3', 'B4', 'B5', 'B6', 'B7', 'pixel_qa'],
['Blue', 'Green', 'Red', 'NIR', 'SWIR1', 'SWIR2', 'pixel_qa']);
}
// Function to get and rename bands of interest from ETM+.
function renameEtm(img) {
return img.select(
['B1', 'B2', 'B3', 'B4', 'B5', 'B7', 'pixel_qa'],
['Blue', 'Green', 'Red', 'NIR', 'SWIR1', 'SWIR2', 'pixel_qa']);
}
در نهایت، تابع تبدیل را تعریف کنید که مدل خطی را بر روی داده‌های ETM+ اعمال می‌کند، نوع داده را برای سازگاری با OLI به Int16 تبدیل می‌کند و باند pixel_qa را برای استفاده بعدی در پوشش ابر و سایه دوباره متصل می‌کند.
function etmToOli(img) {
return img.select(['Blue', 'Green', 'Red', 'NIR', 'SWIR1', 'SWIR2'])
.multiply(coefficients.slopes)
.add(coefficients.itcps)
.round()
.toShort()
.addBands(img.select('pixel_qa'));
}
پوشش ابر و سایه
داده‌های آماده برای تجزیه و تحلیل باید دارای ابرها و سایه‌های ابر باشند که ماسک شده‌اند. تابع زیر از باند pixel_qa مربوط به CFmask ( Zhu et al., 2015 ) که در هر تصویر بازتاب سطحی Landsat USGS وجود دارد، برای تنظیم پیکسل‌های شناسایی شده به عنوان ابر و سایه ابر به حالت تهی استفاده می‌کند.
function fmask(img) {
var cloudShadowBitMask = 1 << 3;
var cloudsBitMask = 1 << 5;
var qa = img.select('pixel_qa');
var mask = qa.bitwiseAnd(cloudShadowBitMask)
.eq(0)
.and(qa.bitwiseAnd(cloudsBitMask).eq(0));
return img.updateMask(mask);
}
محاسبه شاخص طیفی
برنامه‌ی کاربردیِ پیشِ رو از شاخص طیفیِ نسبت سوختگیِ نرمال‌شده (NBR) برای نمایش تاریخچه‌ی طیفیِ یک پیکسل جنگلیِ تحت تأثیر آتش‌سوزیِ جنگلی استفاده می‌کند. NBR به این دلیل استفاده می‌شود که کوهن و همکاران (۲۰۱۸) دریافتند که در میان ۱۳ شاخص/باند طیفی، NBR بیشترین نسبت سیگنال به نویز را در رابطه با اختلال (سیگنال) جنگل در سراسر ایالات متحده دارد. این شاخص با تابع زیر محاسبه می‌شود.
function calcNbr(img) {
return img.normalizedDifference(['NIR', 'SWIR2']).rename('NBR');
}
در کاربرد خودتان، می‌توانید از شاخص طیفی متفاوتی استفاده کنید. اینجا جایی است که می‌توانید یک شاخص اضافی را تغییر دهید یا اضافه کنید.
ترکیب توابع
برای هر مجموعه داده، یک تابع پوششی تعریف کنید که تمام توابع فوق را برای راحتی در اعمال آنها به مجموعه تصاویر مربوطه، تجمیع کند.
// Define function to prepare OLI images.
function prepOli(img) {
var orig = img;
img = renameOli(img);
img = fmask(img);
img = calcNbr(img);
return ee.Image(img.copyProperties(orig, orig.propertyNames()));
}
// Define function to prepare ETM+ images.
function prepEtm(img) {
var orig = img;
img = renameEtm(img);
img = fmask(img);
img = etmToOli(img);
img = calcNbr(img);
return ee.Image(img.copyProperties(orig, orig.propertyNames()));
}
در برنامه خود، ممکن است تصمیم بگیرید که توابع را شامل یا حذف کنید. این توابع را در صورت نیاز تغییر دهید.
مثال سری زمانی
توابع پوششی prepOli و prepEtm فوق را می‌توان روی مجموعه‌های بازتاب سطحی لندست نگاشت کرد تا داده‌های آماده برای تجزیه و تحلیل بین حسگری ایجاد شود و گاهشماری طیفی یک پیکسل یا ناحیه‌ای از پیکسل‌ها را تجسم کند. در این مثال، شما یک سری زمانی بیش از 35 سال ایجاد خواهید کرد و تاریخچه طیفی یک پیکسل واحد را نمایش خواهید داد. این پیکسل خاص، تاریخچه اخیر یک قطعه جنگل مخروطی بالغ شمال غربی اقیانوس آرام (شکل 1) را که در دهه‌های 1980 و 90 دچار برخی اختلالات و در سال 2011 دچار یک سوختگی با شدت بالا شده است، مرتبط می‌کند.
شکل ۱. موقعیت مکانی و ویژگی‌های سایت به عنوان مثال منطقه مورد نظر. جنگل مخروطی بالغ شمال غربی اقیانوس آرام در دامنه شمالی کوه هود، اورگان، ایالات متحده آمریکا. تصاویر با اجازه از: گوگل ارث، سرویس جنگلداری USDA، لندست و کوپرنیک.
تعریف حوزه مورد علاقه (AOI)
نتیجه این برنامه یک سری زمانی از مشاهدات لندست برای یک پیکسل است. یک شیء ee.Geometry.Point برای تعریف مکان پیکسل استفاده می‌شود.
var aoi = ee.Geometry.Point([-121.70938, 45.43185]);
در برنامه خود، می‌توانید پیکسل دیگری را انتخاب کنید. این کار را با جایگزینی مختصات طول و عرض جغرافیایی فوق انجام دهید. به عنوان یک جایگزین، می‌توانید تاریخچه طیفی گروهی از پیکسل‌ها را با استفاده از تعاریف دیگر شیء ee.Geometry ، مانند ee.Geometry.Polygon() و ee.Geometry.Rectangle() خلاصه کنید. لطفاً برای اطلاعات بیشتر به بخش هندسه‌ها در راهنمای توسعه‌دهندگان مراجعه کنید.
بازیابی مجموعه‌های حسگر لندست
مجموعه‌های بازتاب سطحی Landsat USGS را برای OLI ، ETM+ و TM دریافت کنید.
برای کسب اطلاعات بیشتر در مورد هر مجموعه داده، به لینک‌ها مراجعه کنید.
var oliCol = ee.ImageCollection('LANDSAT/LC08/C01/T1_SR');
var etmCol = ee.ImageCollection('LANDSAT/LE07/C01/T1_SR');
var tmCol = ee.ImageCollection('LANDSAT/LT05/C01/T1_SR');
تعریف فیلتر برای مجموعه تصاویر
یک فیلتر تعریف کنید تا مجموعه تصاویر را بر اساس مرز مکانی منطقه مورد نظر، فصل اوج فتوسنتز و کیفیت محدود کند.
var colFilter = ee.Filter.and(
ee.Filter.bounds(aoi), ee.Filter.calendarRange(182, 244, 'day_of_year'),
ee.Filter.lt('CLOUD_COVER', 50), ee.Filter.lt('GEOMETRIC_RMSE_MODEL', 10),
ee.Filter.or(
ee.Filter.eq('IMAGE_QUALITY', 9),
ee.Filter.eq('IMAGE_QUALITY_OLI', 9)));
این معیارهای فیلترینگ را به دلخواه خود تغییر دهید. سایر ویژگی‌هایی را که می‌خواهید فیلتر کنید، در زیر برگه «ویژگی‌های تصاویر» در صفحات توضیحات داده که در بالا پیوند داده شده‌اند، مشخص کنید.
مجموعه‌ها را آماده کنید
مجموعه‌ها را ادغام کنید، فیلتر را اعمال کنید و تابع prepImg را روی همه تصاویر نگاشت کنید. نتیجه قطعه کد زیر یک ee.ImageCollection واحد است که شامل تصاویری از مجموعه‌های حسگر OLI، ETM+ و TM است که معیارهای فیلتر را برآورده می‌کنند و برای NBR آماده تجزیه و تحلیل پردازش می‌شوند.
// Filter collections and prepare them for merging.
oliCol = oliCol.filter(colFilter).map(prepOli);
etmCol = etmCol.filter(colFilter).map(prepEtm);
tmCol = tmCol.filter(colFilter).map(prepEtm);
// Merge the collections.
var col = oliCol.merge(etmCol).merge(tmCol);
یک نمودار سری زمانی ایجاد کنید که تمام مشاهدات را نشان دهد
هماهنگی بین حسگرها را می‌توان به صورت کیفی با رسم تمام مشاهدات به عنوان یک نمودار پراکندگی که در آن حسگر با رنگ تشخیص داده می‌شود، ارزیابی کرد. تابع ui.Chart.feature.groups در Earth Engine این ابزار را فراهم می‌کند، اما ابتدا باید مقدار NBR مشاهده شده به عنوان یک ویژگی به هر تصویر اضافه شود. یک تابع کاهش ناحیه را روی مجموعه تصویر نگاشت کنید که میانه تمام پیکسل‌های متقاطع با AOI را محاسبه کرده و نتیجه را به عنوان یک ویژگی تصویر تنظیم می‌کند.
var allObs = col.map(function(img) {
var obs = img.reduceRegion(
{geometry: aoi, reducer: ee.Reducer.median(), scale: 30});
return img.set('NBR', obs.get('NBR'));
});
در برنامه خودتان، اگر AOI شما یک منطقه بزرگ است، می‌توانید افزایش scale و تعیین آرگومان‌های bestEffort ، maxPixels و tileScale برای تابع reduceRegion را در نظر بگیرید تا مطمئن شوید که عملیات از محدودیت‌های حداکثر پیکسل، حافظه یا زمان‌بندی تجاوز نمی‌کند. همچنین می‌توانید کاهنده میانه را با آمار مورد نظر خود جایگزین کنید. برای اطلاعات بیشتر به بخش آمار یک منطقه تصویر در راهنمای توسعه‌دهندگان مراجعه کنید.
اکنون این مجموعه می‌تواند توسط تابع ui.Chart.feature.groups پذیرفته شود. این تابع انتظار دارد که یک مجموعه و نام ویژگی‌های شیء مجموعه به نمودار نگاشت شوند. در اینجا، ویژگی 'system:time_start' به عنوان متغیر محور x و 'NBR' به عنوان متغیر محور y عمل می‌کند. توجه داشته باشید که ویژگی NBR توسط تابع reduceRegion در بالا بر اساس نام باند در تصویر در حال کاهش نامگذاری شده است. اگر از باند متفاوتی (نه 'NBR') استفاده می‌کنید، آرگومان نام ویژگی محور y را بر این اساس تغییر دهید. در نهایت، متغیر گروه‌بندی (series) را روی 'SATELLITE' تنظیم کنید که به آن رنگ‌های متمایز اختصاص داده می‌شود.
var chartAllObs =
ui.Chart.feature.groups(allObs, 'system:time_start', 'NBR', 'SATELLITE')
.setChartType('ScatterChart')
.setSeriesNames(['TM', 'ETM+', 'OLI'])
.setOptions({
title: 'All Observations',
colors: ['f8766d', '00ba38', '619cff'],
hAxis: {title: 'Date'},
vAxis: {title: 'NBR'},
pointSize: 6,
dataOpacity: 0.5
});
print(chartAllObs);
پس از مدتی پردازش، نموداری مشابه شکل ۲ در کنسول ظاهر می‌شود. به چند نکته توجه کنید:
- چندین مشاهده در سال وجود دارد.
- سری زمانی از سه حسگر مختلف تشکیل شده است.
- هیچ گونه سوگیری حسگریِ قابل حذفی وجود ندارد.
- فراوانی مشاهدات سالانه متفاوت است.
- مقداری واریانس درون سالانه وجود دارد که در طول سال‌ها ثابت و نسبتاً کم است. امتحان کردن این کار با NDVI (img.normalizedDifference(['SWIR1', 'NIR']).rename('NDVI');) تغییرپذیری پاسخ درون و بین سالانه بسیار بیشتری را نشان می‌دهد.
- NBR برای اکثر سری‌های زمانی با برخی اختلالات جزئی همچنان بالا است.
- کاهش سریع و زیاد در پاسخ NBR ناشی از آتش‌سوزی جنگل است که در ... مشهود است.
- با این حال، توجه داشته باشید که آتش‌سوزی ( دریاچه دلار ) در سپتامبر ۲۰۱۱ رخ داده است. خطا در سال تشخیص به دلیل محدوده تاریخ ترکیبی سالانه از ژوئیه تا آگوست است؛ تغییراتی که پس از این دوره رخ می‌دهند تا مشاهده غیر تهی بعدی که می‌تواند سال بعد یا بعد باشد، انتخاب نمی‌شوند.
- بازیابی پوشش گیاهی (افزایش پاسخ NBR) دو سال پس از شناسایی از دست رفتن عمده NBR آغاز می‌شود.
شکل 2. نمودار سری زمانی پاسخ طیفی برای یک پیکسل واحد با نمایش از Landsat TM، ETM+ و OLI. تصاویر TM و ETM+ با تبدیل خطی با OLI هماهنگ شده‌اند. تصاویر از ماه جولای تا آگوست هستند و برای کیفیت بالا فیلتر شده‌اند.
یک نمودار سری زمانی رسم کنید که میانگین سالانه را نشان دهد
برای ساده‌سازی و حذف نویز از سری زمانی، می‌توان از کاهش مشاهدات درون‌سالانه استفاده کرد. در اینجا از میانه استفاده شده است.
اولین قدم، گروه‌بندی تصاویر بر اساس سال است. با نگاشت تنظیمات مجموعه 'year' از ee.Image.Date هر تصویر، یک ویژگی 'year' جدید به هر تصویر اضافه کنید.
var col = col.map(function(img) {
return img.set('year', img.date().get('year'));
});
ویژگی جدید 'year' می‌تواند برای گروه‌بندی تصاویر از یک سال با استفاده از عمل پیوند استفاده شود. پیوند بین یک مجموعه تصویر مجزا از سال و مجموعه کامل، لیست‌های سال مشابه را ارائه می‌دهد که از طریق آنها می‌توان کاهش میانه را انجام داد. مجموعه کامل را به مجموعه‌ای از نمایندگان سال مجزا تقسیم کنید.
var distinctYearCol = col.distinct('year');
یک فیلتر و یک اتصال تعریف کنید که سال‌های منطبق بین مجموعه سال مجزا ( distinctYearCol ) و مجموعه کامل ( col ) را شناسایی کند. تطابق‌ها در یک ویژگی به نام 'year_matches' ذخیره می‌شوند.
var filter = ee.Filter.equals({leftField: 'year', rightField: 'year'});
var join = ee.Join.saveAll('year_matches');
اتصال را اعمال کنید و FeatureCollection حاصل را به ImageCollection تبدیل کنید.
var joinCol = ee.ImageCollection(join.apply(distinctYearCol, col, filter));
// Apply median reduction among matching year collections.
var medianComp = joinCol.map(function(img) {
var yearCol = ee.ImageCollection.fromImages(img.get('year_matches'));
return yearCol.reduce(ee.Reducer.median())
.set('system:time_start', img.date().update({month: 8, day: 1}));
});
اگر ترجیح می‌دهید آماری غیر از میانه برای نمایش مجموعه تصاویر در یک سال مشخص استفاده شود، ee.Reducer بالا را تغییر دهید.
در نهایت، نموداری ایجاد کنید که مقادیر میانه NBR را از مجموعه مشاهدات تابستانی درون سالانه نشان دهد. آماره کاهش منطقه را به دلخواه تغییر دهید.
var chartMedianComp = ui.Chart.image
.series({
imageCollection: medianComp,
region: aoi,
reducer: ee.Reducer.median(),
scale: 30,
xProperty: 'system:time_start',
})
.setSeriesNames(['NBR Median'])
.setOptions({
title: 'Intra-annual Median',
colors: ['619cff'],
hAxis: {title: 'Date'},
vAxis: {title: 'NBR'},
lineWidth: 6
});
print(chartMedianComp);
پس از مدتی پردازش، نموداری مشابه شکل ۳ در کنسول ظاهر می‌شود.
شکل 3. نمودار سری زمانی پاسخ طیفی میانه تابستانی برای یک پیکسل واحد با نمایش از Landsat TM، ETM+ و OLI. تصاویر TM و ETM+ با یک تبدیل خطی با OLI هماهنگ شده‌اند. تصاویر از ماه جولای تا آگوست هستند و برای کیفیت بالا فیلتر شده‌اند.
اطلاعات تکمیلی
توابع تبدیل جایگزین
روی و همکاران (۲۰۱۶) جدول ۲ ضرایب رگرسیون OLS و RMA را برای تبدیل بازتاب سطحی ETM+ به بازتاب سطحی OLI و برعکس ارائه می‌دهد. آموزش فوق فقط تبدیل ETM+ به OLI را با استفاده از رگرسیون OLS نشان می‌دهد. توابع مربوط به همه گزینه‌های ترجمه را می‌توان در اسکریپت ویرایشگر کد یا اسکریپت منبع GitHub یافت.
آتش‌سوزی دریاچه دلار
آتش‌سوزی ذکر شده در این آموزش، آتش‌سوزی دریاچه دلار است که در اواسط سپتامبر ۲۰۱۱ در دامنه‌های شمالی کوه هود، اورگان، ایالات متحده آمریکا رخ داد. برای اطلاعات بیشتر در مورد آن، به این پست وبلاگ مراجعه کنید و به مقاله وارنر و همکاران (۲۰۱۵) مراجعه کنید.