65.9K
CodeProject 正在变化。 阅读更多。
Home

衍射光栅计算器

federico.strati
starIconstarIconstarIconstarIcon
emptyStarIcon
starIcon

4.40/5 (5投票s)

2010年10月27日

GPL3

8分钟阅读

viewsIcon

43085

downloadIcon

690

一个使用VTK库的与凹面衍射光栅相关的MFC Windows程序。

引言

这是关于一款与物理学相关的软件：光栅计算器程序的简短介绍。

光栅计算器程序是一款免费软件，根据GPL GNU许可证版本3.0在SourceForge上提供。本文介绍了该软件所解决的物理问题，以及程序如何在MFC中使用VTK（可视化工具包）库进行物理结果的3D图形表示。VTK库是建立在OpenGL图形基础库之上的、用于科学可视化的3D渲染引擎。所使用的VTK库版本是5.4.2，VTK主页。该程序还使用了Cedric Moonen的图表控件：图表控件 ver. 3.0，可在CodeProject上找到。

程序所解决的问题如下：凹面镜面（即反射式）衍射光栅用于研究中，以获得与衍射元件相关的聚焦特性。这种光栅反射的光在光谱色散的同时还被聚焦。这使得单个元件即可成为完整的光学仪器，将光衍射并聚焦到平面或曲面上进行记录。例如，弯曲的Rowland光栅是最简单的结合了光谱衍射和聚焦特性的光学元件。安装在所谓的Rowland光栅架（即装置）上的光栅具有将光谱色散光聚焦在圆上的特性，该圆的直径等于凹面光栅的半径。如果光源位于该圆上（称为Rowland圆），则衍射后光会聚焦到同一圆上。

最近，物理学家开始使用所谓的“平场光谱仪”。这类光谱仪使用凹面反射光栅，其线密度并非严格周期性的。这类光栅的刻槽周期性从一端到另一端以非周期方式变化。这类光栅的刻划是通过数值控制来实现的，以非常精确地遵循预设的模式。叠加在主要周期性刻划之上的这种非周期性使得这类光栅能够将光谱色散光聚焦到大致平坦的表面上，即平面。这很有用，因为目前大多数可用的记录传感器都有平面（例如，与胶片相对的CCD - 电荷耦合器件）。

本程序能够求解通用光栅的光线追迹方程，给出所求范围内每种色散光波长的聚焦位置。光栅通过其非周期性刻划来表征，可以从简单的Rowland光栅架（无非周期性）到最现代的平场光栅（受控非周期性刻划）。该程序已在装置参数中识别出两种在现代软X射线光谱仪研究中使用的非周期性刻划光栅。这些参数源自参考文献[1,2,3]的研究，其中介绍了当前日期的掠入射平场光谱仪的概念（另见参考文献[5]）。

免责声明：本文档对程序所解决的问题进行了深入的物理描述，并且唯一提供的代码是关于在显示科学数据时使用VTK 3D库的内容。

物理背景

本节将从物理学的角度说明程序所要解决的物理问题。

软X射线领域的平场光谱仪是成熟的光谱学研究技术。提供增强光谱图像聚焦特性的凹面和平面变线距（VLS）光栅已广泛应用于同步辐射仪器[1]、等离子体诊断[2]和空间天体物理学[3]。在设计凹面光栅及其光谱图像聚焦特性的分析中，通常采用一种使用光程函数的解析方法[4]。

如果P点是光栅上的任意一点，具有其表面的局部坐标u、w、l（其中u是w、l的函数），A和B分别是波长为λ、在m阶被光栅衍射的光谱的源点和名义焦点，则光程函数定义为

F = AP + BP + n*m*lambda

其中n是P点处从光栅上参考点O（我们可以假定它是光栅毛坯的中心）开始计算的刻槽数。

GratingComp.GIF

计算的光栅装置

根据费马原理，为了使图中的B点成为A点在波长λ下的像，要求以下微分条件同时满足

dF/dw = 0, dF/dl = 0.

