2010年1月21日星期四

统计笔记(一)

  1. 协方差(covariance)定义:
    Cov( Y1 , Y2 )=E[( Y1 - μ1 )( Y2 - μ2 )]

    其中 E()表示期望值, μ1 =E( Y1 ) μ2 =E( Y2 )
  2. 相关系数(coefficient of correlation)定义:
    ρ= Cov( Y1 , Y2 ) σ1 σ2

    其中 σ1 σ2 分别是 Y1 Y2 的标准差。
  3. 任意分布样本的平均值
    Y1 Y2 ,..., Yn E( Yi )=μ V( Yi )= σ2 的独立随机变量,则平均值
    Y ¯ = 1 n i=1 n Yi

    的期望值 E( Y ¯ )=μ,方差 V( Y ¯ )= σ2 /n
  4. 正态分布样本的平均值
    Y1 Y2 ,..., Yn 是来自正态分布(期望值为 μ,方差为 σ2 )的样本,其平均值
    Y ¯ = 1 n i=1 n Yi

    也呈正态分布,期望值 μ Y ¯ =μ,方差 σ Y ¯ 2 = σ2 /n
  5. Y1 Y2 ,..., Yn 同上,那么
    i=1 n Zi 2 = i=1 n ( Yi -μ σ )2

    呈自由度为 n χ2 分布。
  6. Y1 Y2 ,..., Yn 同上,那么
    (n-1) S2 σ2 = 1 σ2 i=1 n( Yi - Y ¯ )2

    呈自由度为 n-1 χ2 分布。其中,
    S2 = 1 n-1 i=1 n( Yi - Y ¯ )2

    是样本方差。 S2 Y ¯ 为独立随机变量。
  7. 学生 t分布(简称 t分布)定义:
    如果 Z是正态分布随机变量, W是自由度为 ν χ2 分布随机变量,如果 Z W互相独立,则
    T= Z W/ν

    呈自由度为 ν t分布。
  8. 在上面的定义中,如果让 Z=n( Y ¯ -μ)/σ W=(n-1) S2 / σ2 ,则由上面第4和第6条结论可以得出
    n( Y ¯ -μ S )

    呈自由度为 n-1 t分布。
  9. F分布定义:
    如果 W1 W2 为独立 χ2 分布变量,自由度分别为 ν1 ν2 ,则
    F= W1 / ν1 W2 / ν2

    被称作具有 ν1 分子自由度和 ν2 分母自由度的 F分布。
  10. 中心极限定理(Central Limit Theorem)
    Y1 Y2 ,..., Yn 是独立的具有同样分布的随机变量,期望值 E( Yi )=μ,方差 V( Yi )= σ2 <。定义
    Un =n( Y ¯ -μ σ )

    其中
    Y ¯ = 1 n i=1 n Yi

    那么,当 n时, Un 趋近于正态分布。

2010年1月18日星期一

写了一份Maxima的中文文档

以前我为计算机代数系统Maxima写了一些简介,这次抽时间写了一个比较完整的中文文档,希望能为开放源代码软件的推广做点贡献。

《Maxima快速参考手册》 写的仓促,欢迎反馈。

2009年5月11日星期一

在eps图片上添加LaTeX公式

利用MetaPostexteps宏包,可以在metapost文件中加入eps图片,同时利用metapost的label命令在图片任何地方加入LaTeX公式。下面是一个使用exteps的例子:
input exteps
verbatimtex
\documentclass{article}
\begin{document}
etex
prologues:=2;
beginfig(1);
begineps "figure.eps";
endeps;
label.rt(btex $f(x) = \sqrt x$ etex, (2cm,3cm));
endfig;
end
为了使用LaTeX命令,需要设置环境变量
TEX=latex

2009年4月2日星期四

GnuPlot在PostScript终端下使用任意字体

有时候我们在GnuPlot中想使用一些其他字体,例如LaTeX字体或者一些艺术字体,这时候我们可以用fontfile选项把自己的字体内嵌进生成的PostScript文件:
set term postscript fontfile "myfont.ext" "FontName"
这里的字体可以是Type 1字体,也可以是TrueType字体。命令必须包含完整字体文件名,包括扩展名。字体文件可以在当前目录,也可以用GNUPLOT_FONTPATH环境变量控制字体路径。FontName(字体名)不是字体文件名,通常不那么容易直接得到,尤其是TrueType字体。最简单的办法是在GnuPlot里执行
set term postscript fontfile "myfont.ext"
GnuPlot会自动给出字体名。

