Skip to content

harmonics module

addHarmonicTerms(image)

add Time bands to image

Source code in src\rlcms\harmonics.py 
 6
 7
 8
 9
10
11
def addHarmonicTerms(image):
    """add Time bands to image"""
    timeRadians = image.select("t").multiply(2 * math.pi)
    return image.addBands(timeRadians.cos().rename("cos")).addBands(
        timeRadians.sin().rename("sin")
    )

addTimeConstant(imageCollection, timeField)

Add time constant to images in an ImageCollection timeField: time stamp property name (typically is the 'system:time_start' property of an image)

Source code in src\rlcms\harmonics.py 
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
def addTimeConstant(imageCollection: ee.ImageCollection, timeField: str):
    """
    Add time constant to images in an ImageCollection
    timeField: time stamp property name (typically is the 'system:time_start' property of an image)
    """
    def _(image, timeField):
        # // Compute time in fractional years since the epoch.
        date = ee.Date(image.get(timeField))
        years = date.difference(ee.Date("1970-01-01"), "year")
        # // Return the image with the added bands.
        return image.addBands(ee.Image(years).rename("t").float()).addBands(
            ee.Image.constant(1)
        )

    return imageCollection.map(lambda i: _(i, timeField))

calculateHarmonic(imageCollection, dependent)

Calculate harmonic coefficients (phase and amplitude) off of an ImageCollection dependent: band that you fit the harmonic model for, must be contained in ImageCollection

Source code in src\rlcms\harmonics.py 
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
def calculateHarmonic(imageCollection: ee.ImageCollection, dependent: ee.String):
    """
    Calculate harmonic coefficients (phase and amplitude) off of an ImageCollection
    dependent: band that you fit the harmonic model for, must be contained in ImageCollection
    """
    harmonicIndependents = ee.List(["constant", "t", "cos", "sin"])
    #  Add harmonic terms as new image bands.
    harmonicLandsat = imageCollection.map(addHarmonicTerms)
    # The output of the regression reduction is a 4x1 array image.
    harmonicTrend = harmonicLandsat.select(harmonicIndependents.add(dependent)).reduce(
        ee.Reducer.linearRegression(harmonicIndependents.length(), 1)
    )

    # Turn the array image into a multi-band image of coefficients.
    harmonicTrendCoefficients = (
        harmonicTrend.select("coefficients")
        .arrayProject([0])
        .arrayFlatten([harmonicIndependents])
    )

    # // Compute phase and amplitude.
    phase = (
        harmonicTrendCoefficients.select("cos")
        .atan2(harmonicTrendCoefficients.select("sin"))
        .rename(dependent.cat("_phase"))
    )

    amplitude = (
        harmonicTrendCoefficients.select("cos")
        .hypot(harmonicTrendCoefficients.select("sin"))
        .rename(dependent.cat("_amplitude"))
    )
    return ee.Image.cat(phase, amplitude)

doHarmonicsFromOptions(imgColl, **kwargs)

calculateHarmonic function band-wise

kwargs: harmonicsOptions (dict): which band(s) and the DOY start and end date to compute harmonics on formatted like: { 'red':{'start':1,'end':365}, 'blue':{'start':1,'end':365} } returns: ee.Image containing bands [band_phase, band_amplitude]

Source code in src\rlcms\harmonics.py 
 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
def doHarmonicsFromOptions(imgColl:ee.ImageCollection,**kwargs):
    """
    calculateHarmonic function band-wise

    args:
        imgColl (ee.ImageCollection)
    kwargs:
        harmonicsOptions (dict): which band(s) and the DOY start and end date to compute harmonics on
            formatted like: 
            {
                'red':{'start':1,'end':365},
                'blue':{'start':1,'end':365}
                }
    returns:
        ee.Image containing bands [band_phase, band_amplitude]
    """
    imgColl = ee.ImageCollection(imgColl)
    harmonicsOptions = kwargs['harmonicsOptions']

    # get harmonicsOptions dictionary
    if isinstance(harmonicsOptions,dict):

        # get band keys as list
        bands = ee.Dictionary(harmonicsOptions).keys()

        def harmonicByBand(band):
            band = ee.String(band)
            # get the params for that band
            bandwiseParams = ee.Dictionary(harmonicsOptions).get(band)

            # get the start and end DOY parameters
            start = ee.Dictionary(bandwiseParams).get('start')
            end = ee.Dictionary(bandwiseParams).get('end')

            # create temporal filtered imgColl for that band
            imgCollByBand = (ee.ImageCollection(imgColl)
                                .select(band)
                                .filter(ee.Filter.dayOfYear(start,end)))
            # add time bands
            timeField = "system:time_start"
            timeCollection = addTimeConstant(imgCollByBand, timeField)

            return ee.Image(calculateHarmonic(timeCollection,band))
    else:
        raise TypeError(f"harmonicsOptions expects dict type, got: {type(harmonicsOptions)}")

    # do harmonics by band key in model_inputs dictionary
    listOfImages = ee.Image.cat(ee.List(bands).map(harmonicByBand))
    bandStack = ee.Image(ee.ImageCollection.fromImages(listOfImages).toBands())

    # to remove srcImg band name indexing resulting from .toBands() 
    # (i.e. [0_swir1_phase, 0_swir1_amplitude] -> [swir1_phase, swir1_amplitude] )
    bandNames = bandStack.bandNames()
    fixedBandNames = bandNames.map(lambda e: ee.String(e).split("_").slice(-2).join("_"))
    return bandStack.rename(fixedBandNames)

harmonicRGB(harmonics)

Use the HSV to RGB transform to display phase and amplitude

Source code in src\rlcms\harmonics.py 
48
49
50
51
52
53
54
55
56
57
58
59
def harmonicRGB(harmonics: ee.Image):
    """Use the HSV to RGB transform to display phase and amplitude"""
    amplitude = harmonics.select(".*amplitude")
    phase = harmonics.select(".*phase")

    rgb = (
        phase.unitScale(-math.pi, math.pi)
        .addBands(amplitude.multiply(2.5))
        .addBands(ee.Image(1))
        .hsvToRgb()
    )
    return rgb