再次参考图示，我们引入极坐标系并表示距离为

AP = r2 + u2 + l2 - 2*u*r*cos(alpha) - 2*w*r*sin(alpha)
BP = r'2 + u2 + l2 - 2*u*r'*cos(beta) - 2*w*r'*sin(beta)

其中α是入射角，β是从O处出射到光栅的衍射角。

对于具有球形基底的凹面光栅，点P由方程定义

(u-R)2 + w2 + l2 = R2

其中R是光栅的半径。将其展开为幂级数，得到

u = (w2+l2)/2R + (w2+l2)2/8R3 + ...

变线距光栅（非周期性刻划）上的刻槽数n定义为

n = (rho)0 * ( w + b2/R * w2 + b3/R2 * w3 + b4/R3 * w4 + ... )

其中（ρ）0是光栅中心O处的线密度（因此它是线密度的周期性分量），而b2、b3和b4是控制光栅刻划非周期部分的数值参数。

利用上述方程，我们可以将光程函数表示为w和l的幂级数，如下所示

F = r + r' + w*F10 + w2*F20 + l2*F02 
        + w3*F30 + w*l2*F12 + w4*F40 
        + w2*l2*F22 + l4*F04 + O(w5).

在此方程中，每个Fij项都可以表示为

Fij = Cij + m*lambda*rho0*Mij

关于系数Cij和Mij的完整列表，请参见参考文献[3]。

当为成像光谱仪设计非周期性刻划光栅时，光学布局和刻划参数的选择是为了最小化光谱方向和空间方向上的像差。像点成像要求Fij=0直到尽可能高的阶数。

F10=0的条件给出光栅方程

sin(alpha) + sin(beta) = m*lambda*rho0

而沿着光谱色散方向的焦距曲线则由F20=0得到，并与光栅方程联用

r' = r * R * cos2(beta) / 
    ( r*(cos(alpha)+cos(beta)-2*m*lambda*rho0*b2) - R*cos2(alpha))

这表明，给定任何光学布局，只有线性刻划参数b2决定了焦平面的位置。另外两个刻划参数b3、b4通常用于在给定感兴趣的波长下最小化光谱方向上的彗差和球差（分别为F30=0、F40=0）。也可以研究特定的装置来减小空间方向上的像差[3]。

为了估算这些光谱仪的光谱聚焦特性，可以使用光线追迹代码代替Fij=0的方程。

代码通过追迹大量光线穿过光栅并跟踪它们直到找到最小散布来确定聚焦位置。对于反射面的像点成像，光线应在从源发出相同传播时间后收敛到焦点。对于光栅也是如此，如果考虑其元件引入的用于相干干涉的光程，则也是如此。这是代码的核心：对于每条光线，光程会通过一个量进行校正，该量恰好对应于光栅在色散光波长中引入的衍射。对于m阶衍射，此量对应于

m*lambda*(Number of ruled lines between the point P on grating and reference point O)

即，恰好是定义为的光程函数右侧的项

F = AP + BP + n*m*lambda

由于像点不是像点的，因此焦点周围光线的最小几何散布提供了光谱分辨率的估计。这种几何光学预测的分辨率当然不能超过波动光学极限。代码已经过分析公式的测试，在预测焦平面位置方面与之吻合良好。

代码中的VTK渲染

程序中包含的对CodeProject读者最有价值的软件部分是用于将程序与3D渲染库VTK对接的类。该类在此处重新提供，并附带有关其所有成员的注释，描述了VTK 3D图形库的使用。在MFC应用程序中使用VTK时，我们从类CVTKView派生一个CView类型的对象，该类本身派生自CView。这允许VTK渲染引擎无缝集成到MFC视图中。

类的头文件如下

// CVTKView view

#if !defined( __VTK_VIEW__ )
#define __VTK_VIEW__

