-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathREADME.Rmd
90 lines (65 loc) · 3.73 KB
/
README.Rmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
---
output:
github_document: default
---
<!-- README.md is generated from README.Rmd. Please edit that file -->
```{r, echo = FALSE}
knitr::opts_chunk$set(
collapse = TRUE,
comment = "#>",
fig.path = "man/figures/"
)
```
# mincountr
The goal of mincountr is to provide a simple point-counting mechanism for microscope thin section images of rocks
## Installation
You can install mincountr from github with:
```{r gh-installation, eval = FALSE}
# install.packages("devtools")
devtools::install_github("Ignimbrit/mincountr")
```
## Example
Say you have collected an image from a mineral-bearing thin section e.g. in an electron microprobe session and wonder about the relative share of the different phases in your sample. Your image might look something like this:
```{r load image}
library(mincountr)
myimage <- mcr_load_image(
system.file("extdata", "testim.png", package = "mincountr")
)
plot(myimage)
```
As you can see, the different levels of brightness in the image allows the observer to distinguish between several distinct phases. There is a large, very bright mineral, some lightgrey minerals, a darkgrey, glassy matrix, some black holes and so on. With `mincountr` we are able to translate this qualitative optical assessment into practical numbers. First let's have a look at the density-distribution (like a continuous histogram) of the images brightness.
```{r plot brightness}
mcr_inspect_phases(myimage)
```
Here we can see that the four phases we distinguished in the image above show up in the density distribution as distinct peaks of specific brightness values. We can use this to assign a brightness-range to certain phases.
The peak on the far left (most dark) corresponds likely to the hole we've seen in our thin section image. The peak's value range lays between value ~0-0.05. Other peaks (from left to right) range from 0.3-0.45 (glassy matrix?), from 0.5-0.65 (light gray minerals?) and from 0.92-1 (bright minerals?). Let me illustrate what I mean:
```{r illustrate peakborders}
mcr_inspect_phases(myimage) +
ggplot2::geom_vline(
xintercept = c(0, 0.05, 0.3, 0.45, 0.5, 0.65, 0.92, 1),
color = "red"
)
```
Now in this example we just chose the borders of the peak by hand. This is probably the safest method of constraining your brightness-levels, as it allows you to chip in your personal mineralogical expertise. However, sometimes you have a large stack of images you want to work with and probably not the time to constrain peak-ranges by hand every single time. This is why mincountr comes with an automatic mechanism to generate those numbers for you.
```{r automatic constrains}
myconstrains <- mcr_autoconstrain(myimage)
print(myconstrains)
```
The `mcr_autoconstrain` function automatically detects peaks and notes their position (`peakpos`) and then goes on and calculates their borders both on the left-hand-side (`x1`) and on the right hand side (`x2`). Under the hood, `mcr_autoconstrain` identifies turning points in the brightness-"spectra", cuts it into pieces, one piece per peak, and then loops over the single-peak spectra-pieces to calculate the half-height-width.
Now we have all the information we need to assign certain areas of our original image to distinct phases. The mincountr-package comes with a function that lets you inspect what this assignment looks like.
```{r check assignment}
mcr_inspect_assignment(
myimage,
lhs = myconstrains$x1,
rhs = myconstrains$x2
)
```
As every pixel in the original image was now assigned to one of the 7 levels shown in the picture above, we can go ahead and just count the pixels and then calculate the relative share of each group.
```{r sum up}
myresult <- mcr_herd_minerals(
myimage,
lhs = myconstrains$x1,
rhs = myconstrains$x2
)
print(myresult)
```