;+
; NAME:
;	ark_duty_cycle
; PURPOSE:
;	Makes graph of Ark duty cycle, and timeseries graph for sizes,
;	trace count, and cycle number
; CALLING SEQUENCE:
	PRO ark_duty_cycle, data_file, $
		file		= file			, $
		xysize		= xysize		, $
		include 	= include		, $
		silent		= silent		, $
		time_ago	= time_ago		, $
		start_time	= start_time	, $
		stop_time	= stop_time		, $
		bin_width	= bin_width		, $
		ipv4_probing= ipv4_probing  , $
		ipv6_probing= ipv6_probing  , $
		pfx_probing	= pfx_probing	, $
		bdr_probing	= bdr_probing	, $
		loss_probing= loss_probing	, $
		all_probing	= all_probing	, $
		sizes		= sizes			, $
		cycles		= cycles		, $
		traces		= traces		, $
		charsize	= charsize		, $
		kibibyte	= kibibyte		, $
		unlabeled	= unlabeled		, $
		refresh		= refresh		, $
		stack_sizes	= stack_sizes	, $
		show_gaps	= show_gaps		, $
		cdf			= cdf			, $
		monitor		= monitor		, $
		_extra		= _extra
; INPUT:
;	data_file		fully-qualified name with the list of Ark file
;					names on which the plot are based.
;					Default: ark_duty_cycle.txt
;	file			filename on which names of graphics output
;					files are based.
; OPTIONAL INPUTS:
;	include=include	comma-separated list of monitors to be processed
;					default: all the monitors in the data file
;	start_time=start_time
;					start time for graphs; default: earliest file
;	stop_time=stop_time
;					stop time of graphs; default: latest file
;	time_ago=time_ago
;					instead of start_time, specify a time range:
;					start_time then becomes stop_time-time_ago
;
;	/sizes			plots growth of IPv4 and IPv4 datasets as
;					function of time. By default, decimal units are used
;					(i.e. k=1000 bytes, etc.)
;	/kibibyte		if set then binary units are used
;					(i.e. ki=1024 bytes, etc.)
;
;	/cycles			plots graph of cycle numbers as function of time
;
;	/unlabeled		omit labeling
;	/show_gaps		used for duty cycle plot to emphasize times were
;					no data are present: these time ranges are explicitly
;					erased with a final polyfill. This can have the
;					unwelcome side effect that very small gaps (less
;					then one pixel) are shown as a full pixel wide.
;	charsize=charsize
;					character size for labeling
; INCLUDE:
	@compile_opt.pro
; CALLS:
;	InitVar, TimeUnit, GetFileSpec, txt_read
;	TimeGet, TimeSet, strposn, IsType, IsTime
;	set_page, get_page, plotcurves
; MODIFICATION HISTORY:
;	MAR-2017, started adding documentation
;-

InitVar, bin_width		, 1
InitVar, silent	 		, 0
InitVar, charsize		, 1
InitVar, kibibyte		, /key
InitVar, unlabeled		, /key
InitVar, refresh		, /key
InitVar, sizes   		, /key
InitVar, ipv4_probing	, /key
InitVar, ipv6_probing	, /key
InitVar, pfx_probing	, /key
InitVar, bdr_probing    , /key
InitVar, loss_probing   , /key
InitVar, all_probing	, /key
InitVar, cycles			, /key
InitVar, traces			, /key
InitVar, stack_sizes	, /key
InitVar, show_gaps		, /key
InitVar, cdf			, /key


format = 'YYYYMNDD'
uday   = TimeUnit(/day)
usec   = TimeUnit(/second)
sec_per_day = 86400LL

InitVar, data_file, who_am_i(/procedure)+'.txt'
save_file = GetFileSpec(data_file,upto='name')
IF GetFileSpec(data_file,part='type') EQ '.gz' THEN save_file = GetFileSpec(save_file,upto='name')
save_file += '.idl.sav'

; Check for the save file. If not present, then make sure
; refresh=0, forcing (re)reading of list of files.

IF NOT refresh THEN BEGIN
	IF (file_search(save_file))[0] EQ '' THEN BEGIN
		IF silent LT 1 THEN message, /info, "IDL save file not found, '"+save_file+"'"
		refresh = 1
	ENDIF
ENDIF

IF refresh THEN BEGIN

	IF silent LT 1 THEN message, /info, data_file
	IF NOT txt_read(data_file,file_name,silent=silent) THEN RETURN

	nfile = n_elements(file_name)
	compressed_size   = lon64arr(nfile)
	uncompressed_size = lon64arr(nfile)
	trace_count		  = lon64arr(nfile)

	FOR i=0L,nfile-1 DO BEGIN					; Loop over all Ark filenames
		tmp = strtok(file_name[i],/extract)
		file_name        [i] = tmp[0]
		compressed_size  [i] = tmp[1]
		uncompressed_size[i] = tmp[2]
		IF n_elements(tmp) GT 3 THEN trace_count[i] = tmp[3]
	ENDFOR

	; Team-probing:		daily.l7.t1.c000027.20070913.amw-us.warts.gz
	;					aal-dk.team-probing.c006951.20180916.warts.gz
	; Topo-v6:			topo-v6.l8.20170311.1489241615.iad-us.warts.gz
	; Prefix-probing:	ams2-nl.20151209.1449626400.warts.gz
    ; Bordermap:        cmn2-ma.1468113863.bdrmap.warts.gz
    ; Loss data:        acc-gh.loss.20160930.1475256817.warts.gz

	file_name = GetFileSpec(file_name,part='name')
	save, nfile, file_name, compressed_size, uncompressed_size, trace_count, filename=save_file

