Linear interpolation of spectral response function

0 Preface

Domestic satellite data , The official spectral response function is generally interval 1nm The data of , But some atmospheric correction models may require data at other intervals , such as 6s need 2.5nm Interval data , Therefore, it is necessary to interpolate the data . Of course, there are thousands of ways to realize this function , This paper tries to use python Linear interpolation processing , There are other ways to add .

1 Content

1.1 Spectral response function

Most of the spectral response functions of domestic satellites can be downloaded through China Resources Satellite Center , Download address

1.2 np.interp() Interpolation processing

1.2.1 np.interp() brief introduction

numpy There's a function in interp Linear interpolation can be realized , The main parameters are briefly introduced below , Please refer to Official information

numpy.interp(x, xp, fp, left=None, right=None, period=None)
The main parameters are introduced ：
x For the data to be inserted x Axis ,xp For the original data x Axis data ,fp For the original data y Axis data

1.2.2 Code

import numpy as np
import matplotlib.pyplot as plt
import pandas as pd

wavelength = pd.read_excel('GF.xlsx', 'PMS1')
wavelength.head()

length_interp_bin = np.arange(450, 891, 2.5)  # Intercept segment 
x = wavelength.iloc[:,0].tolist()
wavelength_interp = []
for i in range(2,6):
    wavelength_interp_temp = np.interp(length_interp_bin, x, wavelength.iloc[:,i])
    wavelength_interp.append(wavelength_interp_temp)
    
wavelength_interp = np.array(wavelength_interp)

plt.figure(figsize=(20,10))
plt.plot(x, wavelength.iloc[:,2:6], '.')
plt.plot(length_interp_bin, wavelength_interp.T, 'x')
plt.show()

#  Read the interpolation data within each band respectively 
b1 = wavelength_interp[0, np.argwhere(length_interp_bin==450)[0][0]: np.argwhere(length_interp_bin==520)[0][0]+1]
b2 = wavelength_interp[1, np.argwhere(length_interp_bin==520)[0][0]: np.argwhere(length_interp_bin==590)[0][0]+1]
b3 = wavelength_interp[2, np.argwhere(length_interp_bin==630)[0][0]: np.argwhere(length_interp_bin==690)[0][0]+1]
b4 = wavelength_interp[3, np.argwhere(length_interp_bin==770)[0][0]: np.argwhere(length_interp_bin==890)[0][0]+1]

1.3 Use scipy interpolation

1.3.1 Code

Top up

#  Import package 
from scipy import interpolate

plt.figure(figsize=(15,8))
for kind in ["nearest", "zero", "slinear", "quadratic", "cubic"]:
    # "nearest", "zero" Interpolate for steps 
    # slinear  linear interpolation 
    # "quadratic", "cubic"  by 2 rank 、3 rank B Spline interpolation 
    for i in range(2,6):
        f = interpolate.interp1d(x,  wavelength.iloc[:,i], kind=kind)
        inter_wavelength = f(length_interp_bin)
        plt.plot(length_interp_bin, inter_wavelength,'x', label=str(kind+'_band'+str(i)))
plt.legend()
plt.show()

end