这样,由于可以使用LaTeX字体,配合PostScript终端的enhanced选项,我们几乎可以在GnuPlot中输入任意特殊字符。另外,为使用TrueType字体,需要安装ttf2pt1工具进行字体自动转化。

2009年3月16日星期一

sed编辑器的使用

sed被称为流编辑器(stream editor)。它和一般写字板之类的文本编辑器不同,它通过逐行读入文本文件,按照用户指定的方式对文本进行过滤。用户通过脚本命令告诉sed该如何处理文本文件,所以sed可以看作一种程序语言。这看起来没有使用鼠标和键盘编辑文件的编辑器方便,但是当我们需要自动处理大量文本文件时,这会非常方便。 例如,我需要把很多数据文件的某些行注释掉,也就是在所有数据文件的某几行开始处加上字符“# ”。如果一个个文件编辑,工作量就太大了。下面这个命令可以很方便的实现这个功能:
sed -i '1,5 s/^/# /' *.txt
这个命令把每一个.txt作为后缀的文件前5行开始处加上“# ”。 下面把常用的一些sed命令总结一下。
  1. 替换一段文字
  2. 这个恐怕是sed最常用的方式了。
    sed 's/a/A/' <old >new
    这个命令把文件old中的小写字幕a换为大写字幕A,并且把新的文件存为new。但是这个命令只能替换第一个a,如果想要把文件old中所有的a都替换掉,需要加上一个参数:
    sed 's/a/A/g' <old >new
    上面命令中的斜杠/是分隔符,但是sed没有要求一定用斜杠作为分隔符,其实可以用任何字符(只要不和搜索字符矛盾)作为分隔符,例如冒号等等。
  3. 执行多个命令
  4. 上面的例子中只执行了一个命令,如果需要同时执行多个命令,需要用到参数-e
    sed -e 's/a/A/' -e 's/b/B/' <old >new
  5. 指定对某些行操作
  6. 如果只对某些行进行操作,可以像一开始举的例子那样,使用行号指定范围:
    sed '1,10 s/A/a/'
另外,这里有一个很有用的sed常用单行命令参考

2009年1月9日星期五

GnuPlot使用png terminal时找不到字体的问题

GnuPlot不同的终端(terminal)使用不同的处理字体的方法。png/jpeg/gif使用libgd,所以字体目录由GDFONTPATH这个环境变量控制,可以直接使用TrueType字体。例如把需要用到的字体放在$HOME/fonts目录下,在.bashrc中添加:
export GDFONTPATH=$HOME/fonts
当路径设置好后,可以用GNUPLOT_DEFAULT_GDFONT这个环境变量设置默认字体。

2009年1月2日星期五

在USB drive(U盘)上安装Linux

UNetbootin可以很方便的把Linux装到U盘上,只要主板支持,就可以从U盘启动Linux。由于U盘可读写,比Live CD更方便。目前UNetbootin支持Fedora,Debian,Ubuntu等各种主流Linux版本。如果你的U盘比较小,可以试试Puppy Linux,挺小巧的一个Linux版本,只有100M左右大小。

2008年11月21日星期五

Mie散射计算器

Mie散射理论为球形颗粒对于电磁波的散射提供了一套严格解,其数值计算十分复杂。Bohren和Huffman的书里提供了用于数值计算的Fortran代码。另外,Oregon Medical Laser Center的网页提供了一个在线的Mie散射计算器,可以快速的计算出一些常用参数,并且给出一些图线。

2008年11月19日星期三

GNU screen的几个常用命令

GNU screen是命令行界面下的窗口管理器,允许用户在命令行界面下同时在几个虚拟窗口运行不同的程序。GNU screen有一个很有用的功能:如果用户登录到远程机器上利用screen运行程序,该程序在用户切断连接之后可以继续运行,并且下次用户可以登录继续上次的操作。下面是一些常用的命令。