ENDIF ELSE BEGIN

	IF silent LT 1 THEN message, /info, save_file
	restore, save_file

ENDELSE

IF IsType(monitor,/defined) THEN BEGIN
	i = where( strpos( filename, '.'+monitor+'.' ) NE -1 OR strpos( filename, monitor+'.' EQ 0 ), nfile )
	IF nfile EQ 0 THEN message, "no files for monitor '"+monitor+"'"
	file_name			= file_name			[i]
	compressed_size		= compressed_size	[i]
	uncompressed_size	= uncompressed_size	[i]
	trace_count			= trace_count		[i]
ENDIF

times	 = lon64arr(nfile)
names	 = strarr  (nfile)
activity = intarr  (nfile)-1  ; 0: team-probing, 1: ipv6, 2: prefix; 3: bordermap; 4: loss; 5: all

IF cycles THEN cycle_id = strarr(nfile)

ACTIVITY_TEAM_PROBING	= 0
ACTIVITY_TOPO_V6		= 1
ACTIVITY_PREFIX_PROBING = 2
ACTIVITY_BDR			= 3
ACTIVITY_LOSS			= 4
ACTIVITY_ALL			= 5

ii_file = where( strpos(file_name,'daily.'  ) NE -1 OR strpos(file_name,'.team-probing.') NE -1, n4_file )

IF n4_file NE 0 THEN BEGIN	; daily.l7.t1.c000027.20070913.amw-us.warts or aal-dk.team-probing.c006951.20180916.warts
	f = file_name[ii_file]
	i = where( strpos(f,'daily.') NE -1, complement=j )
	IF i[0] NE -1 THEN  BEGIN
		g = f[i]
		n = where( strlen(g) EQ 41 )
		IF n[0] NE -1 THEN g[n] = strmid(g[n],0,35)+' '+strmid(g[n],35)

		i = ii_file[i]
		times   [i] = TimeGet(TimeSet(strmid(g,20,8),format=format),/unix)
		names   [i] = strcompress(strmid(g,29,7),/rem)
		activity[i] = ACTIVITY_TEAM_PROBING
		IF cycles THEN cycle_id[i] = strmid(g,13,6)
	ENDIF

	IF j[0] NE -1 THEN  BEGIN
		g = f[j]
		n = where( strlen(g) EQ 42 )
		IF n[0] NE -1 THEN g[n] = strmid(g[n],0,6)+' '+strmid(g[n],6)

		j = ii_file[j]
		times   [j] = TimeGet(TimeSet(strmid(g,29,8),format=format),/unix)
		names   [j] = strcompress(strmid(g,0,7),/rem)
		activity[j] = ACTIVITY_TEAM_PROBING
		IF cycles THEN cycle_id[j] = strmid(g,22,6)
	ENDIF

ENDIF

ii_file = where( strpos(file_name,'topo-v6.') NE -1, n6_file )

IF n6_file NE 0 THEN BEGIN	; topo-v6.l8.20170311.1489241615.iad-us.warts
	f = file_name[ii_file]
	n = where( strlen(f) EQ 43)
	IF n[0] NE -1 THEN f[n] = strmid(f[n],0,37)+' '+strmid(f[n],37)
	times   [ii_file] = long64(strmid(f,20,10))
	names   [ii_file] = strcompress(strmid(f,31,7),/rem)
	activity[ii_file] = ACTIVITY_TOPO_V6
ENDIF

ii_file = where( strpos(file_name,'daily.'          ) EQ -1 AND   $
                 strpos(file_name,'.team-probing.'  ) EQ -1 AND   $
                 strpos(file_name,'topo-v6.'        ) EQ -1 AND   $
                 strpos(file_name,'.bdrmap.'        ) EQ -1 AND   $
                 strpos(file_name,'.loss.'          ) EQ -1, np_file )

IF np_file NE 0 THEN BEGIN	; ams-nl.20151209.1449626400.warts
	f = file_name[ii_file]
	n = where( strlen(f) EQ 32 )
	IF n[0] NE -1 THEN f[n] = strmid(f[n],0,6)+' '+strmid(f[n],6)
	times   [ii_file] = long64(strmid(f,17,10))
	names   [ii_file] = strcompress(strmid(f,0,7),/rem)
	activity[ii_file] = ACTIVITY_PREFIX_PROBING
ENDIF

ii_file = where( strpos(file_name,'.bdrmap.') NE -1, nb_file )

IF nb_file NE 0 THEN BEGIN	; cmn2-ma.1468113863.bdrmap.warts
	f = file_name[ii_file]  ; OW04A151A3102E.1465600279.bdrmap.warts
	n = where( strlen(f) EQ 38 )
	IF n[0] NE -1 THEN BEGIN
		times[ii_file[n]] = long64(strmid(f[n],15,10))
		names[ii_file[n]] = strcompress(strmid(f[n],0,14),/rem)
	ENDIF
	n = where( strlen(f) EQ 30)
	IF n[0] NE -1 THEN BEGIN
		times[ii_file[n]] = long64(strmid(f[n],7,10))
		names[ii_file[n]] = strcompress(strmid(f[n],0,6),/rem)
	ENDIF
	n = where( strlen(f) EQ 31)
	IF n[0] NE -1 THEN BEGIN
		times[ii_file[n]] = long64(strmid(f[n],8,10))
		names[ii_file[n]] = strcompress(strmid(f[n],0,7),/rem)
	ENDIF

	n = where ( strlen(f) NE 38 AND strlen(f) NE 30 AND strlen(f) NE 31 )
	IF n[0] NE -1 THEN message, 'oops nb_file'

	activity[ii_file] = ACTIVITY_BDR
