2011年5月27日星期五

数码单反(DSLR)相机raw格式图片的读取和处理

随着消费用数码单反(DSLR)相机成本的下降和性能的提升,在很多要求不太高的科研场合DSLR相机可以用作昂贵的高档CCD数码相机的替代品。现在最新的数码相机可以提供14位A/D转换,而动态范围可以到达11档光圈以上【参考链接:12】,很多时候这已经可以满足要求了。一般摄影爱好者使用的JPEG图片格式只能存储8位的层次,而且JPEG图片压缩过程中会有失真,所以作为科研使用的DSLR相机必须使用raw图片格式。

一般各个消费用相机厂商都使用自己专有的raw格式,所以处理起来不像科研级相机那么直接。dcraw是一个Linux平台上的能将各主要品牌相机的raw格式图片转为Netpbm格式图片的小工具。Netpbm提供了一些工具软件,另外很多图像处理软件也包含对Netpbm格式的支持,所以可以在此基础上对图片进行进一步处理。

很多时候我们可能希望自己写代码来处理图片,因为:
  1. 很多图像处理软件只支持8位(或24位RGB)而不支持16位图片;
  2. 科研上用到的算法大都比较复杂,现有的软件不满足需求。
Netpbm可以把图片存成二进制文件或者纯文本(ascii)文件。无论哪种方式,开始的几行都是标识文件格式的纯文本信息,后面跟着是文件内容。

纯文本文件的好处是数值大小一目了然,但是文件比较大。二进制文件比较小,但是处理的时候需要注意:Netpbm存储文件时使用的是Big-endian方式,而一般x86机器都是使用Little-endian方式,所以在读取数据时要进行转换。glibc的endian.h提供这样的转换函数。二进制文件存储顺序是从上到下按行存储,每行内按从左到右顺序存储。

dcraw转换图片时默认使用8位方式,如果要使用16位方式,需要使用参数 -6

2011年5月24日星期二

《Maxima快速参考手册1.1》

这是之前写的一份Maxima中文文档的更新版本。

2011年4月15日星期五

Vim笔记(二)

  1. 分割窗口
    :split file2    把窗口分割成上下两半,在新开窗口打开文件file2
    :new    上下分割窗口,在新开窗口打开新文件
    CTRL-W w    切换窗口
    CTRL-W +    增大窗口尺寸
    CTRL-W -    缩小窗口尺寸
    :vsplit file2    左右分割窗口,在新开窗口打开文件file2
    CTRL-W h    切换到左边窗口,其他方向按j, k, l类推
    CTRL-W H    移动窗口到最左边,其他方向按J, K, L类推
    :qall    退出所有窗口,类似的有:wall
  2. 反复执行复杂命令
    q{register}    开启寄存器,开始记录命令,register可以是az中任何一个
    q    结束记录命令
    @{register}    执行寄存器{register}内所有命令
  3. 查找替换一段字符
    :[range]s/from/to/[flags]    将from替换为to[range]控制替换范围,[flags]表示一些控制参数
    常用[range]
    • .    当前行
    • $    文件最后一行
    • 1,5    第一行到第五行
    • %    全部文件
    • .+3    当前行之下第三行
    常用[flags]
    • g    替换所有匹配字符串,否则只替换第一个匹配
    • c    执行每一个替换之前寻求确认
  4. 可视化块模式
    CTRL-V    开启可视化块模式,用h,j,k,l选择块
    I{string}<Esc>    在块的左边每一行插入一段字符串
    A{string}<Esc>    在块的右边每一行加入一段字符串
    c{string}<Esc>    替换一个块的文字
    U    换为大写字母
    u    换为小写字母
    ~    大小写互换
  5. 数行合并
    J    数行合并,去除换行符

2011年4月14日星期四

Vim笔记(一)