// VTK headers to include
#include "vtkActor.h"
#include "vtkCamera.h"
#include "vtkProperty.h"
#include "vtkProperty2D.h"
#include "vtkPropCollection.h"
#include "vtkRenderer.h"
#include "vtkRenderWindow.h"
#include "vtkRenderWindowInteractor.h"
#include "vtkTextProperty.h"

#include "vtkMFCWindow.h"
#include "vtkWin32OpenGLRenderWindow.h"

class CVTKView : public CView
{
    DECLARE_DYNCREATE(CVTKView)
protected:
    CVTKView();   // protected constructor used by dynamic creation
    virtual ~CVTKView();

protected:
   vtkRenderer*   pvtkRenderer;
   vtkMFCWindow*  pvtkMFCWindow;
   bool m_bInitCreate; // to check to avoid duplicating work on New/Open Documents
   bool m_bInitUpdate; // to check to avoid duplicating work on New/Open Documents

public:
   vtkRenderer*  GetRenderer()  { ASSERT(pvtkRenderer); return pvtkRenderer; }
   vtkMFCWindow* GetMFCWindow() { ASSERT(pvtkMFCWindow); return pvtkMFCWindow; };
   // Render the view
   void Render();
   // Override to avoid duplicating work OnCreate on New/Open Documents
   // Call base class first !, use m_bInitCreate
   virtual void OnViewCreate();
   // Override to avoid duplicating work OnInitialUpdate on New/Open Documents
   // Call base class first !, use m_bInitUpdate
   virtual void OnViewInitUpdate();

public:
    virtual void OnDraw(CDC* pDC);      // overridden to draw this view
#ifdef _DEBUG
    virtual void AssertValid() const;
#ifndef _WIN32_WCE
    virtual void Dump(CDumpContext& dc) const;
#endif
#endif

protected:
    DECLARE_MESSAGE_MAP()
public:
   afx_msg int OnCreate(LPCREATESTRUCT lpCreateStruct);
   afx_msg void OnDestroy();
   afx_msg BOOL OnEraseBkgnd(CDC* pDC);
   afx_msg void OnSize(UINT nType, int cx, int cy);
   virtual void OnInitialUpdate();
};

#endif // !defined( __VTK_VIEW__ )

类的实现部分如下

// VTKView.cpp : implementation file

#include "stdafx.h"
#include "VTKView.h"

// CVTKView

IMPLEMENT_DYNCREATE(CVTKView, CView)

CVTKView::CVTKView()
{
   pvtkMFCWindow = NULL;
   pvtkRenderer  = NULL;

   m_bInitCreate = false;
   m_bInitUpdate = false;
}

CVTKView::~CVTKView()
{
}

BEGIN_MESSAGE_MAP(CVTKView, CView)
   ON_WM_CREATE()
   ON_WM_DESTROY()
   ON_WM_ERASEBKGND()
   ON_WM_SIZE()
END_MESSAGE_MAP()

// CVTKView drawing

void CVTKView::OnDraw(CDC* pDC)
{
    CDocument* pDoc = GetDocument();
   if (pvtkMFCWindow)
   {
      // VTK printing support
      if (pDC->IsPrinting())
         pvtkMFCWindow->DrawDC(pDC);
   }
}

// CVTKView diagnostics

#ifdef _DEBUG
void CVTKView::AssertValid() const
{
    CView::AssertValid();
}

#ifndef _WIN32_WCE
void CVTKView::Dump(CDumpContext& dc) const
{
    CView::Dump(dc);
}
#endif
#endif //_DEBUG

// CVTKView message handlers

int CVTKView::OnCreate(LPCREATESTRUCT lpCreateStruct)
{
   if (CView::OnCreate(lpCreateStruct) == -1)
      return -1;

   // call our virtual create view method
   OnViewCreate();

   m_bInitCreate = true;

   return 0;
}

/*
 * Override to avoid duplicating work OnCreate on New/Open Documents
 * Call base class first !, use m_bInitCreate
 */