ENDIF

ii_file = where( strpos(file_name,'.loss.'  ) NE -1, nl_file )

IF nl_file NE 0 THEN BEGIN	; acc-gh.loss.20160930.1475256817.warts
	f = file_name[ii_file]
	n = where( strlen(f) EQ 37)
	IF n[0] NE -1 THEN f[n] = strmid(f[n],0,6)+' '+strmid(f[n],6)
	times   [ii_file] = long64(strmid(f,22,10))
	names   [ii_file] = strcompress(strmid(f,0,7),/rem)
	activity[ii_file] = ACTIVITY_LOSS
ENDIF

destroyvar, ii_file

i = where( activity EQ -1 )
IF i[0] NE -1 THEN BEGIN
	message, /info, 'files without assigned activities encountered'
	RETURN
ENDIF

IF ipv4_probing+ipv6_probing+pfx_probing+bdr_probing+loss_probing+all_probing EQ 0 THEN BEGIN
	IF n4_file+n6_file+nb_file+nl_file+np_file EQ 0 THEN message, 'no data files available'
	IF sizes+cycles+traces GT 0 THEN message, 'set /ipv4_probing, /ipv6_probing, /pfx_probing, /bdr_probing, /loss_probing and/or /all_probing'
ENDIF ELSE BEGIN
	IF cycles THEN BEGIN
		ipv4_probing = 1
		ipv6_probing = 0
		pfx_probing  = 0
		all_probing  = 0
	ENDIF ELSE BEGIN
		IF sizes+traces EQ 0 THEN message, 'set one of /sizes or /traces'
	ENDELSE
ENDELSE

IF ipv4_probing AND n4_file EQ 0 THEN BEGIN
	message, /info, 'no team-probing files available .. ignore'
	IF n6_file+np_file+nb_file+nl_file EQ 0 THEN message, 'no data files left to plot'
	ipv4_probing = 0
ENDIF
IF ipv6_probing AND n6_file EQ 0 THEN BEGIN
	message, /info, 'no IPv6 files available .. ignore'
	IF n4_file+np_file+nb_file+nl_file EQ 0 THEN message, 'no data files left to plot'
	ipv6_probing = 0
ENDIF
IF pfx_probing AND np_file EQ 0 THEN BEGIN
	message, /info, 'no prefix-probing files available .. ignore'
	IF n4_file+n6_file+nb_file+nl_file EQ 0 THEN message, 'no data files left to plot'
	pfx_probing = 0
ENDIF
IF bdr_probing AND nb_file EQ 0 THEN BEGIN
	message, /info, 'no bdr-probing files available .. ignore'
	IF n4_file+n6_file+np_file+nl_file EQ 0 THEN message, 'no data files left to plot'
	bdr_probing = 0
ENDIF
IF loss_probing AND nl_file EQ 0 THEN BEGIN
	message, /info, 'no loss-probing files available .. ignore'
	IF n4_file+n6_file+np_file+nb_file EQ 0 THEN message, 'no data files left to plot'
	loss_probing = 0
ENDIF
IF all_probing THEN BEGIN
	IF n_elements(file_name) EQ 0 THEN message, 'no data files available'
ENDIF

destroyvar, n4_file, n6_file, np_file, nb_file, nl_file

InitVar, include, strjoin(names[ uniq(names,sort(names)) ],',')

ark_name = strtrim(strtok(include,',',/extract),2)
nark = n_elements(ark_name)

min_time  = lonarr(nark)		; Earliest file
max_time  = lonarr(nark)		; Latest file
ark_cycle = lonarr(nark)		; Cycle

cnt_ipv4  = lonarr(nark)		; Nr of IPv4 files
cnt_ipv6  = lonarr(nark)		; Nr of IPv6 files
cnt_pfx   = lonarr(nark)		; Nr of prefix probing files
cnt_bdr   = lonarr(nark)		; Nr of border map files
cnt_loss  = lonarr(nark)		; Nr of loss files
cnt_all   = lonarr(nark)		; Nr of all files

min_ipv4  = lonarr(nark)		; Earliest IPv4 file
min_ipv6  = lonarr(nark)		; Earliest IPv6 file
min_pfx   = lonarr(nark)		; Earliest prefix probing file
min_bdr   = lonarr(nark)		; Earliest border map file
min_loss  = lonarr(nark)		; Earliest loss file
min_all   = lonarr(nark)		; Earliest file

in_include  = bytarr(nfile)
latest_time = max(times)

