Structure function revisited

I have revisited the usage of the structure function (SF, as defined in Zhang's thesis [hereafter Zhang (2000)]) to characterize the variability of a source, testing a simplified version of it in conjunction with some simple simulated light curves.

I started from the original Paltani-Pian-Zhang Fortran code and made a simplified IDL version which uses the same structure with a triple loop (on tau, i and j) as a reference.
I then made a faster IDL version which removes one loop exploiting IDL array notation and runs 16 times faster with the same results.

The procedure is called as :
sfu4,t,lc,f,nk,binsize,binwidth,tmin,tmax,tt
where :

The simplification with respect to Zhang's full code is that

Also I do not fully understand the roles of binsize and binwidth which apparently can be assigned arbitrarily irrespectively of the bin size of the light curve.
In practice I tried varying them (one or both) to values different from the actual bin size of the light curve, with confusing results, so all simulations are reported with binsize=binwidth=light curve binsize = 600 s.


I can simulate simple light curves via other simple IDL procedures like :

It is advisable to enlarge your Netscape window in x and y so that you can see better the plots below


Let us first assess the effect of the fluctuations (since we do not compute errors on single points we can't call it the noise) generating uniform light curves with the same plateau=1000 and different level (since the same random seed is used, the fluctuations are exactly the same, just scaled differently).
generated with makecurve,100,1000.,x*sqrt(1000.),600.,t,lc,lastt,a2,sdev Structure function generated with sfu4,t,lc,f,nk,600.,600.,0.,lastt,tt
plot,tt,f/nk
Reference light curve A level=sqrt(1000) Comparison of the 4 structure functions The colours correspond to the 4 light curves plotted on the left.

The SF scales exactly as the square of the fluctuations and is of course flat, pure white noise, so that one can reasonably compute its mean.

plateaulevelfactormean sdev mean SFfactor
1000 15.811 1015.089.617184.5531
1000 31.622 1030.1619.23738.2114
1000 63.254 1060.3138.472952.8416
1000 07.910.5 1007.544.80846.13820.25
double fluctuation amplitude B : level=2*sqrt(1000)
quadruple fluctuation amplitude C : level=4*sqrt(1000)
half fluctuation amplitude D : level=0.5*sqrt(1000)

Let us then assess the effect of intensity generating uniform light curves with different plateau (1000 and 2000) and the same percentage level of fluctuations.
Reference light curve I : plateau=1000. level=30 Comparison of the 2 structure functions The colours correspond to the 2 light curves plotted on the left. The SF scales this time exactly as the square of intensity

plateaulevelfactormean sdev mean SFfactor
1000 30 1 1014.309.123166.0971
2000 60 2 2028.6118.25664.3904
double intensity J : plateau=2000 level=60

Let us then examine the SF of a periodic (sinusoidal) modulation, superimposed on the previous white noise light curve.
Sinusoid of amplitude 100, P=3230 s superimposed on (white) reference light curve A (lc)
makesin,t,cc,100.,3230.
lcs=lc+cc

The SF of a periodically modulated white noise light curve (in red on the right) is also periodic, with equal height maxima at P/2+nP and minima tangent to the SF of the reference light curve
The blue vertical line indicates P/2+nP

The example shown is with a period reasonably large and not multiple of the bin size (600 s), i.e. not so well sampled. Tests have been made with other periods, with similar results (although for periods smaller than 600 s the periodicity is poorly determined, and for submultiples of 600 s (or 600 s itself) there is no contribution to the SF (the modulation adds to fixed bins).

Comparison of the structure functions sfu4,t,lc+cc,fs,nks,600.,600.,0.,lastt,tt

The next case is the addition to the reference white noise light curve of finite wavetrains derived from the same sinusoidal light curves (for a non-integer small number of periods)
generated with u=findgen(i+n)+n
lcs=lc
lcs(u)=lcs(u)+cc(u)
Structure function generated with sfu4,t,lcs,fs,nks,600.,600.,0.,lastt,tt
plot,tt,fs/nks
wavetrain i=57 n=39 (3.2P from 23400 to 33600 s) Comparison of the 4 structure functions with reference

The colours correspond to the 4 light curves plotted on the left, while the white line is the SF of the reference light curve.

The first three light curves (red, magenta and yellow) are all wavetrains ending at the same moment but covering a different number of cycles (the number of peaks in the SF match the number of cycles).
The fourth light curve (cyan) is a wavetrain of the same duration as the third (yellow) one, but located elsewhere.

All SF show similar features, that is, a number of peaks of decreasing amplitude (matching the number of cycles in the wavetrain) and a sharp drop which connects exactly with the reference light curve.

The location of the drop seems somewhat correlated with the duration of the wavetrain.

The cyan SF however is strange since the drop is very close to the end of the tau axis, and there is an intermediate dip in the SF.
However one shall also note that the wavetrain in the cyan light curve is located close to the end of the time interval covered, and anyhow well past its middle. The simplest examples of single flares shown below could help to understand the shape of the fourth SF.

wavetrain i=57 n=15 (7.6P from 9000 to 33600 s)
wavetrain i=57 n=43 (2.4P from 25800 to 33600 s)
wavetrain i=90 n=76 (2.4P from 45600 to 53400 s)

Since sinusoidal variations are quite unlikely in the sources we are interested in, a simpler example superimposes onto the white noise light curve a single flare derived from the positive part of a single sinusoidal cycle (with the period and resolution chosen this is just less than 4 time bins).
generated with flare=[cc(0:2),0.0]
lcs=lc
lcs(i+3)=lcs(i+3)+flare
Structure function generated with sfu4,t,lcs,fs,nks,600.,600.,0.,lastt,tt
plot,tt,fs/nks
Example of 5 single flares

The colours correspond to the 5 SF plotted on the right, while the white line is the reference light curve (or its SF). Each analysed light curve differs from the reference only in the single flare highlighted in colour.

Each flare is located at bin i and lasts 4 bins (the last one is forced to be zero, i.e. equal to the reference curve, instead of the negative value in the original sinusoid, compare IDL expression above), i.e. 1800 s.
The range of time bins interested by the flares is printed above each SF.

The two blue L-shaped lines in the SF plots extend respectively :

  • from the left end of the tau range to the end time of the flare tend=t(i+3)
  • from the right end of the tau range (lastt) back to lastt-tend

The features of the SF can be described as follows :

For the first two light curves (red and magenta), which correspond to flares located in the first half of the light curve, the SF shows a shallow peak, followed by a little drop in correspondence of tend, followed by a plateau or shallow dip, followed by a sharp drop in correspondence of lastt-tend after which one reaches exactly the reference SF.

For a flare located in the second half of the light curve (like the third, in yellow) the location of the two drops (the minor one and the sharp one) is reversed. The terminal, sharpest drop corresoponds to the end time of the flare tend
This is true, although less legible, also for a flare located close to the end of the light curve (like the fifth, in green, SF not shown)

A flare located at the very beginning of the light curve (fourth, in cyan) shows a steady rise without final drop, and also the minor drop at tend is almost invisible.
A flare located at the very end of the light curve (not shown) gives a SF at all similar.

i=20, flare from 12000 to 13800 s
i=40, flare from 24000 to 25800 s
i=75, flare from 45000 to 46800 s
i=0, flare from 0 to 1800 s

Let us examine the effect of the intensity of a single flare.
The figure on the right compares the SF of two light curves with a single flare at the same time (the light curve is not shown).
The red curve is the same reported above for a flare at i=20 (12000-13800 s), while the green curve corresponds to a flare of double intensity occurring at the same time.

One can see that the two SF have approximately the same shape, with the same features at the same positions, but with the one corresponding to the more intense flare being scaled up.


Let us examine the effect of multiple flares.
The figure on the right compares the SF of three light curves, two with a single flare at two different times, and the third (not shown) with both flares present.
The red and magenta curves are the same reported above for a flare at i=20 and i=40, while the green curve corresponds to a light curve with two flares of the same intensity occurring at both i=20 and i=40.

One can see that the green SF stays up until the time of the drop of the red curve (end of first flare) where it has a dip. then stays up again mimicking the magenta curve until the time of the final drop of the magenta curve : at this point it drops until it superimposes exactly to the red curve, which then follows until the red (hidden) and green drop together to the reference white noise value.

The figure on the right compares the SF of four light curves, two with a single flare at two different times, and two with two flares, one of which is of double intensity with respecto to the other. The red and magenta dottedcurves are the same reported above for a flare at i=20 and for a flare at i=40, while the yellow and cyan curves corresponds to a light curve with two flares occurring at both i=20 and i=40 : for the yellow curve the first flare is the most intense, while for the cyan curve it is the second flare to be stronger. (the light curves can be viewed clicking here

One can see that the main features occur as sharp drops at characteristics times, connecting otherwise almost flat portions of the SF curve.


I suspended writing these pages on Oct 15. I restart now on Oct 30 adding some results obtained 2 weeks ago and not yet inserted here.
Let us examine now the behaviour of red noise (or something similar to it), which I generate co-adding several sinusoids of different periods and amplitude proportional to a power of period using the following IDL code :
lcr1=lc
lcc1=lcr1-lcr1+1000
seed=123456789L
for i=P1,P2 do begin
    makesin,t+randomu(seed)*i,cc,float(i),float(i)
    lcr1=lcr1+cc/(P2-P1)
endfor
seed=123456789L
for i=P1,P2 do begin
    makesin,t+randomu(seed)*i,cc,float(i),float(i)
    lcc1=lcc1+cc/(P2-P1)
endfor
use white noise light curve A (lc) as base for the "r" light curve
use a flat light curve at intensity 1000 as base for the "c" light curve
assign a random seed
loop from period P1 to period P2
generate a sinusoid of amplitude proportional to period
and co-add it to the "r" light curve, normalize by number of periods used

reassign same random seed to have same noise
loop from period P1 to period P2
generate same sinusoid as above
and co-add it to the "c" light curve, normalize by number of periods used

Note the trick by which each sinusoid is random shifted by a fraction of its period randomu(seed)*i in order to break coherence.
The "r" light curve is the superposition of white noise with red noise from P1 to P2, while the "c" light curve should be pure red noise.

"r" (red) and "c" (green) light curves for P1=3000, P2=10000 s

On the right their SF compared with the white noise one (in white).
The dashed line is a tau2 power law normalized to the first point.
Since the "c" curves are very similar to the "r" curves, just smoother, we limit ourselves to the pure red noise "c" curves in what follows.

One can note that the SF of the "c" curve is smoother, and closer to the tau2 slope in its beginning.
The other feature is that the first peak is at P1, and ends approximately before P2


"c" red noise light curves for the following period ranges (same start)
  • P1=3000, P2=10000 (green, same as above)
  • P1=3000, P2=20000 (cyan)
  • P1=3000, P2=30000 (blue)


"c" red noise light curves for the following period ranges (another start, but same ends as the previous family)

  • P1=7000, P2=10000 (yellow, compare with green above)
  • P1=7000, P2=20000 (magenta, compare with cyan)
  • P1=7000, P2=30000 (red, compare with blue)

On the right their SF with the respective tau2 (dashed) trends which is confirmed.
However the relationship between the location of the peak and the range P1-P2 where the noise is defined is unclear. The situation is even more complex when P2 is 30000 s, i.e. half of the total interval.

I am however not sure at all that the "red noise" simulations are by any mean realistic (but a periodic function is simpler than shot noise).


I have also tried a different "power index " in the red noise using an amplitude proportional to the square of the period, replacing the statement above for the generation of sinusoidal components with the one below (where the prescaling by 10-4 is used to preserve the same approximate total amplitude of the resulting light curve.
makesin,t+randomu(seed)*i,cc,1.E-4*float(i)^2,float(i)

The red curve below has been generated with such formula in the period range i=P1,P2 with P1=3000 and P2=10000 (same as the reference green curve, which has been generated with the default formula (amplitude proportional to period).

The relevant SFs are reported on the right. They scale in a similar way, with just minor shifts of the peak, and the initial rise in both cases follows a tau2 trend (despite the different power index used in the simulation)



I stop here since it is not clear to me which kind of further simulations can be considered realistic. However I note that the interpretation of the peaks in the SF is either not fully obvious, or strictly related to a specific noise model different from the oversimplistic cases I considered.