常用的Vim命令小结,最基本的就不说了,例如如何安装,如何启动和退出,两种基本模式(正常模式和插入模式),用h,j,k,l移动光标,简单的复制粘贴等等。
  1. 移动光标
    1. 以字符为单位
      f    右移至指定字符,例如fh移动至右边第一个h
      F    左移至指定字符
    2. 以词为单位
      w    右移至词首
      b    左移至词首
      e    右移到词尾
      ge    左移到词尾
    3. 以行为单位
      0    左移至行首
      ^    左移至行首第一个非空字符
      $    右移至行尾
      gg    移动至第1行
      G    移动至最后1行
      nG    移动至第n行
      n%    移动至n%位置
    4. 以屏为单位
      H    移动至屏幕顶端
      M    移动至屏幕中间
      L    移动至屏幕底端
      CTRL-U    向上滚半屏
      CTRL-D    向下滚半屏
      CTRL-F    向下滚一整屏
      CTRL-B    向上滚一整屏
  2. 简单搜索
    /string    向后搜索字符串
    ?string    向前搜索字符串
    n    下一个匹配字符串
    N    上一个匹配字符串
    /\<word\>    匹配整个单词,\<\>分别匹配词首和词尾
    %    寻找匹配的括号
  3. 书签
    ma    标记此处为a,共可使用az二十六个书签
    `a    回到标记为a

2011年3月4日星期五

几个关于Gnuplot小技巧的网页

Gnuplot手册不能告诉你的小技巧。Gnuplot也可以画出很惊艳的图片。
  1. not so Frequently Asked Questions
    http://t16web.lanl.gov/Kawano/gnuplot/index-e.html
  2. Gnuplot tricks blog
    http://gnuplot-tricks.blogspot.com/
  3. Impossible gnuplot graphs
    http://www.phyast.pitt.edu/~zov1/gnuplot/html/intro.html
  4. Gnuplot surprising
    http://gnuplot-surprising.blogspot.com/
  5. Gnuplotting
    http://www.gnuplotting.org/

2011年1月12日星期三

Octave笔记

Octave可以看作是开放源代码版本的Matlab,它的命令几乎和Matlab完全一样,因此可以用Matlab的书籍和文档来学习Octave。另外,绝大多数Matlab的代码也可以不经修改在Octave中运行。对于初学者来说,下面这些笔记可以帮助快速上手Octave。更多的功能,可以参考Octave文档或者其他Matlab的资料。
  1. 脚本文件,以下面一行开始:
    #!/usr/bin/octave -qf
  2. 注释,以%或者#开始
  3. 虚数单位:i, j, I, J都可以 两种方法定义复数:3 + 4i或者complex(3, 4)
  4. 定义矩阵:
    octave:1> a = [1, 2; 3, 4]
    a =
    
    1   2
    3   4
  5. 随机数:rand()或者rand(m, n)
  6. 定义一段连续序列:
    octave:1> 1:5
    ans =
    
    1   2   3   4   5
    还可以指定步长:
    octave:2> 1:2:10
    ans =
    
    1    3    5    7    9
  7. 字符串:单引号或者双引号
  8. 定义Data Structure:
    octave:1> x.name = ”Yusuf”
    octave:2> x.age = 33
    octave:3> x.phone = “46103”
    octave:4> x
    x =
    {
    name = Yusuf
    age = 33
    phone = 46103
    }
  9. 定义Cell Array:
    octave1:> a = {”random matrix”, rand(2)}
    a =
    
    {
    [1,1] = random matrix
    [1,2] =
    
    0.417140   0.460380
    0.337105   0.078835
    
    }
    
    octave2:> a{2}(1,2)
    ans =  0.460380
  10. 声明全局变量:
    global c = 3e8
  11. 算术符号
    x + y 加法
    x - y 减法
    x * y 矩阵乘法
    x .* y 元素相乘
    x / y 右除,相当于(inverse (y’) * x’)’
    x ./ y 元素右除
    x \ y 左除,相当于inverse (x) * y
    x .\ y 元素左除,每个y元素被相应的x元素除
    x ^ yx ** y 乘方
    x .^ yx .** y 元素乘方
    x’ 转置复共轭
    x.’ 转置
  12. 比较符号
    x < y 小于
    x <= y 小于等于
    x == y 等于
    x >= y 大于等于
    x > y 大于
    x != yx ~= y 不等于
  13. if语句
    • 形式一:
      if (condition)
        then-body
      endif
    • 形式二:
      if (condition)
        then-body
      else
        else-body
      endif
    • 形式三:
      if (condition)
        then-body
      elseif (condition)
        elseif-body
      else
        else-body
      endif
  14. switch语句
    switch expression
      case label
        command_list
      case label
        command_list
      …
      otherwise
        command_list
    endswitch
  15. while语句
    while (condition)
      body
    endwhile
  16. do-until语句
    do
      body
    until (condition)
  17. for语句
    for var = expression
      body
    endfor
  18. 定义函数
    • 形式一:无自变量 function name body endfunction
    • 形式二:带自变量 function name (arg-list) body endfunction
    • 形式三:单个返回值 function ret-var = name (arg-list) body endfunction
    • 形式四:多个返回值 function [ret-list] = name (arg-list) body endfunction
  19. 可变自变量长度函数 function val = smallest (varargin) val = min ([varargin{:}]); endfunction
  20. 函数Handles(指向函数的指针) @function-name
  21. 匿名函数 @(argument-list) expression
  22. 终端输出 disp (x)
  23. 写读包含分隔符的数据文件 dlmwrite (file, a) dlmwrite (file, a, delim, r, c) dlmwrite (file, a, key, val ...) dlmwrite (file, a, "-append", ...) a是待写矩阵,delim是分隔符,默认为逗号。r, c分别为起始的空行和空列。key包含下列参数:
    • "append" 即上面的"-append","on"或者"off"
    • "delimiter" 即上面的delim
    • "newline" 换行符
    • "roffset" 即上面的r
    • "coffset" 即上面的c
    • "precision" 精确度,可以用类似fprintf的格式字符串,或者有效数字位数
    data = dlmread (file) data = dlmread (file, sep, range) range为四元素的矢量:[R0, C0, R1, C1],为数据块左上角和右下角的标记,从0开始
  24. Octave可以使用大多数类似C的I/O函数(fprintf,fopen等)
  25. 常用数学函数
    exp (x)指数函数sinh (x) 双曲正弦
    log (x)自然对数cosh (x) 双曲余弦
    log10 (x)以10为底的对数 tanh (x) 双曲正切
    log2 (x)以2为底的对数 coth (x) 双曲余切
    sqrt (x) 平方根 sech (x) 双曲正割
    abs (z) 绝对值 csch (x) 双曲余割
    arg (z) 辐角 asinh (x) 反双曲正弦
    conj (z) 复共轭 acosh (x) 反双曲余弦
    real (z) 复数实部 atanh (x) 反双曲正切
    imag (z) 复数虚部 acoth (x) 反双曲余切
    sin (x) 正弦 asech (x) 反双曲正割
    cos (x) 余弦 acsch (x) 反双曲余割
    tan (x) 正切 sum (x) 求和
    cot (x) 余切 prod (x) 求积u
    sec (x) 正割 ceil (x) 不小于x的最小整数
    csc (x) 余割 floor (x) 不大于x的最大整数
    asin (x) 反正弦 factorial (n) 阶乘
    acos (x) 反余弦 max (x) 最大值
    atan (x) 反正切 min (x) 最小值
    acot (x) 反余切 round (x) 四舍五入
    asec (x) 反正割 sign (x) 符号函数
    acsc (x) 反余割
    atan2 (y, x) 计算atan (y / x)
    sind (x),… 以度数为单位的三角函数
  26. 特殊函数
    [a, ierr] = airy (k, z, opt) Airy函数
    [j, ierr] = besselj (alpha, x, opt) 第一类Bessel函数
    [y, ierr] = bessely (alpha, x, opt) 第二类Bessel函数
    [i, ierr] = besseli (alpha, x, opt) 变形第一类Bessel函数
    [k, ierr] = besselk (alpha, x, opt) 变形第二类Bessel函数
    [h, ierr] = besselh (alpha, k, x, opt) 第一类(k=1)和第二类(k=2)Hankel函数
    beta (a, b) Beta函数
    gamma (z) Gamma函数
    bincoeff (n, k) 二项式系数
    erf (z) 误差函数
    legendre (n, x, normalization) n阶Legendre函数,m = 0…n
  27. 线性代数
    det (a) 行列式
    eig (a) 本征值
    dot (x, y) 点乘
    inv (a) 逆矩阵
  28. 解非线性方程 [x, fval, info] = fsolve (fcn, x0, options) fcn是包含待解方程组的函数,x0是自变量初始值,x为返回的数值解,fval为返回的函数值,如果结果收敛,info为1,否则为其他值。例子: function y = f (x) y(1) = -2*x(1)^2 + 3*x(1)*x(2) + 4*sin(x(2)) - 6; y(2) = 3*x(1)^2 - 2*x(1)*x(2)^2 + 3*cos(x(1)) + 4; endfunction [x, fval, info] = fsolve (@f, [1; 2]) 结果: x = 0.57983 2.54621 fval = 6.1872e-08 -3.2708e-07 info = 1
  29. 数值积分(各函数对应不同算法) [v, ier, nfun, err] = quad (f, a, b, tol, sing) quadl (f, a, b, tol) quadgk (f, a, b, prop, val, ...) quadl (f, a, b, tol)