FOR i=0,nark-1 DO BEGIN			; All monitors in include list
	name = ark_name[i]

	n = where(name EQ names)
	IF n[0] EQ -1 THEN BEGIN

		IF silent LT 1 THEN message, /info, "no files for monitor '"+name+"'"
		min_time[i] = latest_time
		max_time[i] = latest_time

	ENDIF ELSE BEGIN

		in_include[n] = 1B

		; IPv4 files cover one day,  IPv6 files cover two days.
		; If the latest file is an IPv4 file add one day to max_time;
		; if it is an IPv6 file then add two days.

		min_time[i] = min( times[n]   )
		max_time[i] = max( times[n], m) + sec_per_day*( activity[n[m]] EQ ACTIVITY_TOPO_V6 ? 2 : 1 )

		; IPv4 files have a team number; for IPv6 files the team nr is null string.

		i4 = where( activity[n] EQ ACTIVITY_TEAM_PROBING  , n4)
		i6 = where( activity[n] EQ ACTIVITY_TOPO_V6       , n6)
		ip = where( activity[n] EQ ACTIVITY_PREFIX_PROBING, np)
		ib = where( activity[n] EQ ACTIVITY_BDR           , nb)
		il = where( activity[n] EQ ACTIVITY_LOSS          , nl)
		na = n4+n6+np+nb+nl

		cnt_ipv4[i] = n4		; # IPv4 files
		cnt_ipv6[i] = n6		; # IPv6 files
		cnt_pfx [i] = np		; # Prefix-probing files
		cnt_bdr [i] = nb		; # Bordermap files
		cnt_loss[i] = nl		; # Loss files
		cnt_all [i] = n_elements(n)

		IF n4 NE 0 THEN min_ipv4[i] = min( times[n[i4]] )
		IF n6 NE 0 THEN min_ipv6[i] = min( times[n[i6]] )
		IF np NE 0 THEN min_pfx [i] = min( times[n[ip]] )
		IF nb NE 0 THEN min_bdr [i] = min( times[n[ib]] )
		IF nl NE 0 THEN min_loss[i] = min( times[n[il]] )
		min_all [i] = min( times[n    ] )

	ENDELSE

ENDFOR

; For some monitors no files may be available.
; Typically, these are new monitors that have just been activated,
; but for which no data are available yet.

i = where(in_include, nfile)
IF nfile EQ 0 THEN BEGIN
	message, /info, 'no files for this selection of monitors'
	RETURN
ENDIF

; At this point we have selected all the monitors we want to process
; (in ark_name). Throw away all data for monitors not selected.

file_name			= file_name			[i]
times 				= times				[i]
names 				= names				[i]
activity			= activity			[i]
compressed_size		= compressed_size	[i]
uncompressed_size	= uncompressed_size	[i]
trace_count			= trace_count		[i]

IF cycles THEN cycle_id = cycle_id[i]

i = sort(min_time)		; Sort monitors on time of 1st datafile.

ark_name  = ark_name [i]
ark_cycle = ark_cycle[i]
min_time  = min_time [i]
max_time  = max_time [i]

cnt_ipv4  = cnt_ipv4 [i]
cnt_ipv6  = cnt_ipv6 [i]
cnt_pfx   = cnt_pfx  [i]
cnt_bdr   = cnt_bdr  [i]
cnt_loss  = cnt_loss [i]
cnt_all   = cnt_all  [i]

min_ipv4  = min_ipv4 [i]
min_ipv6  = min_ipv6 [i]
min_pfx   = min_pfx  [i]
min_bdr   = min_bdr  [i]
min_loss  = min_loss [i]
min_all   = min_all  [i]

first_time = min(times)
last_time  = max(times,m) + sec_per_day*(activity[m] EQ ACTIVITY_TOPO_V6 ? 2 : 1)

n = where(cnt_all NE 0)

IF IsType(time_ago,/defined) THEN BEGIN
	InitVar, stop_time, max( max_time[n] )
	IF IsType(stop_time,/string) THEN stop_time = TimeSet(stop_time)
	IF IsTime(stop_time) THEN stop_time = TimeGet(stop_time,/unix,/scalar)
	IF IsTime(time_ago) THEN time_ago = TimeGet(time_ago,/unix,/scalar) ELSE time_ago *= sec_per_day
	start_time = stop_time-time_ago
ENDIF ELSE IF IsType(start_time,/defined) THEN BEGIN
	IF IsType(start_time,/string) THEN start_time = TimeSet(start_time)
	IF IsTime(start_time) THEN start_time = TimeGet(start_time,/unix,/scalar)
	InitVar, stop_time , max( max_time[n] )
	IF IsType(stop_time,/string) THEN stop_time = TimeSet(stop_time)
	IF IsTime(stop_time) THEN stop_time = TimeGet(stop_time,/unix,/scalar)
ENDIF ELSE IF IsType(stop_time,/defined) THEN BEGIN
	IF IsType(stop_time,/string) THEN stop_time = TimeSet(stop_time)
	IF IsTime(stop_time) THEN stop_time = TimeGet(stop_time,/unix,/scalar)
	InitVar, start_time, min( min_time[n] )
	IF IsType(start_time,/string) THEN start_time = TimeSet(start_time)
	IF IsTime(start_time) THEN start_time = TimeGet(start_time,/unix,/scalar)
ENDIF ELSE BEGIN
	start_time = first_time
	stop_time  = last_time
ENDELSE

trange = [start_time,stop_time]+0.025*(stop_time-start_time)*[-1,1]