首先在终端执行screen命令进入程序。执行之后虽然看不到任何变化,其实screen已经启动,我们看到的就是一个虚拟终端。
  1. 建立一个新的窗口(虚拟终端):C-a c
    也就是按Ctrl+a然后再按c
  2. 切换到之后的窗口:C-a n
  3. 切换到之前的窗口:C-a p
  4. 切换的第N个窗口:C-a N
  5. 显示现有窗口列表,用上下箭头(或者类似于vi中用j和k)移动光标,按Enter选择:C-a "
  6. 改变当前窗口的名字:C-a A
  7. 关闭目前的screen但是保留其中运行的程序(detach):C-a d
    如果直接关闭窗口也可以有同样作用。
  8. 运行screen继续上次在screen里运行的程序(reattach):screen -r
  9. 进入拷贝模式:C-a [
    在screen里不能向上翻屏,进入拷贝模式可以使用光标上下移动翻屏。用空格键可以选择要拷贝的文字,第一次按空格选择开始的字符,第二次选择结束的字符。
  10. 粘贴选择的文字:C-a ]

2008年11月16日星期日

Linux下用sopcast看电视

Sopcast是一个基于P2P技术的网络电视软件,在Linux和Windows上都可以运行。Sopcast网站提供的Linux版本倒挺简单,直接下载解压就可以用了。但是这只是一个命令行的软件。有几个图形界面可以用,觉得gsopcast不错,但是源代码在Fedora和Debian上编译都有问题,header文件需要改动几处。可能因为大家用的编译器不一样?编译好的文件在这里,是我在Debian下用gcc4.3编译的,gsopcast版本是0.4.0,可以直接下载下来用,应该不依赖特定的发行版。

如果自己在桌面建立一个gsopcast快捷方式图标的话,似乎当从图标启动的时候,$PATH变量不起作用,所以gsopcast可能找不到sopcast的执行文件。另外,新的sopcast执行文件是sp-sc-auth,而不是以前的sp-sc。要保证正常使用,最好在/usr/local/bin目录建立一个符号链接,取名作sp-sc。

2008年11月10日星期一

光圈档数(aperture stop)的计算方法