void CVTKView::OnViewCreate()
{
   if ( !m_bInitCreate )
   {
      if (pvtkMFCWindow)
      {
         delete pvtkMFCWindow;
         pvtkMFCWindow = NULL;
      }
      // Create the new hosted MFC window
      pvtkMFCWindow = new vtkMFCWindow(this);
      // Create the renderer and add it to window's view
      pvtkRenderer  = vtkRenderer::New();
      if (pvtkMFCWindow)
         pvtkMFCWindow->GetRenderWindow()->AddRenderer(pvtkRenderer);
   }
}

void CVTKView::OnDestroy()
{
   CView::OnDestroy();

   if (pvtkMFCWindow)
   {
      delete pvtkMFCWindow;
      pvtkMFCWindow = NULL;
   }
   if (pvtkRenderer)
   {
      pvtkRenderer->Delete();
   }
}

// Avoid possible flickering
BOOL CVTKView::OnEraseBkgnd(CDC* pDC)
{
   // return CView::OnEraseBkgnd(pDC);
   return TRUE;
}

void CVTKView::OnSize(UINT nType, int cx, int cy)
{
   CView::OnSize(nType, cx, cy);

   // resize our view
   if (pvtkMFCWindow)
      pvtkMFCWindow->MoveWindow(0, 0, cx, cy);
}

void CVTKView::OnInitialUpdate()
{
   CView::OnInitialUpdate();

   // Call our virtual method for initial update
   OnViewInitUpdate();
   Render();

   m_bInitUpdate = true;
}

/*
 * Override to avoid duplicating work OnInitialUpdate on New/Open Documents
 * Call base class first !, use m_bInitUpdate
 */
void CVTKView::OnViewInitUpdate()
{
   if ( !m_bInitUpdate )
   {
      // Add a default interactor to our view
      vtkWin32OpenGLRenderWindow* pRenderWindow = NULL;
      pRenderWindow = pvtkMFCWindow->GetRenderWindow();
      if (pRenderWindow)
      {
         pRenderWindow->GetInteractor()->Initialize();
         pRenderWindow->GetInteractor()->Start();
      }
   }
}

/*
 * Render the view
 */
void CVTKView::Render()
{
   if (!pvtkMFCWindow)
      return;

   vtkWin32OpenGLRenderWindow* pRenderWindow = NULL;
   pRenderWindow = pvtkMFCWindow->GetRenderWindow();
   if (pRenderWindow)
   {
      pvtkRenderer->Render();
      pRenderWindow->Render();
   }
}

VTK围绕图形管线的概念工作。数据被处理，从源对象开始，然后进行过滤（即处理），最后映射到视图中的一个Actor。

Source objects -> Filters -> Mappers -> Actors -> Insertion of actors in the Scene

在从CVTKView派生的类中，使用一个单独的方法添加用于在VTK派生视图中显示我们数据的图形管线