IF silent LT 2 THEN BEGIN
	IF silent LT 1 THEN BEGIN
		tmin = TimeGet(/ymd,TimeSet(unix=min_time),upto=uday)
		tmax = TimeGet(/ymd,TimeSet(unix=max_time),upto=uday)
		FOR i=0,nark-1 DO BEGIN
			IF ark_name[i] EQ '' THEN continue
			message, /info, string(format='(A8,2(1X,A),6(A,I5))', ark_name[i], tmin[i], tmax[i], $
				' IPv4 ' , cnt_ipv4[i],	$
				' IPv6 ' , cnt_ipv6[i], $
				' Pfx  ' , cnt_pfx [i], $
                ' Bdr  ' , cnt_bdr [i], $
                ' Loss ' , cnt_loss[i], $
				' All  ' , cnt_all [i])
		ENDFOR
		destroyvar, tmin, tmax
	ENDIF
	message,/info,'Data time range '+strjoin(TimeGet(/ymd,TimeSet(unix=[first_time,last_time]), upto=usec),' to ')
	message,/info,'Axis time range '+strjoin(TimeGet(/ymd,TimeSet(unix=trange                ), upto=usec),' to ')
	message,/info,'Plot time range '+strjoin(TimeGet(/ymd,TimeSet(unix=[start_time,stop_time]), upto=usec),' to ')
ENDIF

ctable = 12
color = [14,8,15]*16
InitVar, file, 'ark_duty_cycle.eps'
is_eps = GetFileSpec(file,part='type') EQ '.eps'

title = TimeGet(TimeSet(unix=start_time),format='YYYY Mon DD')+' to '+TimeGet(TimeSet(unix=stop_time),format='YYYY Mon DD')

igraph = [ipv4_probing,ipv6_probing,pfx_probing,bdr_probing,loss_probing,all_probing]
nprobe = n_elements(igraph)
ngraph = fix(total(igraph))

IF cycles EQ 1 THEN BEGIN				; Graph of cycle numbers

	ndays = 1LL+(stop_time-start_time)/(bin_width*sec_per_day)
	tunix = start_time+bin_width*sec_per_day*l64indgen(ndays)
	tplot = TimeSet(unix=tunix)

	current_cycle_id = lon64arr(ndays)
	FOR i=0L,nfile-1 DO BEGIN
		tmp = min(abs(tunix-times[i]),n)
		IF tmp LE 43200*bin_width THEN BEGIN
			cycle = cycle_id[i]
			IF cycle NE '' THEN current_cycle_id[n] >= long64(cycle)
		ENDIF
	ENDFOR

	i = where( current_cycle_id NE 0 )
	IF i[0] EQ -1 THEN BEGIN
		message, /info, 'no cycles for this selection of monitors'
		RETURN
	ENDIF

	tunix = tunix[i]
	tplot = tplot[i]
	current_cycle_id = current_cycle_id[i]

	InitVar, xysize, [8,6]*80, set=xysize_
	xysize_ *= ngraph EQ 1 ? [1,1] : ngraph GT 2 ? [2,2] : stack_sizes ? [1,2] : [2,1]
	xysize_ = long( xysize_*(is_eps ? 1.0/100.0*3 : 1) )
	scale_factor = is_eps ? 10 : 1

	set_page, $
		_extra			= _extra	, $
		file			= file		, $
		xysize			= xysize_	, $
		ctable 		 	= ctable	, $
		/new						, $
		silent			= silent	, $
		scale_factor	= scale_factor

	plotcurves, tplot, current_cycle_id, $
		;xrange	= trange		, $
		/noxtitle				, $
		silent	= silent		, $
		charsize= charsize		, $
		ymargin = [3, 5]		, $
		ytitle	= "IPv4 Cycle number", $
		/exact					, $
		_extra	= _extra

	IF NOT unlabeled THEN BEGIN
		x0 = !x.crange
		y0 = !y.crange
		dyval = 0.05*(y0[1]-y0[0])
		xyouts, mean(!x.crange), !y.crange[1]+2.0*dyval, align=0.5, 'IPv4 Ark Cycles '+title, charsize=charsize, _extra=_extra
	ENDIF

	get_page, $
		silent	= silent	, $
		_extra	= _extra

