本文共 2410 字，大约阅读时间需要 8 分钟。

IDL制作卫星天顶角和方位角数据集，对应NC文件中的SAA、SAZ

PRO temp  COMPILE_OPT idl2  e=envi(/headless)    locRaster=e.OpenRaster('**\loc.dat')  locData=locRaster.GetData()  pi=3.14159265358979323846  d2r=0.01745329252;pi/180  r2d=57.295779513;180/pi  cols=(locdata.dim)[0]  lines=(locdata.dim)[1]  outdata=make_array(cols,lines,2,/FLOAT,VALUE=0)  lon=locData[*,*,0]  lat=locData[*,*,1]  alt=locData[*,*,2]#单位：km  ;计算卫星天顶角、方位角  satz=obslook_sat(lon, lat, alt)  ;输出角度数据集  outraster=ENVIRaster(outdata,URI='*\satz.dat')  outraster.Save  locRaster.CLose  outraster.Close     e.CloseEND

代码中的’**\loc.dat’是影像的像元对应的经纬度数据，可以查看https://blog.csdn.net/qq_33339770/article/details/102957857

角度数据集生存代码参看：https://download.csdn.net/download/urbancorbie/10343539?ops_request_misc=%257B%2522request%255Fid%2522%253A%2522158769212919726869052219%2522%252C%2522scm%2522%253A%252220140713.130102334.pc%255Fall.57649%2522%257D&request_id=158769212919726869052219&biz_id=1&utm_source=distribute.pc_search_result.none-task-download-2_allfirst_rank_v2~rank_v25-1

;卫星天顶角、方位角Function obslook_sat,lon,lat,alt  COMPILE_OPT idl2  pi=3.14159265358979323846  d2r=0.01745329252;pi/180  r2d=57.295779513;180/pi   cols=(lon.dim)[0]  rows=(lon.dim)[1]    clat=atan(0.993305616*tan(lat*d2r));'Rpol^2 / Req^2'=0.993305616  clon=make_array(cols,rows,2,/FLOAT,value=d2r*lon);clon=d2r*lon  slon=make_array(cols,rows,2,/FLOAT,value=d2r*140.7);'sub_lon'=140.7,slon=d2r*140.7  ;计算卫星与地面距离  ;  'Earth's polar radius'=6356.7523,'(Req^2 / Rpol^2) / Req^2'=0.00669438444  rp=6356.7523/(sqrt(1-0.00669438444*cos(clat)*cos(clat)))+alt  ;计算地面坐标  x0=rp*cos(clat)*cos(clon)  y0=rp*cos(clat)*sin(clon)  z0=rp*sin(clat)  ;计算卫星坐标  ;'Distance from Earth's center to virtual satellite'=42164km  xs=42164*cos(slon)  ys=42164*sin(slon)  zs=make_array(cols,rows,2,/FLOAT,value=0);zs=0  ;计算差向量  rx=xs-x0  ry=ys-y0  rz=zs-z0  ;计算地心坐标系下的极坐标形式的差向量  ;南方,极坐标下：x=r*sin(el)*cos(az)  ls=sin(lat*d2r)*cos(clon)*rx+sin(lat*d2r)*sin(clon)*ry-cos(lat*d2r)*rz  ;西方,极坐标下：y=r*sin(el)*sin(az),计算的是东方  lse=-sin(clon)*rx+cos(clon)*ry  ;天顶,极坐标下：z=r*cos(el)  lz=cos(lat*d2r)*cos(clon)*rx+cos(lat*d2r)*sin(clon)*ry+sin(lat*d2r)*rz  ;计算方位角：az=atan(y/x)  az=atan(-lse/ls)  az[where(ls gt 0)] += pi  az[where(az lt 0)] += pi*2  ;计算高度角  rg=sqrt(rx*rx+ry*ry+rz*rz)  el=asin(lz/rg)  el=el*r2d  az=az*r2d  ;方位角修正为[-180.0,180.0]  az[where(az GT 180.0)] -= 360.0  ;天顶角为高度角的余角  return,[az,90.0-el]END

转载地址：https://blog.csdn.net/qq_33339770/article/details/103023725 如侵犯您的版权，请留言回复原文章的地址，我们会给您删除此文章，给您带来不便请您谅解！