1

OSX v10.10.5 和 Gnuplot v5.0

我有一个包含三列数字的数据文件，我读取了存储在其中的值以进行一些计算。但是很费时间！

这是我到目前为止所拥有的：

#user defined function to read data in a file
#see stackoverflow: "Reading dataset value into a gnuplot variable (start of X series)"
at(file, row, col) = system( sprintf("awk -v row=%d -v col=%d 'NR == row {print $col}' %s", row, col, file) )
file="myFile"

do for [k=1:10] { #we read line by line and we want the ratio between column 2/1 and 3/1
f(k) = at(file,k,2)/at(file,k,1)
g(k) = at(file,k,3)/at(file,k,1)

# example of calculation: least square to find the best "i"
do for [i=1:10] {
    f1(i) = (a*i**2 + b*i + c) #function for the least square. a,b,c: floats
    g1(i) = (d*i**2 + e*i + f) #d,e,f: floats
    h(i) = sqrt( (f1(i)-f(k))**2 + (g1(i)-g(k))**2 )
    if (h(i)<hMin) {
        hMin=h(i)
                   }
    else {}
            } #end loop i
print i," ",hMin
} #end loop k

它有效，但正如我所说，这需要时间（每 k 大约 2 分钟）。当我不进行任何计算而只询问 print f(k),g(k) 时，它是 << 1 秒。我怀疑除法可能会导致数字过多和计算效率低下。我使用 round2 函数首先保留 n=4 ：

#see stackoverflow: How to use floor function in gnuplot
round(x) = x - floor(x) < 0.5 ? floor(x) : ceil(x)
round2(x, n) = round(x*10**n)*10.0**(-n)
f(k) = round2((at(file,k,2)/at(file,k,1)),4)
g(k) = round2((at(file,k,3)/at(file,k,1)),4)

但它并没有改变所需的时间。知道发生了什么吗？

4

1 回答 1

0

您没有发布完整的代码（缺少 a、b、...、f 的定义）。但是在您发布的代码部分中，我认为您可以避免awk经常调用。您可以用简单的变量和替换函数f(k)和，因为实际上它们在每个 k 迭代中都是常数。似乎没有必要在每次迭代中重新计算它们。g(k)fkgk

#user defined function to read data in a file
#see stackoverflow: "Reading dataset value into a gnuplot variable (start of X series)"
at(file, row, col) = system( sprintf("awk -v row=%d -v col=%d 'NR == row {print $col}' %s", row, col, file) )
file="myFile"

do for [k=1:10] { #we read line by line and we want the ratio between column 2/1 and 3/1
    at1 = at(file,k,1)
    fk = at(file,k,2)/at1
    gk = at(file,k,3)/at1

    # example of calculation: least square to find the best "i"
    do for [i=1:10] {
        f1i = (a*i**2 + b*i + c) #function for the least square. a,b,c: floats
        g1i = (d*i**2 + e*i + f) #d,e,f: floats
        hi = sqrt( (f1i-fk)**2 + (g1i-gk)**2 )
        if (hi<hMin) {
            hMin=hi   
        } else {
        }
    } #end loop i 
    print i," ",hMin
} #end loop k

但是在缺失的代码中可能有更多有趣的细节会抑制这个解决方案。

于 2017-03-16T20:44:10.163 回答