ENDIF ELSE IF traces THEN BEGIN

	ndays = 1LL+(stop_time-start_time)/(bin_width*sec_per_day)
	tunix = start_time+bin_width*sec_per_day*l64indgen(ndays)
	tplot = TimeSet(unix=tunix)

	total_trace_count = lon64arr(ndays,nprobe)

	FOR i=0L,nfile-1 DO BEGIN
		tmp = min(abs(tunix-times[i]),n)
		IF tmp LE 43200*bin_width THEN BEGIN
			j = activity[i]
			total_trace_count[n,j] += trace_count [i]
			total_trace_count[n,5] += trace_count [i]
		ENDIF
	ENDFOR

	InitVar, xysize, [8,6]*80, set=xysize_
	xysize_ *= ngraph EQ 1 ? [1,1] : ngraph GT 2 ? [2,2] : stack_sizes ? [1,2] : [2,1]
	xysize_ = long( xysize_*(is_eps ? 1.0/100.0*3 : 1) )
	scale_factor = is_eps ? 10 : 1

	set_page, $
		_extra			= _extra	, $
		file			= file		, $
		xysize			= xysize_	, $
		ctable  		= ctable	, $
		/new						, $
		silent			= silent	, $
		scale_factor	= scale_factor

	IF ngraph GT 1 THEN BEGIN
		p_multi = !p.multi
		!p.multi = ngraph GT 4 ? [0,2,3] : ngraph GT 2 ? [0,2,2] : stack_sizes ? [0,1,2] : [0,2,1]
	ENDIF

	FOR i=0,nprobe-1 DO BEGIN
		IF igraph[i] EQ 0 THEN continue

		y = cdf ? total( total_trace_count[*,i] ) : max( total_trace_count[*,i] )
		unit = floor(alog10(y))

		kunit = unit/3
		runit = unit-3*kunit
		post = (['','10 ','100 '])[runit]+(['','k','M','G','T','P'])[kunit]

		unit = 10.0^unit
		yrange = [-0.05,1]*ceil( y/unit )

		ytitle = (['IPv4', 'IPv6', 'IPv4','Ark probing'])[i]+', '+post+(cdf ? '' : ' per day')

		TimeXAxis, TimeSet(unix=trange), $
			/noxtitle			, $
			silent	= silent	, $
			charsize= charsize	, $
			ymargin = [3, 5]	, $
			/exact				, $
			_extra	= _extra

		plotbars, tplot, total_trace_count[*,i]/unit , $
			cdf		= cdf		, $
			/oplot				, $
			/newyaxis			, $
			yrange	= yrange	, $
			low_margin= 0		, $
			ystyle	= 1			, $
			ytitle  = ytitle	, $
			silent	= silent	, $
			shade	= color[0]	, $
			/noline				, $
			charsize= charsize	, $
			_extra	= _extra

		cnt_ip = i EQ 0 ? cnt_ipv4 : i EQ 1 ? cnt_ipv6 : i EQ 2 ? cnt_pfx : cnt_all
		mons = where( ark_name NE '' AND cnt_ip NE 0, nmon )

		y0 = !y.crange
		dyval = 0.080*(y0[1]-y0[0])

		IF NOT unlabeled THEN BEGIN
			s = (['Ark IPv4 team probing','Ark IPv6 probing','Ark IPv4 prefix probing','Ark bordermap probing','Ark loss probing','Ark probing'])[i]
			xyouts, mean(!x.crange), !y.crange[1]+2.0*dyval, align=0.5, s+' Traces '+title, charsize=charsize, _extra=_extra
			xyouts, mean(!x.crange), !y.crange[1]+1.2*dyval, align=0.5, strcompress(nmon,/rem)+' contributing monitors', charsize=charsize, _extra=_extra
		ENDIF

	ENDFOR

	get_page, $
		silent	= silent	, $
		_extra	= _extra

	IF ngraph GT 1 THEN !p.multi = p_multi

