数据介绍：Lunar Reconnaissance Orbiter Camera (LROC) Narrow Angle Cameras (NACs)——月球轨道窄角相机，是月球勘测轨道器(LRO)上的一对高分辨率推扫式(pushbroom)相机。两台相机（NAC-Left，NAC-Right）略有重叠地平行扫描月面，并同时获取图像。典型的LROC NAC图像的像素分辨率在0.5-2.0 m/px之间，部分图像的空间分辨率可达0.25-0.3 m/px，是目前分辨率最高的月球影像数据集。一个典型的NAC成对图像大约有10,000个像素宽（东西向）和52,000个像素长（南北向），对应的正常覆盖范围大小为5x26公里到20x100公里左右。NAC数据集由NASA's Planetary Data System（PDS）提供，可供公众下载。

数据下载网址：

LROC Image Search (asu.edu)
LROC | PDS Archive (asu.edu)
pds-imaging.jpl.nasa.gov （需要能够访问外网）
QuickMap
Lunar Orbital Data Explorer - Home Page (wustl.edu)

注意：获取的数据为.IMG格式，为了读取和提取图像的地理信息，需要将其转换为(.geotiff or .gtif)格式，该数据一般包括：EDR/CDR/RDR，我们下载EDR数据，CDR数据的处理比EDR更复杂。

EDR：未经过处理的原始NAC数据；
CDR：经过光度校正但无地理投影的NAC数据；
RDR：经过几何校正和地理投影的NAC数据

打开第一个网址后如图所示，可根据NAC的（轨道号、相机角度、入射角和发射角、分辨率和图像ID）查询数据，Product ID（eg.M1172873803LE，M1172873803RE）下载对应的EDR数据。

2.利用ISIS处理LROC-NAC(EDR)数据

我们使用美国地质调查局United States Geological Survey (USGS)开发的成像仪和光谱仪集成软件Integrated Software for Imagers and Spectrometers (ISIS)处理数据。

ISIS的下载安装参考：深空探测处理软件USGS-ISIS的下载安装_qq_31721595的博客-CSDN博客

2.1 将NAC-EDR转换为ISIS的cube格式

（1）首先找到你所下载的数据文件的位置，然后打开终端，利用conda进入ISIS，如图：

（2）将M1172873803LE.IMG，M1172873803RE.IMG转换为ISIS的cube格式，命令如下：

lronac2isis from=M1172873803LE.IMG to=M1172873803L.raw.cub

2.2 附加相机位置信息

ISIS使用Navigation and Ancillary Information Facility SPICE system在观测过程中确定相机位置，以便将图像正确投影到地面上。该步只要在光度校正前执行即可，开始时执行此步骤便于立即查询图像的几何信息，或在进行额外处理前确定是否有所需要的kernels可用。命令如下：

spiceinit from=M1172873803L.raw.cub spksmithed=true web=true

spiceinit from=M1172873803R.raw.cub spksmithed=true web=true

对M1172873803R做同样处理，如图：

2.3 辐射校正

lronaccal from=M1172873803L.raw.cub to=M1172873803L.cal.cub

lronaccal from=M1172873803R.raw.cub to=M1172873803R.cal.cub

2.4 回波校正

lronacecho from=M1172873803L.cal.cub to=M1172873803L.echo.cub

lronacecho from=M1172873803R.cal.cub to=M1172873803R.echo.cub

注：回波校正后会生成伪影，利用trim去除。

如果对图像求和(在EDR标签中，查找“CROSSTRACK_SUMMING = 2”，或者检查图像宽度是否为2532像素)，将这些数字除以2:修剪量应该是23和13 ：

NAC-L：

trim from=M1172873803L.echo.cub to=M1172873803L.tr.cub left=46 right=26

NAC-R：

trim from=M1172873803R.echo.cub to=M1172873803R.tr.cub left=26 right=46

2.5 光度校正（可选）

这一步不是必须的

photomet from=M1172873803L.tr.cub to=M1172873803L.pho.cub frompvl=basicpho.pvl

photomet from=M1172873803R.tr.cub to=M1172873803R.pho.cub frompvl=basicpho.pvl

不过，若要进行光度校正首先需要新建一个basicpho.pvl文件，该文件内容为：

Object = NormalizationModel
	Group = Algorithm
 		Name = Mixed
 		Incmat=87.0
 		Albedo=0.05
 		Incref=0.0
 		Thresh=10E30
 	EndGroup