摄影的时候计算曝光量经常用到几“档”曝光这个术语。这个曝光的“档数”(stop)也经常被称作EV(exposure value)。这个“档数”到底怎么算出来的? 曝光量是由光圈和快门共同决定的。曝光量增加一倍,就叫做曝光增加一档。当光圈一定时,曝光量和快门时间(也就是曝光时间)成正比,例如1/30秒就比1/60秒曝光增加一倍,也就是增加一档。这比较容易理解。比较容易迷惑的是光圈。 光圈一般用F数(F number或者F/#)来表示,其定义为:
f/#=fD
其中f是镜头焦距,D是镜头入射孔径。如果拿起一个镜头,从前面看进去,所看到的那个圆形小孔的直径,就是镜头的入射孔径了。对于一般相机镜头,这个孔径大小是可调的,由f/#的设置来控制。有些人可能会想,既然计算曝光量,为什么不直接拿孔径大小D来算,还要定义一个什么f/#。这是因为,感光器件(胶片或者CCD)的曝光量不仅和孔径大小有关,还和焦距有关。用f/#能更准确的表示曝光量。但是这有一个前提,就是物体离镜头比较远,如果物体离镜头很近(例如显微镜),就不能用f/#来算曝光量了。这时候一个更常用的量是数值孔径(numerical aperture或NA)。 因为孔径D在分母上,所以f/#数值越小,光圈越大。例如光圈f/2比f/4要大。但是f/2并非比f/4大一档光圈,而是二档。因为镜头通光量和孔径面积成正比,也就是和孔径大小的平方成正比。所以,比f/4大一档的光圈是
422.8
实际上还会经常碰到分数表示,例如半档,1/3档。对于普遍情况下,可以用下面的公式计算两个光圈究竟差几档:
N=2log2f/#2f/#1
例如,f/1.4和f/1.8相差
N=2log2f/1.8f/1.423

2008年11月4日星期二

LaTeX CJK生成PDF文件的中文搜索以及复制粘贴

以前的blog介绍过如何在Fedora下安装配置LaTeX CJK,其他的Linux发行版也可以用类似的方法。但是如果直接用ps2pdf或者dvipdfm生成PDF文件,这样的PDF文件没办法搜索中文,也不能复制粘贴中文内容。可以用dvipdfmx将生成的dvi文件转换为PDF,这样生成的PDF可以完成上述功能。 另外,我把自己的LaTeX CJK的字体和配置文件放到了网上,地址是http://www.mediafire.com/file/dwutwxjjmdn/texmf.tar.gz,有需要的朋友可以去下载。这个文件包含CJK包和三个字体,分别是Fedora自带的文鼎宋体,楷体,还有文泉驿的黑体,应该满足大多数应用。下载后解压到$HOME目录,执行texhash $HOME/texmf。不需要进行任何其他配置。下面的文件可以用来测试安装结果:
\documentclass{article}
\usepackage{CJKutf8}
\begin{document}
\begin{CJK}{UTF8}{song}
CJK 测试。

\emph{CJK 测试。}

\textbf{CJK 测试。}
\end{CJK}
\end{document}

2008年10月31日星期五

用WebPlayer在网页里播放mp3

WebPlayer利用Flash技术,允许你在你的网页里直接播放mp3文件,而不需要任何其他软件。WebPlayer使用非常简单,只要下载一个很小的文件放在你的网站上,把要播放的mp3文件也放上去,再在你的网页的HTML文件里加几行代码,就可以了。 下面是我用口琴演奏的《天鹅》。
祝老婆纪念日快乐!

2008年10月29日星期三

光束质量因子M2(M squared)的测量

各种文献经常会用M2(M squared)因子来评价光束的质量。到底什么是M2?一般文献会告诉你,M2就是一个光束的腰半径和发散角的乘积跟TEM00高斯光束相比的比值。M2≥1,因为TEM00高斯光束的该乘积具有最小值。要想从实验得到M2的值,不可避免的要涉及光束半径测量的问题。
  1. 光束半径
  2. 光束没有明显的边界,因此它的半径有很多不同定义方法。对于高斯光束,最常见的定义就是1/e2半径,也就是把光强下降到最大值的1/e2的位置定义为半径,通常在公式里以W表示。在统计学里,我们对高斯分布的宽度(半径)使用的是方差,也就是说,对于高斯分布函数
    f(x)=ae-x22σ2
    它的方差(二次矩)半径为σ。把这个公式和一维高斯光束光强公式
    I(x)=e-2x2wx2
    比较,就会发现W=2σ。无论是W还是σ,实验上都没办法直接测量。但是如下定义的宽度实验上是可以测量的:用一个刀片,横向扫描过光束,这样一部分光束被刀片遮挡,另外一部分会通过。如果通过的光强等于总光强的90%时刀片坐标为x1,通过的光强等于总光强的10%时刀片坐标为x2,这样光束的宽度可以用
    D=|x2-x1|
    表示。当然,这个比例不局限于10%和90%,也可能是其他值。这种实验方法称为“刀锋(knife edge)法”。这样的话,我们需要找出一个D和W(或者σ)的关系。可惜,这个关系依赖于光束的光强分布,所以不存在一种普遍适用的关系。但是对于TEM00高斯光束,这个关系不难得到。
  3. TEM00高斯光束
  4. 因为刀片在一维方向移动,所以我们只要对光强进行一维积分就可以了。高斯函数的积分为误差函数。例如我们计算宽度D内通过的光强占总光强的比例:
    -D2D2e-2x2w2dx-e-2x2w2dx=erf(2D2w)
    按照上面的例子,这个比值应该是80%。通过查表或者数值计算,我们可以得到
    D=1.2816w=2.5631σ
    这是TEM00高斯光束的结果。对于其他分布种类的光束,这个结果会有变化。
  5. M2因子
  6. 对于TEM00高斯光束来讲,光束半径在光轴方向(z方向)的变化可以表示为
    w2(z)=w02+(λπw0)2(z-z0)2
    其中z0表示TEM00高斯模式腰的位置。对于一般光束,我们已经提到过,光束发散程度一定超过TEM00高斯光束,所以可以引入一个M2因子,这样一般光束半径可以表示为
    w2(z)=w02+M4×(λπw0)2(z-z0)2
    上述各式均是在某单一横向方向(例如x),对于一般光束,x和y方向可能有不同的半径和M2(例如半导体激光),甚至各个方向都不相同,这时候就需要更多方向上的测量才能确定光束分布。M2可以作为光束质量的一个量度。M2越接近1,光束质量越好,发散越小,聚焦时形成的斑点越小。
  7. 刀锋法测量
  8. 上面的例子中,我们用10%和90%作为光束宽度的测量标准,那么可不可以用其他数值呢?另外,上面得到的D和W之间的关系只对于TEM00高斯光束适用,那么一般的光束半径如何测量呢?Siegman的一篇论文阐述了这个问题,结论是对于绝大多数激光光束,如果阈值选择在8.5%到11.6%之间,TEM00高斯光束的结果都可以作为很好的近似,10%无疑是一个很好的折衷。其实对于大多数激光光束而言,TEM00都是主要成分,这样的近似在大多数时候都是没什么问题的。 测量时,可以在不同的位置z多测几组,然后进行拟合,得到光束的M2和其他各个参数。

2008年10月28日星期二

让blogger支持MathML公式输出

我想在我的blog里显示公式怎么办?我可以把每个公式变成图片文件,然后把公式作为图片显示。我不需要截屏软件,OpenOffice可以很方便的把包含公式的文档转换为包含图片的HTML文件。但是这么做太笨了。更好的办法是用MathML在网页上显示公式。 MathMLXML的一个应用,是W3C推荐使用的网页数学公式显示方案。遗憾的是,目前不是所有的浏览器都支持MathMLFirefox本身包含对MathML的支持,但是你可能需要安装显示数学公式需要的字体。如果你使用Linux,你的发行版很可能包含(甚至已经安装)了该字体,请以mathml为关键词搜索以确定是否已经安装。如果你不幸使用Windows,请到这里按照说明自行安装字体。如果你更不幸使用IE浏览器,那只好自生自灭了。(开个玩笑,请安装MathPlayer插件)。 现在你可以到这里测试一下你的浏览器是不是可以正确显示MathML公式了。能显示,很好。可是,blogger并不能自动支持MathML的显示,因为blogger的网页扩展名都是html,浏览器不会自动按照xml来显示。这个问题可以这么解决:到这里下载一个小javascript文件,把这个文件直接放到你的blogger template文件的head里面,或者放到互联网上某个地方,在template的head里面加入下面一行:
<script type="text/javascript" src="javascript文件的链接"></script>
这样blogger就可以显示你的MathML公式了。 MathML语法非常复杂,要想直接用MathML写出一个公式并不容易。TtM是个把LaTeX公式转换成MathML代码的小工具。这下就方便多了。把转换好的MathML代码复制到你的blog就可以了。 下面是一个MathML显示公式的例子:
U(r)=A0W0W(z)exp[-ρ2W2(z)]exp[-jkz-jkρ22R(z)+jζ(z)]

2008年10月21日星期二

不用重启切换runlevel

Linux下telinit命令允许不重启的情况下切换runlevel。

2008年10月19日星期日

魔术组合键SysRq

很多用户不熟悉SysRq键到底有什么用,甚至很多用户根本没注意过电脑键盘上还有这么一个键存在。在Linux系统出现问题时,只要内核没有完全崩溃,SysRq就可以为我们提供一个诊断和控制内核的途径。SysRq键的使用方法是同时按下:
Alt + SysRq + CommandKey
有些键盘无法同时处理那么多键,例如一些笔记本的键盘,很可能还需要按下额外的Fn键。这时候可以按照下面的顺序:
按住Alt不放 => 按下并放开SysRq => 按下CommandKey => 放开所有键
Fedora默认是不开启SysRq功能的。如果想开启这项功能,可以修改/etc/sysctl.conf文件,把kernel.sysrq设置为1:
kernel.sysrq = 1
常用的CommandKey有如下这些:
CommandKey功能
m目前内存信息
t目前任务信息
s同步所有挂载文件系统
u重新挂载所有文件系统为只读
b立即重启(不sync或者unmount)
r键盘从raw模式切换到XLATE模式
kSecure Access Key杀死现在终端上所有程序
c故意让系统崩溃
e向除init外所有进程发出SIGTERM信号
i向除init外所有进程发出SIGKILL信号
如果syslogd仍然起作用,SysRq触发结果将会保存在/var/log/messages

2008年9月30日星期二

如何使用Maxima(四)

矩阵运算
  1. 输入矩阵
  2. 可以用entermatrix(m,n)输入一个m乘n矩阵。如果m=n,Maxima会自动询问是什么类型的矩阵,以减少输入的工作量。
    (%i1) m: entermatrix(2,2);
    
    Is the matrix  1. Diagonal  2. Symmetric  3. Antisymmetric  4. General
    Answer 1, 2, 3 or 4 :
    4;
    Row 1 Column 1:
    a;
    Row 1 Column 2:
    1;
    Row 2 Column 1:
    0;
    Row 2 Column 2:
    b;
    
    Matrix entered.
                                      [ a  1 ]
    (%o1)                             [      ]
                                      [ 0  b ]
    也可以用matrix以行为单位输入
    (%i2) matrix([a,b,c],[1,2,3],[d,e,f]);
                                      [ a  b  c ]
                                      [         ]
    (%o2)                             [ 1  2  3 ]
                                      [         ]
                                      [ d  e  f ]
  3. 矩阵的转置,逆和行列式
  4. (%i3) transpose(m);
                                       [ a  0 ]
    (%o3)                              [      ]
                                       [ 1  b ]
    (%i4) invert(m);
                                     [ 1     1  ]
                                     [ -  - --- ]
                                     [ a    a b ]
    (%o4)                            [          ]
                                     [      1   ]
                                     [ 0    -   ]
                                     [      b   ]
    (%i5) determinant(m);
    (%o5)                                 a b
  5. 矩阵的代数运算
  6. (%i6) m . %o4;
                                       [ 1  0 ]
    (%o6)                              [      ]
                                       [ 0  1 ]
    (%i7) %o6 + %o4;
                                   [ 1         1  ]
                                   [ - + 1  - --- ]
                                   [ a        a b ]
    (%o7)                          [              ]
                                   [        1     ]
                                   [   0    - + 1 ]
                                   [        b     ]
  7. 提取矩阵的某一行或某一列
  8. (%i8) row(m,1);
    (%o8)                              [ a  1 ]
    (%i9) col(%o7,2);
                                       [    1  ]
                                       [ - --- ]
                                       [   a b ]
    (%o9)                              [       ]
                                       [ 1     ]
                                       [ - + 1 ]
                                       [ b     ]
  9. 为矩阵增加一行(列)或几行(列)
  10. (%i10) addcol(m,[0,1],[c,1]);
                                     [ a  1  0  c ]
    (%o10)                           [            ]
                                     [ 0  b  1  1 ]
    (%i11) addrow(%o10,[1,0,0,1]);
                                    [ a  1  0  c ]
                                    [            ]
    (%o11)                          [ 0  b  1  1 ]
                                    [            ]
                                    [ 1  0  0  1 ]
  11. 另外,+, -, *, /, ^, exp(x)等运算都可以应用于矩阵,这些运算符分别作用于每一个元素。
  12. (%i12) m/c;
                                       [ a  1 ]
                                       [ -  - ]
                                       [ c  c ]
    (%o12)                             [      ]
                                       [    b ]
                                       [ 0  - ]
                                       [    c ]

2008年9月25日星期四

Vim输入特殊字符(例如希腊字母)

Digraphs可以用来在Vim中输入特殊字符。Digraphs使用两个一般字符来表示一个特殊字符,输入之前先按Ctrl+K,然后紧接着键入两个一般字符。 例如要键入拼音字母ü,可以先按Ctrl+K,然后键入u:。如果要键入希腊字母θ,可以先按Ctrl+K,然后键入h*。 这个双字符组合比较有规律,例如第二个字符
  • *表示希腊字母
  • 表示字母上面加两点
  • R表示罗马数字
等等。默认的组合列表可以在Vim里用:dig命令查看。另外,还可以手动添加字符。详细说明可以参考Vim Digraphs文档。

2008年9月22日星期一

改变Linux终端分辨率

grub.conf文件里kernel启动选项加入参数
vga=xxx
其中xxx是代表显示模式的数字,具体多少可以从这里查出来。下面是一些常用的数值:
640x480800x6001024x768
256 color palette769771773
16-bit785788791
24-bit786789792