У меня нет опыта в python, и я конвертирую скрипт NCL в python, в надежде python будет работать намного быстрее. Осматривая, я не нахожу ответа на то, что я считаю простейшими вычислениями в скрипте NCL. Глядя на то, как выполняются более сложные вычисления, я также не нахожу ответа на вопрос, как это можно сделать в python.
. Большая часть вычислений выполняется после преобразования 3-мерных переменных в 1-мерные. переменных, и запрашивая их значения и позиции в пространстве массива. Зная положение переменных t в пространстве массива, мы можем получить значения переменной p, соответствующие целочисленным значениям переменной t.
Вычисления следующие:
- Установить значение в переменная p по умолчанию _FillValue,
- считать число (объем) точек сетки, значение происходит для каждого возможного целого значения в переменной t (сумма во времени и пространстве),
- вычислить индексы времени начала и окончания для каждого из возможных целочисленных значений в переменной t,
- вычислить время продолжительности как разницу (+ 1, потому что числа) между временем окончания и началом в переменной t,
- вычисляет среднюю (пространственно-временную) широту и долготу для каждого из возможных целочисленных значений в переменной t,
- вычисляет площадь (объем / длительность) для каждого возможного целочисленного значения в переменной t ,
- вычисляет среднее значение p из переменной p, где оно соответствует в пространстве-времени каждому возможному целочисленному значению в переменной t, и
- вычисляет процентили p из переменной p, где она соответствует в пространстве-времени каждому возможному целочисленному значению в переменной t.
Все эти вычисления сохраняют значения в 1- размерные массивы, размеры размеров которых равны максимальному целочисленному значению в переменной t. Например, переменная at может иметь целые числа от 0 до 100. Целочисленное значение 0 игнорируется, поэтому каждый из одномерных массивов должен иметь в этом примере 100 значений; (100 томов, 100 начальных времен, 100 конечных времен и т. Д. c.).
Наконец, все одномерные массивы записываются в текстовый файл (с разделителями табуляции), причем каждый столбец является Одномерные массивы.
;===================================================================
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
;===================================================================
begin
;===============================================================
begTime = get_cpu_time()
; Data I/O and data names
; T output, and raw data input
f_t = addfile("t_in.nc","r")
f_p = addfile("p_in.nc","r")
; Data variables
time = f_t->time
lat = f_t->lat
lon = f_t->lon
t_var = f_t->t
p_var = f_p->p
p_fix = p_var
p_fix = where(p_var.eq.9.96921e+36, p_var@_FillValue, p_var)
delete(p_var)
p_var = p_fix
delete(p_fix)
; t = 0, is not measured
; Compute volume, start and end time indices, delta time, centroid lat and centroid lon, area, and percentiles
volume = new(max(t_var)+1, "integer")
start_time = new(max(t_var)+1, "integer")
end_time = new(max(t_var)+1, "integer")
delta_time = new(max(t_var)+1, "integer")
centroid_lat = new(max(t_var)+1, "double")
centroid_lon = new(max(t_var)+1, "double")
area = new(max(t_var)+1, "float", -9999.)
v_av = new(max(t_var)+1, "float", -9999.)
p_10 = new(max(t_var)+1, "float", -9999.)
p_25 = new(max(t_var)+1, "float", -9999.)
p_50 = new(max(t_var)+1, "float", -9999.)
p_75 = new(max(t_var)+1, "float", -9999.)
p_90 = new(max(t_var)+1, "float", -9999.)
t1D = ndtooned(t_var)
p1D = ndtooned(p_var)
dsizes_t = dimsizes(t_var)
do i=1,max(t_var)
indices_t = ind_resolve(ind(t1D.eq.i),dsizes_t)
volume(i) = num(t_var.eq.i)
start_time(i) = indices_t(0,0)
end_time(i) = indices_t(dimsizes(indices_t(:,0))-1,0)
delta_time(i) = 1+end_time(i)-start_time(i)
centroid_lat(i) = avg(lat(indices_t(:,1)))
centroid_lon(i) = avg(lon(indices_t(:,2)))
area(i) = volume(i)/delta_time(i)
v_av(i) = avg(p1D(ind(t1D.eq.i)))
p_10(i) = Percentile(p1D(ind(t1D.eq.i)),10)
p_25(i) = Percentile(p1D(ind(t1D.eq.i)),25)
p_50(i) = Percentile(p1D(ind(t1D.eq.i)),50)
p_75(i) = Percentile(p1D(ind(t1D.eq.i)),75)
p_90(i) = Percentile(p1D(ind(t1D.eq.i)),90)
delete(indices_t)
end do
; Write data as table to text file
r = ispan(1,max(t_var),1)
system("/bin/rm -f var.txt")
fname = "var.txt"
fhead = systemfunc("echo -e tnum $'\t' start $'\t' end $'\t' dt $'\t' c_lat $'\t' c_lon $'\t' vol $'\t' area $'\t' v_avg $'\t' p_10 $'\t' p_25 $'\t' p_50 $'\t' p_75 $'\t' p_90 >> "+fname)
print(fhead)
do i=1,max(t_var)
str_var = sprinti("%8.0i",r(i-1))+"$'\t'"+sprinti("%4.0i",start_time(i))+"$'\t'"+sprinti("%4.0i",end_time(i))+"$'\t'"+sprinti("%4.0i",delta_time(i))+"$'\t'"+\
sprintf("%2.2f",centroid_lat(i))+"$'\t'"+sprintf("%3.2f",centroid_lon(i))+"$'\t'"+\
sprinti("%10.0i",volume(i))+"$'\t'"+sprintf("%8.2f",area(i))+"$'\t'"+sprintf("%3.2f",v_av(i))+"$'\t'"+\
sprintf("%3.2f",p_10(i))+"$'\t'"+sprintf("%3.2f",p_25(i))+"$'\t'"+\
sprintf("%3.2f",p_50(i))+"$'\t'"+sprintf("%3.2f",p_75(i))+"$'\t'"+\
sprintf("%3.2f",p_90(i))
cmd = systemfunc("echo -e " + str_var + " >> "+fname)
print(cmd)
end do
print("Total run time: " + (get_cpu_time() - begTime))
end