ENDIF ELSE IF sizes THEN BEGIN		; Dataset size graph for IPv4 and IPv6

	ndays = 1LL+(stop_time-start_time)/(bin_width*sec_per_day)
	tunix = start_time+bin_width*sec_per_day*l64indgen(ndays)
	tplot = TimeSet(unix=tunix)

	total_compressed_size   = lon64arr(ndays,nprobe)
	total_uncompressed_size = lon64arr(ndays,nprobe)

	FOR i=0L,nfile-1 DO BEGIN
		tmp = min(abs(tunix-times[i]),n)
		IF tmp LE 43200*bin_width THEN BEGIN
			j = activity[i]
			total_compressed_size  [n,j] += compressed_size  [i]
			total_uncompressed_size[n,j] += uncompressed_size[i]
			total_compressed_size  [n,5] += compressed_size  [i]
			total_uncompressed_size[n,5] += uncompressed_size[i]
		ENDIF
	ENDFOR

	InitVar, xysize, [8,6]*80, set=xysize_
	xysize_ *= ngraph EQ 1 ? [1,1] : ngraph GT 2 ? [2,2] : stack_sizes ? [1,2] : [2,1]
	xysize_ = long( xysize_*(is_eps ? 1.0/100.0*3 : 1) )
	scale_factor = is_eps ? 10 : 1

	set_page, $
		_extra			= _extra	, $
		file			= file		, $
		xysize			= xysize_	, $
		ctable  		= ctable	, $
		/new						, $
		silent			= silent	, $
		scale_factor	= scale_factor

	IF ngraph GT 1 THEN BEGIN
		p_multi = !p.multi
		!p.multi = ngraph GT 2 ? [0,2,2] : stack_sizes ? [0,1,2] : [0,2,1]
	ENDIF

	FOR i=0,nprobe-1 DO BEGIN
		IF igraph[i] EQ 0 THEN continue

		y = cdf ? total( total_uncompressed_size[*,i] ) : max( total_uncompressed_size[*,i] )
		unit = floor(alog10(y))

		kunit = unit/3
		runit = unit-3*kunit
		post = (['','10 ','100 '])[runit]+(['','k','M','G','T','P','E','Z','Y'])[kunit]+'B'

		unit = 10.0^unit
		yrange = [-0.05,1]*ceil( y/unit )

		ytitle = (['IPv4', 'IPv6', 'IPv4', 'Bordermap','Loss','Ark probing'])[i]+', '+post+(cdf ? '' : ' per day')

		TimeXAxis, TimeSet(unix=trange), $
			/noxtitle			, $
			silent	= silent	, $
			charsize= charsize	, $
			ymargin = [3, 5]	, $
			/exact				, $
			_extra	= _extra

		plotbars, tplot, total_uncompressed_size[*,i]/unit , $
			cdf		= cdf		, $
			/oplot				, $
			/newyaxis			, $
			yrange	= yrange	, $
			low_margin= 0		, $
			ystyle	= 1			, $
			ytitle  = ytitle	, $
			silent	= silent	, $
			shade	= color[0]	, $
			/noline				, $
			charsize= charsize	, $
			_extra	= _extra

		plotbars, tplot, total_compressed_size[*,i]/unit, $
			cdf		= cdf		, $
			/oplot				, $
			low_margin= 0		, $
			silent	= silent	, $
			shade	= color[2]	, $
			/noline				, $
			_extra	= _extra

		x0 = !x.crange
		y0 = !y.crange
		dxval = 0.075*(x0[1]-x0[0])
		dyval = 0.080*(y0[1]-y0[0])

		cnt_ip = i EQ 0 ? cnt_ipv4 : i EQ 1 ? cnt_ipv6 : i EQ 2 ? cnt_pfx : cnt_all
		mons = where( ark_name NE '' AND cnt_ip NE 0, nmon )

		IF NOT unlabeled THEN BEGIN
			x = x0[0]+dxval+[0,1,1,0]*dxval
			y = y0[1]+(-0.10+0.02*[0,0,1,1])*(y0[1]-y0[0])
			polyfill, x, y, color=color[0]
			xyouts, x[1]+15,y[1], 'uncompressed', charsize=charsize, _extra=_extra

			y -= 0.05*(y0[1]-y0[0])
			polyfill, x, y, color=color[2]
			xyouts, x[1]+15,y[1], 'compressed'  , charsize=charsize, _extra=_extra

			s = (['Ark IPv4 team probing','Ark IPv6 probing','Ark IPv4 prefix probing','Ark bordermap probing','Ark loss probing','Ark probing'])[i]
			xyouts, mean(!x.crange), !y.crange[1]+2.0*dyval, align=0.5, s+' Sizes '+title, charsize=charsize, _extra=_extra
			xyouts, mean(!x.crange), !y.crange[1]+1.2*dyval, align=0.5, strcompress(nmon,/rem)+' contributing monitors', charsize=charsize, _extra=_extra
		ENDIF

		; Plot # contributing monitors along top of graph

		IF nmon GT 0 THEN BEGIN
			min_ip = i EQ 0 ? min_ipv4 : i EQ 1 ? min_ipv6 : i EQ 2 ? min_pfx : min_all
			on_time = min_ip[mons]
			on_time = on_time[ sort(on_time) ]
			on_tag  = 1+indgen(n_elements(on_time))

			n = where( start_time LE on_time AND on_time LE stop_time )
			IF n[0] NE -1 THEN BEGIN
				on_tag  = on_tag [n]
				on_time = on_time[n]

				n = where ( on_tag/25*25 EQ on_tag, cnt )
				IF cnt LT 3 THEN n = where ( on_tag/ 5* 5 EQ on_tag, cnt )
				on_tag  = on_tag [n]
				on_time = on_time[n]

				on_time = TimeScale( TimeSet(unix=on_time) ) 

				x1 = (convert_coord(!y.crange,/to_normal))[0]
				dx = 0.05
				FOR j=0,n_elements(on_tag)-1 DO BEGIN
					x2 = (convert_coord( [on_time[j],!y.crange[1]], /to_normal ))[0]
					IF x2-x1 GT dx THEN BEGIN
						plots , on_time[j]*[1,1],!y.crange[1]+0.25*[0,dyval]
						xyouts, on_time[j], !y.crange[1]+0.5*dyval, align=0.5, strcompress(on_tag[j],/rem), charsize=charsize, _extra=_extra
						x1 = x2
					ENDIF
				ENDFOR
			ENDIF
		ENDIF

	ENDFOR

	get_page, $
		silent	= silent	, $
		_extra	= _extra

	IF ngraph GT 1 THEN !p.multi = p_multi

