一维统计量，多维统计量，OIF， python程序

• 计算需要用到的库

from PIL import Image
import matplotlib.pyplot as plt
import scipy.stats as sts
import numpy as np
import cv2 as cv
from numpy import cov,corrcoef
from scipy.special import comb, perm
from itertools import combinations, permutations

• 一维统计量

1. mean    均值
2. ptp       极差
3. var       方差
4. std       标准差
5. min  max   极值

代码如下:

img1 = np.array(cv.imread("X:/WDR/ETM/LE71220352017080EDC00_B1.TIF",1))
img2 = ....
img7 = ....


IMGN2P=[img1,img2,img3,img4,img5,img7,img8]   //  二维矩阵数组



img1np=np.array(img1).flatten()   
img2np=.....
img7np=......

IMGNP=[img1np,img2np,img3np,img4np,img5np,img7np,img8np]   //  一维矩阵数组

// 创建 数组

mean=[0.00,1.00,2.00,3.00,4.00,5.00,6.00]   # 均值数组
ptp=range(7) #   极差
var=range(7) #  方差
std=range(7)# 标准差

正文代码：

mean=np.array(mean)          ###  均值矩阵
ptp=np.array(ptp,dtype=float)   ###   极差矩阵
var=np.array(var,dtype=float)   ###   方差矩阵
std=np.array(std,dtype=float)   ###   标准差
for x in range(7):
    mean[x]=IMGN2P[x].mean()
    var[x]=IMGN2P[x].var()
    std[x]=IMGN2P[x].std()
    ptp[x]=IMGN2P[x].ptp()
for x in range(7):      ### 最大值，最小值5

    print("各波段  max: ",IMGNP[x].max())  
    print("各波段  min: ",IMGNP[x].min())

• 多维统计量

1. cov  协方差
2. corrcoef  相关系数
3. OIF   指数因子

先写单一协方差和相关系数 还有OIF :

#cov corrcof   ###### 协方差矩阵           矩阵
 data1 =  np.array([img1np,img2np]) ##波段1.2
 cove = cov(data1,bias=1)

 print("波段1 2  协方差:",cove)
 print("波段1 2 相关系数: ",corrcoef(data1))

###1 2 3 波段OIF
 Sk=std[0]+std[1]+std[2]

 data2= np.array([img1np,img3np])  ## 波段 1  3
 data3=np.array([img2np,img3np])  ##波段 2 3
 R1=abs(corrcoef(data1).min())
 R2=abs(corrcoef(data2).min())
 R3=abs(corrcoef(data3).min())
 R=R1+R2+R3
 OIF=Sk/R
 print("OIF(1 2 3 ) : \n",OIF)

7波段协方差，相关系数和OIF：

# ncb73=comb(7,3)
nstd=list(range(100))
# ncb72=comb(7,2)
ncorrcoef=range(21)

 

k=0
npcrroef=np.zeros((7,7),dtype=float)
for x in range(6):   ###相关系数
    mxk=x
    for y in range(1,6-x):
        mxk=mxk+1
        ndata=np.array([IMGNP[x],IMGNP[mxk]])   ### ju zhen  qiu xiang guan xi shu 相关系数矩阵  12  13 ....
        npcrroef[x][mxk]=abs(corrcoef(ndata).min())
        k=k+1
m2=0
mk=0

for x in  range(5):
    m2=m2+1
    men=m2
    for y in range(5-m2):
        men=men+1
        nstd[mk]=std[x]+std[m2]+std[men]
for x in range(35):
    if x < 5:
        nco = npcrroef[0][x+2]+npcrroef[0][1]+npcrroef[1][x+2]
        OIF = nstd[x]/nco
        print("OIF( 1 2 ",x+3," )",OIF)
    elif x < 9:
        nco = npcrroef[0][2]+npcrroef[0][x-2]+npcrroef[2][x-2]
        OIF = nstd[x]/nco
        print("OIF(1 3 ",x-1," )",OIF)
    elif x < 12:
        nco = npcrroef[0][3]+npcrroef[0][x-5]+npcrroef[3][x-5]
        OIF = nstd[x]/nco
        print("OIF ( 1 4 ",x-4," )",OIF)
    elif x < 14:
        nco = npcrroef[0][4]+npcrroef[0][x-7]+npcrroef[4][x-7]
        OIF= nstd[x]/nco
        print("OIF ( 1 5 ",x-7," )",OIF)
    elif x <15:
        nco = npcrroef[0][5]+npcrroef[0][6]+npcrroef[5][6]
        OIF= nstd[x]/nco
        print("OIF ( 1 6 7)",OIF)
    elif x < 19:
        nco = npcrroef[1][2]+npcrroef[1][x-12]+npcrroef[2][x-12]
        OIF=nstd[x]/nco
        print("OIF ( 2 3 ",x-11," )",OIF)
    elif x < 22:
        nco=npcrroef[1][3]+npcrroef[1][x-15]+npcrroef[3][x-15]
        OIF=nstd[x]/nco
        print("OIF ( 2 4 ",x-14," )",OIF)
    elif x< 24:
        nco=npcrroef[1][4]+npcrroef[1][x-17]+npcrroef[4][x-17]
        OIF=nstd[x]/nco
        print("OIF (2 5 ",x-16," )",OIF)
    elif x <25:
        nco=npcrroef[1][5]+npcrroef[1][x-18]+npcrroef[5][x-18]
        OIF=nstd[x]/nco
        print("OIF (2 6 7)",OIF)
    elif x < 28:
        nco = npcrroef[2][3]+npcrroef[2][x-21]+npcrroef[3][x-21]
        OIF=nstd[x]/nco
        print("OIF (3 4 ",x-20,OIF)
    elif x < 30:
         nco=npcrroef[2][4]+npcrroef[2][x-23]+npcrroef[4][x-23]
         OIF=nstd[x]/nco
         print("OIF (3 5 ",x-22," )",OIF)
    elif x < 31:
        nco=npcrroef[2][5]+npcrroef[2][6]+npcrroef[5][6]
        OIF=nstd[x]/nco
        print("OIF （3 6 7）",OIF)
    elif x < 33:
        nco=npcrroef[3][4]+npcrroef[3][x-26]+npcrroef[4][x-26]
        OIF=nstd[x]/nco
        print("OIF (4 5 ",x-25," )",OIF)
    elif x<34:
        nco=npcrroef[3][5]+npcrroef[3][6]+npcrroef[5][6]
        OIF=nstd[x]/nco
        print("OIF (4 6 7)",OIF)
    elif x <35:
        nco=npcrroef[4][5]+npcrroef[4][6]+npcrroef[5][6]
        OIF=nstd[x]/nco
        print("OIF( 5 6 7)",OIF)

若要运行OIF，需先求得标准差  ！ (std)

PS： 求7波段OIF代码，个人认为非常冗杂，也没想到其他办法 去求，如若有知道的好友，麻烦告知！！！