-
Notifications
You must be signed in to change notification settings - Fork 14
/
MODIS_TFA_Day_or_Night_Export_2_Assets
162 lines (138 loc) · 6.47 KB
/
MODIS_TFA_Day_or_Night_Export_2_Assets
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
// This code uses Temporal Fourier Analysis to build 365 images, for each day of year,
// of the climatological tempratures based on the day (or night) temperatures.
// Choose here either day or night band
var Temperature_Band = 'LST_Day_1km';
// var Temperature_Band = 'LST_Night_1km';
// Insert here the coordinates of the area of interest.
var geometry = ee.Geometry.Rectangle([33.2,29.0,36.6,34.0]);
var geometry_json = geometry.toGeoJSON();
// Change the path to your own assets folder
var Assets_path = 'users/shilosh/MODIS_LST_TFA'
var collection = 'MODIS/006/MYD11A1';
// The number of cycles per year to model.
var harmonics = 3;
// Make a list of harmonic frequencies to model.
// These also serve as band name suffixes.
var harmonicFrequencies = ee.List.sequence(1, harmonics);
// Function to get a sequence of band names for harmonic terms.
var getNames = function(base, list) {
return ee.List(list).map(function(i) {
return ee.String(base).cat(ee.Number(i).int());
});
};
// Construct lists of names for the harmonic terms.
var cosNames = getNames('cos_', harmonicFrequencies);
var sinNames = getNames('sin_', harmonicFrequencies);
var k2celsius = function(image) {
return image.multiply(ee.Image(0.02))
.subtract(ee.Image(273.15)) //convert Kelvin to Celsius
.set('system:time_start', image.get('system:time_start'));
};
// Function to compute the specified number of harmonics
// and add them as bands.
// for p = 0, Nharmonics-1 do begin
// a[p] = a[p] + 2*TS_t*cos(2*!pi*(p+1)*(t+1)/N)
// b[p] = b[p] + 2*TS_t*sin(2*!pi*(p+1)*(t+1)/N)
var addHarmonics = function(freqs) {
return function(image) {
// Make an image of frequencies.
var frequencies = ee.Image.constant(freqs);
// This band should represent lst.
var lst = ee.Image(image).select(Temperature_Band);
// This band should represent days from start.
var time = ee.Image(image).select('t');
// Get the cosine terms.
var cosines = lst.multiply(2.0).multiply((frequencies.multiply(2.0 * Math.PI).multiply(time).divide(365.0)).cos())
.rename(cosNames);
// Get the sin terms.
var sines = lst.multiply(2.0).multiply((frequencies.multiply(2.0 * Math.PI).multiply(time).divide(365.0)).sin())
.rename(sinNames);
return image.addBands(cosines).addBands(sines);
};
};
var LST = ee.ImageCollection(ee.List.sequence(1, 365).map(function (doy){
return ee.ImageCollection(collection)
.select(Temperature_Band)
.filter(ee.Filter.calendarRange(doy, doy, 'day_of_year'))
// .filter(ee.Filter.eq('start_hour', 12))
.map(k2celsius)
.mean()
.set({doy: doy})
.addBands(ee.Image(ee.Number(doy)).rename('t'))
// .map(addHarmonics(harmonicFrequencies));
}))
print('LST = ', LST)
LST = LST.map(addHarmonics(harmonicFrequencies));
// a = a / N_TS
// b = b / N_TS
var num_days = 365;//ee.Date(lastDay).difference(ee.Date(firstDay), 'day');
var a_coef_1 = LST.select('cos_1').sum().divide(num_days);
var a_coef_2 = LST.select('cos_2').sum().divide(num_days);
var a_coef_3 = LST.select('cos_3').sum().divide(num_days);
var b_coef_1 = LST.select('sin_1').sum().divide(num_days);
var b_coef_2 = LST.select('sin_2').sum().divide(num_days);
var b_coef_3 = LST.select('sin_3').sum().divide(num_days);
// ; ========================================
// ; prepare the first 3 harmonics of the TFA
// ; ========================================
// H0 = total(TS) / N_TS ;mean
// FTA = MAKE_ARRAY(N_TS, /FLOAT, VALUE=0.)
// H = MAKE_ARRAY(Nharmonics, N_TS, /FLOAT, VALUE=0.)
// omega = MAKE_ARRAY(Nharmonics, /FLOAT, VALUE=0.)
// Phase = MAKE_ARRAY(Nharmonics, /FLOAT, VALUE=0.)
// Amplitude = MAKE_ARRAY(Nharmonics, /FLOAT, VALUE=0.)
// for p = 0, Nharmonics-1 do begin
// omega[p] = 2. * !pi * (p+1) / N
// Phase[p] = atan( -b[p], a[p] )
// Amplitude[p] = sqrt( a[p]^2 + b[p]^2 )
var H0 = LST.select(Temperature_Band).mean();
// Get the omega terms.
var omega_1 = ee.Image(1.0).multiply(2.0 * Math.PI).divide(365.0);
var omega_2 = ee.Image(2.0).multiply(2.0 * Math.PI).divide(365.0);
var omega_3 = ee.Image(3.0).multiply(2.0 * Math.PI).divide(365.0);
var phase_1 = a_coef_1.atan2(b_coef_1.multiply(ee.Image(-1.0)));
var phase_2 = a_coef_2.atan2(b_coef_2.multiply(ee.Image(-1.0)));
var phase_3 = a_coef_3.atan2(b_coef_3.multiply(ee.Image(-1.0)));
var amplitude_1 = (a_coef_1.pow(2).add(b_coef_1.pow(2))).sqrt();
var amplitude_2 = (a_coef_2.pow(2).add(b_coef_2.pow(2))).sqrt();
var amplitude_3 = (a_coef_3.pow(2).add(b_coef_3.pow(2))).sqrt();
// for t = 0, N_TS-1 do begin
// H[p,t] = Amplitude[p] * cos(omega[p]*(t+1) + Phase[p])
// Function to add a H band.
var addH = function(image) {
var H_1 = amplitude_1.multiply((omega_1.multiply(image.select('t')).add(phase_1)).cos()).rename('H_1');
var H_2 = amplitude_2.multiply((omega_2.multiply(image.select('t')).add(phase_2)).cos()).rename('H_2');
var H_3 = amplitude_3.multiply((omega_3.multiply(image.select('t')).add(phase_3)).cos()).rename('H_3');
// FTA[t] = H0 + H[0,t] + H[1,t] + H[2,t]
var TFA = H0.add(H_1).add(H_2).add(H_3).rename('TFA');
var Anomalies = TFA.subtract(image.select(Temperature_Band)).rename('Anomalies');
return image.addBands(TFA).addBands(Anomalies);
};
var LST_TFA = LST.map(addH);
// Define the visualization parameters.
var vizParams = {
bands: ['mean', 'phase', 'amplitude'],
min: [285.0,2.0,5.0],
max: [295.0,3.0,10.0],
gamma: [1.0, 1.0, 1.0]
};
//display the image.
var img = ee.Image(H0.addBands(phase_1).addBands(amplitude_1).rename('mean', 'phase', 'amplitude'));
// Map.addLayer(img, vizParams, 'RGB');
//Iterating over the image collection using this function....
var TFA_Images = LST_TFA.select('TFA')
.sort('t')
.iterate(function(img, all) {
return ee.Image(all).addBands(img);
}, ee.Image().select());
print('TFA_Images = ', TFA_Images );
// Export one single file with all the bands
Export.image.toAsset({
image: ee.Image(TFA_Images),
assetId: Assets_path + Temperature_Band,
region: geometry_json,
crs: ee.Image('MODIS/006/MYD11A1/2018_01_01').projection().crs().getInfo(),
crsTransform: [926.625433056,0,-20015109.354,0,-926.625433055,10007554.677],
description: 'exportToAsset-MODIS-TFA' + Temperature_Band,
pyramidingPolicy: {'.default': 'sample'}
});