dc_intervals <- function(dat, threshold=0.02/100, price="price", time="time", vol="amount") {
  # candidates for confirmed directional change
  dcc_up_cand <- c(FALSE, c(dat[2:nrow(dat), price]<dat[1:(nrow(dat)-1), price]*(1-threshold)))
  dcc_down_cand <- c(FALSE, c(dat[2:nrow(dat), price]>dat[1:(nrow(dat)-1), price]*(1+threshold)))
  
  # initialize new columns
  dat[,"interval"] <- NA  # 4 types of intervals: down, down_os, up, up_os (os=overshoot, NA means unknown)
  dat[,"point"] <- NA     # 3 types of points: dcc_up, dcc_down,  ep (NA means ignore)
  
  # threshold too large?
  if(all(!dcc_up_cand)) stop("threshold=",threshold, " too large, no price moves above this threshold")
  if(all(!dcc_down_cand)) stop("threshold=",threshold, " too large, no price moves above this threshold")
  
  # figure out if the first change is up or down
  first_up <- which.max(dcc_up_cand) 
  first_down <- which.max(dcc_down_cand)
  
  if(first_up<first_down) {
    dat[first_up,"point"] <- "dcc_up"
    mode <- "dcc_up"
    i <- first_up
  } else {
    dat[first_down, "point"] <- "dcc_down"
    mode <- "dcc_down"
    i <- first_down
  }
  
  # iterate till the end of the data
  while (i<nrow(dat)) {
    if(mode=="dcc_up") {
      j <- which.max(dcc_down_cand[i+1:length(dcc_down_cand)])+i
      if(!dcc_down_cand[j]) break # end of dataset
      dat[i,"point"] <- "dcc_up"
      dat[i,"interval"] <- "uoi"
      if((j-1)-(i+1)>0) { # if not changing direction immediately
        ep_idx <- which.min(dat[(i+1):(j-1), price])+i
        dat[ep_idx, "point"] <- "ep"
        dat[i:ep_idx-1, "interval"] <- "uoi"
        dat[ep_idx:j,"interval"] <- "dei"
      } 
      i <- j+1
      mode <- "dcc_down"
    } else { # mode=="dcc_down"
      j <- which.max(dcc_up_cand[i+1:length(dcc_up_cand)])+i
      if(!dcc_up_cand[j]) break # end of dataset
      dat[i,"point"] <- "dcc_down"
      dat[i,"interval"] <- "doi"
      if((j-1)-(i+1)>0) { # if not changing direction immediately
        ep_idx <- which.max(dat[(i+1):(j-1), price])+i
        dat[ep_idx, "point"] <- "ep"
        dat[i:ep_idx-1, "interval"] <- "doi"
        dat[ep_idx:j,"interval"] <- "uei"
      }
      i <- j+1
      mode <- "dcc_up"      
    }
  }
  
  dat <- dat[!is.na(dat[,"point"]),] # drop all points that are not dcc points
  return(dat)
}

ohlc_summary <- function(dat, browse=TRUE, verbose=TRUE) {
  dat$hour <- strftime(dat$time, format="%Y-%m-%d %H")
  one_hour <- function(dset) {
    i1 <- which.min(dset$time)
    i2 <- which.max(dset$time)
    ans <- data.frame(
      time=strptime(dset[1,"hour"], format="%Y-%m-%d %H"),
      open=dset[i1, "price"],
      hi=max(dset$price),
      lo=min(dset$price),
      close=dset[i2,"price"],
      amount=sum(dset$amount)
    )
    return(ans)
  }
  ret <- do.call(rbind, lapply(split(dat, dat$hour), FUN=one_hour))
  if(verbose) message("ohlc summary built.")
  if(browse) browser()
  return(ret)
}

# begin, end ... "YYYY-MM-DD"
plot_dc <- function(dc_dat=NULL, ohlc=NULL, ticks=NULL, begin=NULL, end=NULL, lwd=3, price="price", time="time", vol="amount", ...) {
  if(is.null(ohlc) & !is.null(ticks)) ohlc <- ohlc_summary(ticks) 
  if(is.null(dc_dat)  & !is.null(ticks)) dc_dat <- dc_intervals(dat=ticks, time=time, price=price, vol=vol)
  if(is.null(dc_dat)) stop("missing dc data.")
  
  if(is.null(end)) end <- max(dc_dat[,time])
  if(is.character(end)) end <- as.POSIXct(strptime(end, "%Y-%m-%d"))
  if(is.null(begin)) begin <- min(dc_dat[,time])
  if(is.character(begin)) begin <- as.POSIXct(strptime(begin, "%Y-%m-%d"))
  
  plot(
    dat=1, col="white", dc_dat=dc_dat, xlab="", 
    xlim=c(min(ohlc$time), max(ohlc$time)), ylim=c(min(ticks$price), max(ticks$price)),
    main=sprintf("%.2f%% Direction Changes for BTC", threshold*100),
    cex.main=1
  )
  
  if(!is.null(ohlc)) {
    for (i in seq(1, nrow(ohlc))) {
      segments(
        x0=ohlc[i,"time"], y0=ohlc[i,"lo"],
        x1=ohlc[i,"time"], y0=ohlc[i, "hi"],
        col="darkgrey", lwd=0.5
      )
    }
  }
  
  dcc_idx <- which(dc_dat$point=="dcc_up"|dc_dat$point=="dcc_down")
  for (i in setdiff(dcc_idx,1)) {
    segments(
      x0=dc_dat[i-1,time], y0=dc_dat[i-1,price],
      x1=dc_dat[i,time], y1=dc_dat[i,price],
      col="red", lwd=lwd
    )
  }
  for (i in setdiff(dcc_idx,nrow(dc_dat))) {
    segments(
      x0=dc_dat[i,time], y0=dc_dat[i,price],
      x1=dc_dat[i+1,time], y1=dc_dat[i+1,price],
      col="darkgreen", lwd=lwd
    )
  }
  axis(1,at=dc_dat[dcc_idx, time],labels=FALSE)
}