pnulls = fltarr(100,50) ;array that will give us the density of absorbing nulls per Mm^3
nnn = findgen(100)/2. ;Height array
nn = findgen(100000)/2000. ;Another height array
denulls = fltarr(100000,50) ;interpolated density ofabsorbing nulls
for j = 0,49 do begin
for i = 3,99 do begin
pnulls[i,j] = pnull(gcz[i,j],q[i])
endfor
denulls[*,j] = interpol(pnulls[*,j], nnn, nn)
endfor
nulls = fltarr(50)
for i = 0,49 do begin
nulls[i]=int_tabulated(zz,denulls[*,i]) ;this integrates the density of nulls over height, producing the density of nulls for a column above the photosphere, for each frequency
endfor
w = findgen(50)/100.+0.01 ;frequency array
ww = findgen(50000)/100000.+0.01 ;another frequency array
wnulls = fltarr(50) ; array giving the normalised kinetic energy spectrum of frequencies
beta = 1.0 ;factor affecting how sharp the tail of the spectrum is, should be between 1 and 2.5
for i = 0,49 do begin
wnulls[i]=nulls[i]*(w[i])^(-beta)
endfor
wnulls = wnulls/(max(wnulls))
wnulls = interpol(wnulls,w,ww,/spline)
end