之前讨论过FITS文件的投影问题( http://blog.sciencenet.cn/blog-117333-604768.html ),这是处理天文观测数据时常见的问题。除了投影,还有一个重要的问题,就是旋转。在观测的时候,探测器的竖直边并不一定总是保持在南北方向,有时会与南北方向有一定夹角,这个时候得到的数据除了有投影的问题,还有旋转的问题。 在没有投影的情况下,旋转的问题是简单的( http://blog.sciencenet.cn/home.php?mod=spaceuid=117333do=blogid=258897 ),因为球面上的旋转是非常明确的,有明确的表达式可以计算。于是,在有投影的情况下,通常需要先将投影过的数据去投影(变换为球面上的数据),在球面上进行旋转,然后再进行投影,这样就可以得到竖直边为南北方向的投影数据。 基本程序如下: theta=90.0-dec_p phi=ra_p tempx=sin(theta*convert)*cos(phi*convert) tempy=sin(theta*convert)*sin(phi*convert) tempz=cos(theta*convert) vect= vect=transpose(vect) theta0=90.0-crvaly phi0=crvalx matrix1= ,$ ,$ ] matrix1t=transpose(matrix1) matrix2= ,$ ,$ ] matrix2t=transpose(matrix2) matrix3= ,$ ,$ ] vect_f=matrix1t##matrix2t##matrix3##matrix2##matrix1##vect tempx=vect_f(0) tempy=vect_f(1) tempz=vect_f(2) if (abs(tempx) lt 1.0e-8) then begin if (tempy ge 0.0) then tempra=!pi/2.0 if (tempy lt 0.0) then tempra=3.0*!pi/2.0 endif else begin if (tempy gt 0.0) then begin if (tempx gt 0.0) then begin tempra=atan(tempy/tempx) endif else begin tempra=atan(tempy/tempx)+!pi endelse endif else begin if (tempx gt 0.0) then begin tempra=atan(tempy/tempx)+2.0*!pi endif else begin tempra=atan(tempy/tempx)+!pi endelse endelse endelse tempdec=asin(tempz) ra=tempra/convert dec=tempdec/convert
use warnings; use strict; open FH,"test.txt"; open RESULT,"result.txt"; open LINE,"line.txt"; open LEN, "length.txt"; while (FH) { #chomp; print RESULT "line $. is $_"; print LINE "line $.\n"; my $leng=length$_; print LEN "The length of line $. is $leng \n"; } close FH; test.txt as follow 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 3 3 3 3 3 3 3 3 3