void CPlotFocalDensity::PlotData()
{
   if ( !m_pFocalSpotObj )
      return;

   CWaitCursor wc;

   int nX, nY; nX = nY = 0;
   m_pFocalSpotObj->GetAttribute(FS_DENSITY_3D_NX)->Get_int(nX);
   m_pFocalSpotObj->GetAttribute(FS_DENSITY_3D_NY)->Get_int(nY);

   CAttribute_Generic *at, *atX, *atY;
   at = atX = atY = NULL;
   at  = m_pFocalSpotObj->GetAttribute(FS_DENSITY_3D);
   atX = m_pFocalSpotObj->GetAttribute(FS_DENSITY_3D_X);
   atY = m_pFocalSpotObj->GetAttribute(FS_DENSITY_3D_Y);

   size_t Naux;
   at->GetSize(Naux);

   // Set-up VTK Source objects
   vtkStructuredGrid* sgrid  = vtkStructuredGrid::New();
   vtkFloatArray*     scalar = vtkFloatArray::New();
   vtkPoints*         points = vtkPoints::New();

   double x, y, z;
   for (size_t j=0; j<Naux; j++)
   {
      atX->GetAt(x,j);
      atY->GetAt(y,j);
      points->InsertPoint(j,x,y,0.0);
       at->GetAt(z,j);
      scalar->InsertValue(j,z);
   }

   sgrid->SetDimensions(nX,nY,1);
   sgrid->SetPoints(points);
   sgrid->GetPointData()->SetScalars(scalar);

   // Create a color lookup table
   vtkLookupTable* lut = vtkLookupTable::New();
   lut->SetNumberOfTableValues(256);
   lut->SetHueRange(0.0,0.667);
   lut->Build();

   // Map the grid datasets
   vtkDataSetMapper*  mapper = vtkDataSetMapper::New();
   mapper->SetInput(sgrid);
   mapper->SetLookupTable(lut);
   double tmp[2];
   sgrid->GetScalarRange(tmp);
   mapper->SetScalarRange(tmp);

   // Apply a contour filter to the data
   vtkContourFilter* contourFilter = vtkContourFilter::New();
   contourFilter->SetInput(sgrid);
   contourFilter->GenerateValues(20,tmp);

   // Map the contour field
   vtkPolyDataMapper* contourMapper = vtkPolyDataMapper::New();
   contourMapper->SetInput(contourFilter->GetOutput());
   sgrid->GetScalarRange(tmp);
   contourMapper->SetScalarRange(tmp);

   // Add the actors in the scene
   if ( carpet )
   {
      pvtkRenderer->RemoveActor(carpet);
      Render();
   }
   carpet = vtkActor::New();
   carpet->SetMapper(mapper);

   if ( contour )
   {
      pvtkRenderer->RemoveActor(contour);
      Render();
   }
   contour = vtkActor::New();
   contour->SetMapper(contourMapper);

   if ( scalarBar )
   {
      pvtkRenderer->RemoveActor(scalarBar);
      Render();
   }
   scalarBar = vtkScalarBarActor::New();
   scalarBar->SetLookupTable(lut);
   scalarBar->SetTitle(_T("Focal Spot Density"));
   scalarBar->GetPositionCoordinate()->SetCoordinateSystemToNormalizedViewport();
   scalarBar->GetPositionCoordinate()->SetValue(0.8,0.2);

   // assign our actors to the renderer
   pvtkRenderer->AddActor(carpet);
   pvtkRenderer->AddActor(contour);
   pvtkRenderer->AddActor(scalarBar);
   pvtkRenderer->SetBackground(1.0,1.0,0.4);

   // draw the resulting scene
   pvtkRenderer->ResetCamera();
   pvtkRenderer->GetActiveCamera()->Zoom(3.0);
   pvtkRenderer->GetActiveCamera()->Elevation(45);
   pvtkRenderer->GetActiveCamera()->Azimuth(30);
   pvtkRenderer->GetActiveCamera()->Dolly(1.0);
   pvtkRenderer->ResetCameraClippingRange();
   Render();
}

关注点

在MFC Visual Studio 2008 SP1构建的应用程序中使用VTK 3D库的示例。

参考文献

  • [1] T. Harada, T. Kita, M. Itou, H. Taira, and A. Mikuni, Nucl. Instrum. Methods A 246, 272 (1986). M. Itou, T. Harada, and T. Kita, Appl. Opt. 28, 146 (1989)。
  • [2] T. Kita, T. Harada, N. Nakano, and H. Kuroda, Appl. Opt. 22, 512 (1983). N. Nakano, H. Kuroda, T. Kita, and T. Harada, Appl. Opt. 23, 2386 (1984)。
  • [3] T. Harada, H. Sakuma, K. Takahashi, T. Watanabe, H. Hara, and T. Kita, Appl. Opt. 37, 6803 (1998)。
  • [4] T. Namioka, J. Opt. Soc. Am. 49, 446 (1959)。
  • [5] Instruments division, Hitachi Ltd, 882 Ichige, Hitachinaka-shi, Ibaragi 312− 0033, Japan。

历史

  • 2010年10月27日 - 初稿。
© . All rights reserved.