如何将 bash 命令集成到 python 代码中

Question

大家，我希望将 bash 命令集成到我的 python 代码中以计算索引。我的问题是，我想为每个计算出的指数生成带波段的输出图像，但我无法通过 bash 命令将这些指数整合到用 [= 创建的 'im_index' 矩阵中16=]代码。我不知道如何 link 他们两个......你有什么想法吗？

import numpy as np
import sys
import os
import spectral as sp
from scipy import ndimage
import pylab as pl
from math import *
import spectral.io.envi as envi

#------------------------------------
def reject_outliers(data, m=1):
    return data[abs(data - np.mean(data)) < m * np.std(data)]

#------------------------------------
def find_nearest(array, value):
    #For a given value, find the nearest value in an array

    array = np.asarray(array)
    idx = (np.abs(array - value)).argmin()

    return idx
#------------------------------------

#Open existing dataset
src_directory = "/d/afavro/Bureau/4_reflectance/"
dossier = os.listdir (src_directory)
print(dossier)

for fichier in dossier:
    print (fichier)

    ssrc_directory = "/d/afavro/Bureau/4_reflectance/" + fichier + "/"
    rasters = os.listdir (ssrc_directory) 
    print(rasters) 

    OUTPUT_FOLDER = "/d/afavro/Bureau/5_indices2/" + 'indices_' + fichier + '/'
    print(OUTPUT_FOLDER)
    if not os.path.exists(OUTPUT_FOLDER):
        os.makedirs(OUTPUT_FOLDER)

    for image in rasters:
        print(image)
            name, ext = os.path.splitext(image)

        if ext == '.hdr':

            img = sp.open_image(ssrc_directory + image)
            print(image)
            im_HS = img[:,:,:]

            cols = im_HS.shape[0]  # Number of column
            rows = im_HS.shape[1]   # Number of lines
            bands = im_HS.shape[2]  # Number of bands
            NbPix = cols * rows   # Number of pixels

            #Get wavelengths from hdr file
            wv = np.asarray(img.bands.centers)

            if len(wv) == 0 :
                print("Wavelengths not defined in the hdr file")
                sys.exit("Try again!")

            if wv[0] > 100:
                wv=wv*0.001 # Convert to micrometers if necessary

            im_HS=im_HS.reshape(NbPix, bands)

 #Compute HC index------------------------------------------------------
            Nind=4 # Number of indice to be computed
            im_index=np.zeros((cols*rows, Nind))
            names = []

 ##NDVI computation
            names.append('NDVI')
            bande_ref=[0.67, 0.8]
            bRef0 = find_nearest(wv,bande_ref[0])
            bRef1 = find_nearest(wv,bande_ref[1])
            #Check if the required specral bands are available
            if (np.abs(wv[bRef0]-bande_ref[0])<=0.1 and np.abs(wv[bRef1]-bande_ref[1])<=0.1):
                b0 = im_HS[:, bRef0]
                b1 = im_HS[:, bRef1]  
                index = (b0 - b1) /  (b0 + b1)
            else:
                index = np.zeros_like(im_HS[:,0])
                print("Wavelengths selection problem, NDVI not computed")

            im_index[:,0]=  index

    # bash command :

            inRaster = ssrc_directory + image
            print(inRaster)
            outRaster = OUTPUT_FOLDER + 'indices_' + image
            print (outRaster)
            cmd = 'otbcli_RadiometricIndices -in inRaster -list Soil:BI Vegetation:MSAVI Vegetation:SAVI -out outRaster'         
            os.system(cmd)

#saving
            im_index=im_index.reshape(cols, rows, Nind)

            file_image = OUTPUT_FOLDER + "indices2_" + fichier 

            header = envi.read_envi_header(ssrc_directory + image)
            header ['description'] = "fichier d'origine " + image
            header ['band names'] = ['NDVI', 'Sober filter', 'NDWI', 'IB(1)', 'IB(2)']
            del header['wavelength units'] 
            del header['wavelength'] 

            sp.envi.save_image(file_image + '.hdr', im_index, metadata=header, force = True, interleave = 'bsq')

Answer 1

我想您可能正在寻找 subprocess 包。一个例子：

>>> import subprocess as sp
>>> output = sp.check_output('echo hello world', shell=True)
>>> print(output)
b'hello world\n'

check_output() 方法可用于从命令收集标准输出。之后您需要解析输出以获得整数索引。

Answer 2

假设这是您实际询问的代码：

inRaster = ssrc_directory + image
print(inRaster)
outRaster = OUTPUT_FOLDER + 'indices_' + image
print (outRaster)
cmd = 'otbcli_RadiometricIndices -in inRaster -list Soil:BI Vegetation:MSAVI Vegetation:SAVI -out outRaster'
os.system(cmd)

当然，单引号内的inRaster只是一个文字串；插入变量的值你可以说

cmd = 'otbcli_RadiometricIndices -in ' + inRaster + \
     ' -list Soil:BI Vegetation:MSAVI Vegetation:SAVI -out ' + \
     outRaster

或

cmd = 'otbcli_RadiometricIndices -in {0} -list Soil:BI Vegetation:MSAVI Vegetation:SAVI -out {1}'.format(
   inRaster, outRaster)

或 Python 中的许多其他字符串插值技术（遗留 % 格式、f 字符串等）。但更好的解决方案是用更灵活和通用的 subprocess 替换 os.system，甚至在 os.system 文档中也有建议。

subprocess.run([
    'otbcli_RadiometricIndices',
    '-in', inRaster,
    '-list', 'Soil:BI', 'Vegetation:MSAVI', 'Vegetation:SAVI',
    '-out', outRaster], check=True)

subprocess.run 是在 Python 3.5 中引入的；如果您需要与旧版本兼容，请尝试 subprocess.check_call 甚至粗略的 subprocess.call.

如何将 bash 命令集成到 python 代码中

How to integrate a bash command into a python code

python

bash

indices