ENDIF ELSE BEGIN				; Duty-cycle graph

	InitVar, xysize, [8,10]*80, set=xysize_
	xysize_ = long( xysize_*(is_eps ? 1.0/100.0*3 : 1) )
	scale_factor = is_eps ? 10 : 1

	set_page, $
		_extra			= _extra	, $
		file			= file		, $
		xysize			= xysize_	, $
		ctable  		= ctable	, $
		silent			= silent	, $
		/new						, $
		scale_factor	= scale_factor

	; Plot x-axis (time axis) and left y-axis (monitor count)

	yrange = [0,nark]+[-1,1]*3
	IF IsType(time_ago,/undefined) AND IsType(start_time,/undefined) AND IsType(stop_time,/undefined) THEN from = TimeUnit(/year)

	IF unlabeled THEN BEGIN
		nmajor = 3
		xmargin = [3,1]
		ymargin = [2,1]
	ENDIF ELSE BEGIN
		xmargin = [8,18]
		ymargin = [3, 6]
		ytitle	= 'nr of active monitors'
	ENDELSE

	; Plot the axes (no data yet)

	plotcurves, TimeSet(unix=trange), yrange	, $
		/exact					, $
		nmajor	= nmajor		, $
		yrange	= yrange		, $
		ystyle	= 1				, $
		ytitle	= ytitle		, $
		/nodata					, $
		xmargin = xmargin		, $
		ymargin = ymargin		, $
		/noxtitle				, $
		silent	= silent		, $
		charsize= charsize		, $
		from	= from			, $
		_extra	= _extra


	; 0: team, 1: ipv6, 2: prefix; 3: bordermap; 4: loss; 5: all

	plot_activity = (where(activity EQ ACTIVITY_PREFIX_PROBING	))[0] NE -1 ? 2 : (where(activity EQ ACTIVITY_BDR))[0] NE -1 ? 3 : 0
	activity_title = (['','IPv6 ','Prefix Probing ', 'Bordermap ', 'Loss '])[plot_activity]

	; Plot title above graph

	IF NOT unlabeled THEN xyouts, /normal, 0.5,0.970, align=0.5, 'Ark '+activity_title+'Coverage, '+title, charsize=charsize, _extra=_extra

	nn = [0,0,0]

	FOR i=0,nark-1 DO BEGIN			; Plot duty cycle for all monitors

		IF ark_name[i] EQ '' THEN continue	; Skip if monitor name not given

		; IPv4 files between start_time and end_time.
		; t4 defines 1-day long block for each file

		i4 = where(start_time LE times AND times LE stop_time AND names EQ ark_name[i] AND activity EQ plot_activity,n4)
		IF n4 EQ 0 THEN destroyvar, t4 ELSE t4 = [1,1]#times[i4]+[0,sec_per_day]#replicate(1,n4)
		i4 = merge_ranges(t4,merge=t4)

		; IPv6 files between start_time and end_time.
		; t6 defines 2-day long block for each file.

		IF plot_activity EQ ACTIVITY_TEAM_PROBING THEN BEGIN
			i6 = where(start_time LE times AND times LE stop_time AND names EQ ark_name[i] AND activity EQ ACTIVITY_TOPO_V6,n6)
			IF n6 EQ 0 THEN destroyvar, t6 ELSE t6 = [1,1]#times[i6]+[0,sec_per_day]#replicate(2,n6)
			i6 = merge_ranges(t6,merge=t6)
		ENDIF ELSE BEGIN
			n6 = 0
			destroyvar, t6
			i6 = merge_ranges(t6,merge=t6)
		ENDELSE

		nb = 1						; Full time range
		tb = [start_time,stop_time]

		; Get time ranges in tb that are not covered by IPv4 and IPv6 files

		IF nb NE 0 AND n4 NE 0 THEN tmp = where_overlap(tb,t4,nn1=nb,rr1=tb)
		IF nb NE 0 AND n6 NE 0 THEN tmp = where_overlap(tb,t6,nn1=nb,rr1=tb)

		; At this point tb (nb pairs) are the time ranges in [start_time,stop_time]
		; that are not covered by IPv4 or IPv4 files (i.e. no data available)

		; Split the ranges with data (IPv4 and/or IPv6) into time ranges
		; where both IPv4 and IPv6 data are available: t46 (n46 pairs),
		; where only IPv4 data are available: tt4 (n4 pairs), and
		; where only IPv6 data are available: tt6 (n6 pairs).

		n46 = where_overlap(t4,t6,overlap=t46,nn1=n4,rr1=tt4,nn2=n6,rr2=tt6)

		; Fill squares for IPv4-only time ranges

		IF n4 NE 0 THEN tt4 = TimeScale(TimeSet(unix=tt4))
		FOR n=0,n4-1 DO polyfill, tt4[[0,1,1,0],n], i+[0,0,1,1], color=color[0]

		; Fill squares for IPv6-only time ranges

		IF n6 NE 0 THEN tt6 = TimeScale(TimeSet(unix=tt6))
		FOR n=0,n6-1 DO polyfill, tt6[[0,1,1,0],n], i+[0,0,1,1], color=color[1]

		; Fill squares with both IPv4 and IPv6 data.

		IF n46 NE 0 THEN t46 = TimeScale(TimeSet(unix=t46))
		FOR n=0,n46-1 DO polyfill, t46[[0,1,1,0],n], i+[0,0,1,1], color=color[2]

		; Explicitly erase squares where there are no data
		; Note that this can become misleading if there are lots of gaps
		; that are smaller than one pixel: if they are shown they will always
		; be one pixel wide.

		IF show_gaps THEN BEGIN
			IF nb NE 0 THEN tb = TimeScale(TimeSet(unix=tb))
			FOR n=0,nb-1 DO polyfill, tb[[0,1,1,0],n], i+[0,0,1,1], color=0
		ENDIF

		nn += [n4,n6,n46]

		IF NOT unlabeled THEN BEGIN	; Print monitor names on right side of graph
			tmp = !x.crange[1]+(0.01+0.13*(i-i/3*3))*(!x.crange[1]-!x.crange[0])
			cname = ark_name[i]
			xyouts, tmp, i, cname, color=!p.color, align=0.0, _extra=_extra	 ;, charsize=charsize
		ENDIF

	ENDFOR

	IF plot_activity EQ ACTIVITY_TEAM_PROBING AND NOT unlabeled THEN BEGIN			; Plot a legend
		FOR i=0,2 DO BEGIN				; IPv4, IPv6, IPv4+IPv6
			IF nn[i] NE 0 THEN BEGIN
				x0 = 0.1
				y0 = 0.94
				x = x0+[0,1,1,0]*0.01+0.25*i
				y = y0+[0,0,1,1]*0.015
				polyfill, /normal, x, y, color=color[i]
				xyouts, /normal, x[1]+0.02,y[1]+0.005, (['IPv4','IPv6','IPv4+IPv6'])[i]+' measurements'+([' only',' only',''])[i], _extra=_extra
			ENDIF
		ENDFOR
	ENDIF

	get_page, $
		silent	= silent	, $
		_extra	= _extra

ENDELSE

RETURN  &  END