EndObject
Object = PhotometricModel
 	Group = Algorithm
 		Name = HapkeHen
 		Wh=0.278906221
 		Hh=0.0438
 		B0=1.0
 		Theta=10.0
 		HG1=-0.229315607
 		HG2=0.0
 	EndGroup
EndObject

2.6 创建一个地图文件mosaic.map

首先将所有想要组合成一个mosaic的图像的路径放入自建的文本文件"mosaic_level_1.lis"中，如图：

然后，利用mosrange命令生成地图文件。此外，如果知道感兴趣区域的确切边界或中心点，还可用maptemplate，但这里用的mosrange命令，使用的是等矩形投影equirectangular，若纬度在70°以上或靠近极点的区域，利用polar stereographic。

注：因为所有发布的LRO产品使用的都是0-360的经度值，而ESRI程序（如ArcMap）经度值在-180-180，所以若要在ArcMap中使用地图投影产品，还需在后面添加命令，londom=180。

mosrange fromlist=mosaic_level_1.lis to=mosaic.map projection=equirectangular londom=180

2.7 地图投影

cam2map from=M1172873803L.pho.cub to=M1172873803L.map.cub map=mosaic.map warpalgorithm=forwardpatch patchsize=50

cam2map from=M1172873803R.pho.cub to=M1172873803R.map.cub map=mosaic.map warpalgorithm=forwardpatch patchsize=50

2.8 镶嵌

ls M1172873803?.map.cub > M1172873803.lis

automos fromlist=M1172873803.lis mos=M1172873803.mosaic.cub

运行完第一句命令后，生成M1172873803.lis文件，如图：

然后运行第二条命令，但是报错了，目前还不清楚是什么原因，如图：

但我们可以将得到的 M1172873803L.map.cub 和 M1172873803R.map.cub 放在ArcMap中进行镶嵌，见第3节。

2.9 配准（可选）

有时一对NAC图像可能存在明显的偏移，需要使用corge来对齐两个图像。

map2map from=M1172873803R.map.cub to=M1172873803R.Lbounds.cub map=M1172873803L.map.cub matchmap=true

下一步，首先需要创建nac_pair_align.def文件，内容如下：

Object = AutoRegistration
 	Group = Algorithm
		Name = MaximumCorrelation
 		Tolerance = 0.7
 		SubpixelAccuracy = True
 	End_Group
 	
	Group = PatternChip
 		Samples = 20
 		Lines = 20
 	End_Group

# If you determine that the vertical or horizontal offset is >50 pixels, adjust the SearchChip 
# size to be at least `20+offset*2`, where 20 is the PatternChip size.
 
	Group = SearchChip
 		Samples = 120
 		Lines = 120
 	End_Group
 
	Group = SurfaceModel
 		DistanceTolerance = 1.5
 		WindowSize = 5
 	End_Group
End_Object
End

然后，运行如下命令，会得到 SampleAverage 和 LineAverage，用于下一步translate命令：

coreg from=M1172873803R.Lbounds.cub match=M1172873803L.map.cub deffile=nac_pair_align.def

利用translate命令纠正偏移，将$sampleaverage 和 $lineaverage换为上图corge输出的值，-0.58和-1.13.

translate from=M1172873803R.map.cub to=M1172873803R.map.shifted.cub strans=$sampleaverage ltrans=$lineaverage

然后再对M1172873803L.map.cub、M1172873803R.map.shifted.cub运行2.8中的镶嵌命令，但依然报错，所以我直接跳过2.8、2.9的命令，利用ArcMap拼接影像。

最后需要的就是这两个数据文件：

3.利用GDAL转换为地理参考格式（.tif）

因为未在虚拟机中安装GDAL，我在本机电脑的Pycharm中利用Python进行转换，前提是也需要安装gdal包，或从osgeo中导入gdal包。

from osgeo import gdal

cubdata = r"D:\Moon_Data\processedData\M1096293859L.map.cub"   # 经上面步骤得到的.cub文件路径
tifdata = r"D:\Moon_Data\processedData\M1096293859L.tif"       # .tif文件路径

gdal.Translate(tifdata, cubdata)

最终得到下图的.tif文件，并在ArcMap中打开它们，可以看到这些.tif图像带有投影信息了：