## 加载运行环境 --------------------------------

library(imager)

library(imagerExtra)

library(tidyverse)

library(cowplot)

# 设置参数及定义函数------------

figPath <- "木质-1%-9000-3.jpg"

tmpfile <- tempfile(fileext = ".tiff")

tiff(

  filename = tmpfile,

  width = 16,

  height = 8.5,

  units = "cm",

  res = 600,

  compression = "lzw"

)

par(

  oma = c(0, 0, 0, 0),

  mfrow = c(2, 4),

  mai = c(0, 0, 0.1, 0),

  mgp = c(0, 0, 0),

  tck = -0.03,

  xaxt="none",

  yaxt="none",

  pty = "s",

  xpd = NA

)

fig_label <- function(text, region="figure", pos="topleft", cex=NULL, ...) {

  region <- match.arg(region, c("figure", "plot", "device"))

  pos <- match.arg(pos, c("topleft", "top", "topright", 

                          "left", "center", "right", 

                          "bottomleft", "bottom", "bottomright"))

  if(region %in% c("figure", "device")) {

    ds <- dev.size("in")

    # xy coordinates of device corners in user coordinates

    x <- grconvertX(c(0, ds[1]), from="in", to="user")

    y <- grconvertY(c(0, ds[2]), from="in", to="user")

    # fragment of the device we use to plot

    if(region == "figure") {

      # account for the fragment of the device that 

      # the figure is using

      fig <- par("fig")

      dx <- (x[2] - x[1])

      dy <- (y[2] - y[1])

      x <- x[1] + dx * fig[1:2]

      y <- y[1] + dy * fig[3:4]

    } 

  }

  # much simpler if in plotting region

  if(region == "plot") {

    u <- par("usr")

    x <- u[1:2]

    y <- u[3:4]

  }

  sw <- strwidth(text, cex=cex) * 60/100

  sh <- strheight(text, cex=cex) * 60/100

  x1 <- switch(pos,

               topleft     =x[1] + sw, 

               left        =x[1] + sw,

               bottomleft  =x[1] + sw,

               top         =(x[1] + x[2])/2,

               center      =(x[1] + x[2])/2,

               bottom      =(x[1] + x[2])/2,

               topright    =x[2] - sw,

               right       =x[2] - sw,

               bottomright =x[2] - sw)

  y1 <- switch(pos,

               topleft     =y[2] - sh,

               top         =y[2] - sh,

               topright    =y[2] - sh,

               left        =(y[1] + y[2])/2,

               center      =(y[1] + y[2])/2,

               right       =(y[1] + y[2])/2,

               bottomleft  =y[1] + sh,

               bottom      =y[1] + sh,

               bottomright =y[1] + sh)

  old.par <- par(xpd=NA)

  on.exit(par(old.par))

  text(x1, y1, text, cex=cex, ...)

  return(invisible(c(x,y)))

}

# 1. 原始图片--------------------------------

IMG1 <- load.image(figPath)

IMG1 %>% plot()

fig_label("A")

mtext(side=3, line=0, "原图片",cex =0.5)

# 2. 灰度化----------------------

IMG2 <- IMG1 %>% imager::grayscale()

IMG2 %>% plot()

fig_label("B")

mtext(side=3, line=0, "灰度化",cex =0.5)

# 3. 提取内部有效部分--------------------

markbyCircle = function(im){

  w <- width(im)/2

  h <- height(im)/2

  im %>% as.data.frame() %>% 

    dplyr::mutate(value=if_else((x-w)^2+(y-h)^2>(h/1.5)^2,0,value)) %>%

    as.cimg(dim=dim(im)) %>%

    crop.borders(nx=h-h/1.5-10,ny=h-h/1.5-10)

}

IMG3 <- markbyCircle(IMG2)

IMG3 %>% plot()

fig_label("C")

mtext(side=3, line=0, "提取有效部分",cex =0.5)

# 4. 颜色平衡------------------

IMG4 <- BalanceSimplest(IMG3, 1, 1,range=c(0,1)) 

IMG4 %>% plot

fig_label("D")

mtext(side=3, line=0, "颜色平衡",cex =0.5)

# 5. 去杂质--------------------

IMG5 <- DenoiseDCT(IMG4, 0.1,flag_dct16x16 = TRUE) 

IMG5 %>% plot

fig_label("E")

mtext(side=3, line=0, "去杂",cex =0.5)

markbyCircle2 = function(im){

  w <- width(im)/2

  h <- height(im)/2

  mds1 <- expand.grid(x = 1:width(im), y = 1:height(im)) %>% 

    add_column(value =1) %>%

    dplyr::filter((x-w)^2+(y-h)^2>(h-20)^2)

  mds2 <- im %>% as.data.frame() %>% 

    dplyr::select(x,y,value=z,cc)

  full_join(mds1,mds2,by=c("x","y")) %>% 

    dplyr::transmute(x,y,value=if_else(is.na(value.x),value.y,value.x),

                     cc=1) %>% 

    as.cimg(dim=c(width(im),height(im),1,1)) 

}

# 6. 二值化--------------------

IMG6 <- ThresholdAdaptive(IMG5, 0.2,  range = c(0,1)) 

IMG6 %>% plot

fig_label("F")

mtext(side=3, line=0, "二值化",cex =0.5)

# 7. 获取裂纹--------------------

IMG7 <- 1-markbyCircle2(IMG6)

IMG7 %>% plot()

fig_label("G")

mtext(side=3, line=0, "获取裂纹部分",cex =0.5)

# 8. 连通域扫描

IMG8.list <- imager::split_connected(IMG7)

IMG8.len <- IMG8.list %>% map_dbl(~sum(.x==1))

IMG8.list[order(IMG8.len,decreasing = TRUE)[1:30]] %>% 

  reduce(.f = `+`) %>% plot()

fig_label("G")

mtext(side=3, line=0, "连通域扫描",cex =0.5)

dev.off()

shell.exec(tmpfile)

运行结果如下：