Please provide python code, thanks!

Construct a piecewise cubic polynomial function that interpolates the N+1 data points (x0,y0),,(xN,yN). In each subinterval Ik=[xk,xk+1], define the interpolant as a cubic polynomial sk(x)=a0(k)+a1(k)x+a2(k)x2+a3(k)x3. The cubic polynomials sk are constructed such that for k=0,1,,N3,sk(x) should interpolate points (xk,yk),(xk+1,yk+1),(xk+2,yk+2),(xk+3,yk+3). Implement a function in Python with the following arguments: - x is a row vector containing the N+1 values x=[x0,,xN]. - y is a row vector containing the N+1 values y=[y0,,yN]. The function should return a (N24) matrix M with the coefficient of each sk on the respective row, i.e., M=a0(0)a0(1)a0(N3)a1(0)a1(1)a1(N3)a2(0)a2(1)a2(N3)a3(0)a3(1)a3(N3). Then, use the code to print the associated matrix M, given an arbitrary input N and data range, from the sixth-order polynomial P(x)=924x62772x5+3150x41680x3+420x242x+1