【问题标题】:How to integrate a bash command into a python code [duplicate]如何将 bash 命令集成到 python 代码中[重复]
【发布时间】:2019-05-06 07:51:14
【问题描述】:

大家, 我希望将 bash 命令集成到我的 python 代码中以计算索引。我的问题是我想为每个计算的索引生成一个带带的输出图像,但是我无法通过 bash 命令将这些索引集成到使用我的 python 代码创建的“im_index”矩阵中。我不知道如何链接它们...你有什么想法吗?

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')

【问题讨论】:

标签: python bash indices


【解决方案1】:

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

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-string 等)。但更好的解决方案是将os.system 替换为更灵活、更通用的subprocess,甚至在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

【讨论】:

  • 详情和更多选项,请见stackoverflow.com/questions/4256107/…
  • 非常感谢,它有效!但是,我的输出文件是 img 和 img.aux.xml 格式,是否可以有一个 .img 和一个 .hdr?我想告诉他保持与另一个相同的标题,但有一些修改
  • 恐怕我不太明白你的后续问题。这听起来有点像一个简单的字符串操作问题; stackoverflow.com/questions/678236/…
  • 我只想输出一个头文件而不是 .img.aux.xml
  • 你在问如何string.replace('.img.aux.xml', '.hdr')
【解决方案2】:

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

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

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

【讨论】:

    猜你喜欢
    • 1970-01-01
    • 2019-03-13
    • 2014-04-08
    • 1970-01-01
    • 2021-04-17
    • 1970-01-01
    • 2014-12-01
    • 2021-11-04
    • 2010-09-06
    相关资源
    最近更新 更多