1 / 50

An Introduction to IDL (The Interactive Data Language) and IDL in Astronomy

Astronomy Seminar Series 11/12/2005. An Introduction to IDL (The Interactive Data Language) and IDL in Astronomy. 1. IDL Features 2. Basics of IDL 3. 1D & 2D Display 4. FITS I/O in IDL 5. Image Processing. cqras@126.com. I.IDL Features.

aaralyn
Télécharger la présentation

An Introduction to IDL (The Interactive Data Language) and IDL in Astronomy

An Image/Link below is provided (as is) to download presentation Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author. Content is provided to you AS IS for your information and personal use only. Download presentation by click this link. While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server. During download, if you can't get a presentation, the file might be deleted by the publisher.

E N D

Presentation Transcript

  1. Astronomy Seminar Series 11/12/2005 An Introduction to IDL (The Interactive Data Language)and IDL in Astronomy 1. IDL Features 2. Basics of IDL 3. 1D & 2D Display 4. FITS I/O in IDL 5. Image Processing cqras@126.com

  2. I.IDL Features http://www.astro.virginia.edu/class/oconnell/astr511/IDLguide.html

  3. IDL vs. Mathematica, Matlab, Maple http://amath.colorado.edu/computing/mmm/brief.html

  4. IDL vs. Traditional Astronomical Software http://www.astro.virginia.edu/class/oconnell/astr511/IDLguide.html

  5. IDL in Astronomy http://idlastro.gsfc.nasa.gov/other_url.html

  6. Very Good Places and Packages IDL Astronomy User's Library http://idlastro.gsfc.nasa.gov/homepage.html http://idlastro.gsfc.nasa.gov/ftp/astron.tar.gz FITShttp://idlastro.gsfc.nasa.gov/fitsio.html Solar Software http://lmsal.com/solarsoft/sswdoc/index_menu.html ftp://sohoftp.nascom.nasa.gov/solarsoft/offline/swmaint/tar/ ssw_ssw_gen.tar.Z Coyote’s Guide to IDL Programming http://www.dfanning.com ftp://ftp.dfanning.com/pub/dfanning/outgoing/coyote2nd/ IDL Newsgroup (comp.lang.idl-pvwave) http://groups.google.com/group/comp.lang.idl-pvwave Markwardt IDL Library (Fitting) http://cow.physics.wisc.edu/~craigm/idl/idl.html JHUAPL IDL Library http://fermi.jhuapl.edu/s1r/idl/s1rlib/local_idl.html IDL + EMACS http://www.idlwave.org/ http://idlwave.org/download/idlwave-help.tar.bz2 IDL Astro Package astron.dir.tar.gz Coyote Program Library coyote2ndfiles.zip SSW General Package ssw_ssw_gen.tar.Z tar xvfz ssw_ssw_gen.tar.Z mv gen /usr/local/rsi/idl/lib/sswgen Others in a similar way; or see Page 8 Important: Practice, Take notes, Google idl + keywords

  7. compile run Working with IDL Menu Toolbars Demo Edit Window Output Log Variable Watch My favorite Linux/Unix Console Command Input Status Bar

  8. Personal SETUP for IDL(recommended) Flashing Colorsin Linux/Unix #edit ~/.Xresources idl.gr_visual: TrueColor   idl.gr_depth: 24   idl.retain: 2   idl.colors: -1 Personal IDL_STARTUP 1.edit ~/.cshrc : setenv IDL_STARTUP ~/.idl_startup.pro 2.edit ~/.idl_startup.pro (or other file names), Device, retain = 2, decomposed = 0 !path = !path + ‘:~/idlpros1:~/idlpros2’ #Windows: File  Preferences  Startup... you can select any file that is similar to that under Linux. Working Directory IDL> cd, ‘myworkdir’ Windows: File  Preferences …  Startup  Startup Memory Limitsin Hubble (Alpha) # edit~/.cshrc limit stacksize unlimited limit datasize unlimited

  9. II. Basics of IDL IDL>PRINT, 3 * 5, [30,5,50] 15 30 5 50 IDL>x = 'Hello! IDL World' & HELP, x X STRING = 'Hello! IDL World' IDL>x = indgen(15) & y = sin(2*!dpi*x/15) X INT = Array[10] Y DOUBLE = Array[10] IDL>FOR i = 0, 15-1 DO PRINT, i, x[i], y[i] IDL>DEVICE, decomposed = 0 IDL>plot, loaddata(1), psym=-4, $ title='plot', xtitle='Month', ytit='Sth‘ ;,$ ;/ylog, yrange=[5e-1,40], ystyle=1 Dynamic Datatype Array Operation Array Zero-Ordered Parameters Direct Graphics

  10. Some Symbols, Definitions and Others ? – online help .run, .compile, .r ; & $ ! ↑↓ http://fermi.jhuapl.edu/s1r/idl/idl_syntx.html

  11. IDL Variables: Dynamic • Scalar, Array (1—8D) • Structure: collection of scalars, arrays, or other structures • System Variables (!) • !dpi (3.1415926) • !p: Display. e.g., !p.font, !p.color • !d: Device. e.g. !d.name • HELP, variable or HELP, variable, /struct • Data type: dynamic

  12. Type Len Creation Array Conversion Byte 1 A=5B Bytarr Byte Integer 2 B=0;b=0S Intarr Fix Uint 2 C=0U Uintarr unit Long 4 D=0L Lonarr Long Ulong 4 E=0UL Ulonarr Ulong Long64 8 F=0LL Long64arr Long64 Ulong64 8 G=0ULL Ulon64arr Ulong64 Float 4 H=0.0 Fltarr Float Double 8 I=0.0D Dblarr Double Complex 8 J=complex(1.0,0.0) Complexarr Complex Dcomplex 16 K=dcomplex(1.0,0.0) Dcomplexarr Dcomplex String ? L=’hello’ Strarr String Pointer 4 M=ptr_new() Ptrarr --- Object 4 N=obj_new() Objarr --- IDL Variables: Basic Datatypes Indexed Array Creator: e.g., findgen()float index generator Table c.f. intenet

  13. Control Statements • IF: conditional if exp then statem if exp then statem1 else statem2 • For Loops: for i = init, limit, step do statem • While Loops: while exp do statem • Repeat Loops: repeat statem until exp • Case: • GOTO: goto, label • Blocks: Begin statem1 …… statemx Endxxx if x lt 0 then begin print,x a=2 endif for i=0, 10 do begin readf, lun, txt print,txt endfor e.g. if x lt 0 then begin print,x & a=2 & endif http://fermi.jhuapl.edu/s1r/idl/idl_syntx.html

  14. Programs: Procedures & Functions IDL: an interactive tool, also a powerful programming language. • Batch files: one or more IDL statements or commands. $ @batchfile: interpreted one by one, exactly as if it was from the keyboard. • Main-level Programs: a series of program statements that are compiled and executed once an END statement is encountered. • Programs: pro name, param1, param2, ... paramx ended with END • Functions: function name, param1, param2, ... Paramx ended with END • Variable Access: Batch and Main (globe), Program & function (local)

  15. Parameters Passing • Actual (caller) and Formal (called) Parameters Correspondence by position or keyword Keyword Inheritance (_Extra) • Passing Mechanism • Expressions, constants, system variables, and subscripted variable references are passed by value. • Variables are passed by reference. • Parameters and Keywords Checking • n_params()：number of parameters in calling an procedure/ function • n_elements()：returns zero for undefined variable. • keyword_set()：check a Boolean keyword parameter. • arg_present()：defined and reference passing?

  16. Working with Arrays! (Essential for IDL) A(n,m): Array with n columns and m rows A[n,m]: Element A[i1, j:*] pos = n*i+j j = pos MOD n i = (pos - j)/n

  17. Array-Oriented Operation • Arrays work the same as Scalars. e.g. 2*A, A + B, A-B, A/B, SQURT(A), …… Try to avoid use of loops (slow!) • Array Creation: xxxArr, xIndGen, Replicate, … • Array Manipulation: [], ARRAY_INDICES, CONGRID, HISTOGRAM, INVERT, MAX, MIN, MEDIAN, N_ELEMENTS, REBIN, REFORM, REVERSE, ROT, ROTATE, SHIFT, SIZE, SORT, TOTAL, MEAN, TRANSPOSE, UNIQ, WHERE, etc.

  18. WHERE much faster than IF Set all values between 5 to 8 equal = 15. Loop Way: For j=0,2 Do Begin For k=0,3 Do Begin IF (array(j,k) GE 5) AND (array(j,k) LE 8) THEN array(j,k) = 15 EndFor EndFor IDL Way: index = Where((array GE 5) AND (array LE 8), count) IF count GT 0 THEN array[index] = 15 http://www.dfanning.com/powerpoint/index.html

  19. Where(): where are they? file= FILEPATH('galaxy.dat', subdir = ['examples', 'data']) imagesize = [256, 256] image = READ_BINARY(file, data_dims =imagesize) DEVICE, decomposed = 0 & LOADCT, 4 WINDOW, 0, xsize = imagesize[0], ys= 2 * imagesize[1] indices = Where((image GE 200) AND $ (image LE 230), count) IF count GT 0 THEN BEGIN result = Array_Indices(image, indices) col = Reform(result [0,*]) row = Reform(result [1,*]) TV, image, 0 & TV, image, 1 PlotS, col, row, /Device, $ Color=FSC_Color(‘white'), psym=1 ENDIF

  20. Mask using Where() file = FILEPATH('worldelv.dat', $ subdir = ['examples', 'data']) file = FILEPATH('worldtmp.png', $ subdir = ['examples', 'demo', 'demodata']) TV, elvImage, 0 & TV, tmpImage, 1 ocean = WHERE(elvImage LT 125) image = tmpImage image[ocean] = elvImage[ocean] TV, image, 2 land = WHERE(elvImage GE 125) image = tmpImage image[land] = elvImage[land] TV, image, 3 Temperature Distribution in Land and Ocean

  21. Index Manipulation Simultaneous look at two images. img1 = BytScl(Loaddata(4), Top=99) img2 = BytScl(Loaddata(5), Top=99)+100B Window, XSize=256, YSize=256*3 LoadCT, 13, NColors=100 & TV, img1, 2 LoadCT, 3, NColors=100, Bottom=100 & TV, img2, 1 LoadCT, 13, NColors=100 LoadCT, 3, NColors=100, Bottom=100 index = $ Where((Indgen(256L*256L) MOD 2) EQ 0) img1[index] = img2[index] TV, img1, 0 http://www.dfanning.com/powerpoint/index.html

  22. Examples: the IDL Way e.g., Upper Triangular Matrix n = 5 i = REBIN(LINDGEN(n), n, n) j = REBIN(TRANSPOSE(LINDGEN(n)), n, n) Print, i Print, j 0 1 2 3 4 0 0 0 0 0 0 1 2 3 4 1 1 1 1 1 0 1 2 3 4 2 2 2 2 2 0 1 2 3 4 3 3 3 3 3 0 1 2 3 4 4 4 4 4 4 mask = (i GE j) Print, mask 1 1 1 1 1 0 1 1 1 1 0 0 1 1 1 0 0 0 1 1 0 0 0 0 1 http://www.dfanning.com/idl_way/

  23. Direct Graphics • Device oriented; display in a specific device • Set_plot: X/WIN/MAC, PS, PRINTER, METAFILE, Z, CGM, PCL, NULL • Device: retain, decomposed, set_character_size, pseudo_color, index_color, true_color • Window Coordinates: Data, Device, Nomal

  24. PS files Creation curName=!D.name ; ‘X’ or ‘Win’ set_plot,'ps' device, file=‘test.ps‘, posi=[1,1,9.5,9]/10, bits_per_pixel=8 ;device, file=‘test.eps‘, posi=[1,1,9.5,9]/10 ,/encapsulated ;device, file=‘test.ps‘, posi=[1,1,9.5,9]/10,/color plot, loaddata(1) device,/close & set_plot,curName PS Layout configurations: http://www.dfanning/documents/programs.html http://cow.physics.wisc.edu/~craigm/idl/printing.html For X/MS windows devices, write_jpeg, ‘test.jpg', tvrd() write_jpeg, ‘test.jpg', tvrd(true=1),true=1

  25. Working with Colors • (R,G,B) triple:any Color decomposed to Red, Green, and Blue components, each with value 0~255 • Indexed Color Model and RGB Color Model • Indexed Color Model (also 8 bit color): Color Lookup Table  Color Palette 2^8 = 256 colors Device, decomposed = 0 • RGB Color Model (also 24 bit color) : specify color values explicitly, using an RGB triple 2^8 * 2^8 * 2^8 = 16777216 colors Device, decomposed = 1 Indexed? Dynamic HELP, /DEVICE LOADCT, TVLCT, XLOADCT, XPALETTE About COLORBAR, FSC_COLOR, etc. http://www.dfanning/documents/programs.html

  26. Comparison: 24-Bit and 8-Bit Color http://www.dfanning.com/powerpoint/index.html

  27. Color Decomposition Data=loaddata(1) Plot, data, Color=255, posi=[0,0,1,1] Help, /Device Device, Decomposed=1 ; ON Plot, data, Color=11829830L,$ Background='00ff00'xL, posi=[0,0,1,1] TVLCT, 70, 130, 180, 240 Device, Decomposed=0; OFF Plot, data, Color=240,$ Background=255,posi=[0,0,1,1]

  28. RGB Image: Ant Nebula device,decom=1 window, 0, xsize=2*info[2], ysize=2*info[3] tv, ant, true=1,0 tvscl,r,channel=1 & tvscl,g,channel=2& tvscl,b,channel=3 tvscl,r,channel=1,1 &tvscl,g,channel=2,2& tvscl,b,channel=3,3 read_jpeg,'Ant.jpg',ant r=reform(ant[0,*,*]) g=reform(ant[1,*,*]) b=reform(ant[2,*,*]) Help, ant Google Images: Ant Nebula

  29. Image Types • Byte: 0~255; otherwise, BYTSCL() • Binary Image: values only with 0 or 1 • Gray-scale Image: B-W LUT • Indexed Images: 0~255, LUT • RGB Images: (R,G,B), each 0~255 COLOR_QUAN() : RGB Indexed Image FILEPATH, QUERY_IMAGE, READ_BINARY, READ_IMAGE READ_JPEG/WRITE_JPEG, WRITE_IMAGE,……

  30. III. 1D and 2D Data Display • plot/oplot, plots, axis, xyouts, ploterr/oploterr/errplot/, ploterror(!!!), vel/velovect/plot_field/flow3 • BTLSCL, REBIN, REFORM • tv/tvscl/bytscl, imdisp(!!!), plot_image, tvimage, contour, image_cont, surface, surf_shade, show3, median, smooth/convol, reberts/sobel, defroi, profiles, etc. • Histogram, plot_hist • box_cursor, plot_box, rdpix, curval • live_tools(live_plot,...), itools(iplot,icontour,...) • iTools (iplot, ……) • !p.multi, position, xrange, xstyle, psym, ……

  31. User PSYM device, decomposed=0 !p.multi=[0,1,2] data=loaddata(1) Plot, data, PSym=-2 ;filled circle phi = Findgen(32) * (!PI * 2 / 32.) phi = [ phi, phi[0] ] UserSym, Cos(phi), Sin(phi), /Fill Plot, data, /NoData ;tvlct, 178,34,34,10 & OPlot, data, Color=10 OPlot, data, Color=FSC_Color('firebrick') Oplot, data, color=FSC_Color('forest green'), $ PSym=8, Symsize=1.5

  32. Plotting Error Bars xtime = indgen(101) & data = loaddata(1) xerr = randomN(seed,101)*2 yerr = randomN(seed,101)*4 device,decomposed=0 & loadct, 0 tvlct, r, g, b, /get & tvlct, 255-[[r],[g],[b]] & tvlct, 0,255,0,100 ploterror, xtime, data, xerr, yerr, psym = -2, xstyle=1, $ xtitle='!7b(!6cm!U-2!N !6s!D-1!N)', ytitle='!6H!7a' oplot,xtime,data,color=100,thick=2 ERRPLOT, X, Y-yerr1, Y+yerr2

  33. TV and TVSCL file = FILEPATH('hurric.dat', subdir = ['examples', 'data']) hurric = READ_BINARY(file, DATA_DIMS = [440, 340])

  34. Pseudo-Color Images in PS file = FILEPATH('worldelv.dat', subdir = ['examples', 'data']) image = READ_BINARY(file, DATA_DIMS = [360, 360]) curName=!d.name & set_plot,'PS' device,/color,file='elev.eps',bits_per_pixel=8,/encap loadct,13 & imdisp, image ; plot_image,image device,/close & set_plot,curName IMDISP better than TV/TVSCL Also try PLOT_IMAGE

  35. True-Color Images in PS file=filepath('rose.jpg',subdir='examples/data') read_jpeg,file,rose & tvlct,r,g,b,/get curName=!d.name & set_plot,'PS‘ device,/color,file='rose.eps',bits_per=8,/encap loadct,0 & imdisp,rose device,/close & set_plot,curName tvlct,r,g,b Color table always active in 24 bits mode, ps device Google Images: RosetteNebula

  36. VI. FITS I/O in IDL • FITS (Flexible Image Transport System) is a standardized data format which is widely used in astronomy. • Briefly, a FITS file consists of a sequence of one or more Header and Data Units (HDUs). A header is composed of ASCII card images that in IDL is usually read into a string array variable. The header describes the content of the associated data unit, which might be a spectrum (IDL vector), an image (IDL array), or tabular data in ASCII or binary format (often read as an IDL structure). Image and vector data can be present in any HDU, but tabular data cannot appear in the first HDU. The HDUs following the first (or primary) HDU are also known as extensions, and thus a FITS file containing tabular data must contain at least one extension. Four Classes of Procedures: MRDFITS()/MWRFITS: READFITS()/WRITEFITS FX* Procedures FITS_* and FTAB_* Procedures FITS I/O in IDLAstro http://idlastro.gsfc.nasa.gov/fitsio.html The FITS Support Office/NASAhttp://fits.gsfc.nasa.gov/

  37. FITS I/O: File Information http://fits.gsfc.nasa.gov/fits_samples.html IDL>images=mrdfits(file,0,head0) MRDFITS: Image array (200,200,4) Type=Real*4 IDL> help,images,head0 IMAGES FLOAT = Array[200, 200, 4] HEAD0 STRING = Array[263] IDL>table=mrdfits(file,0,head1) MRDFITS: Image array (200,200,4) Type=Real*4 IDL> help,table,head1 TABLE FLOAT = Array[200, 200, 4] HEAD1 STRING = Array[263] IDL> file =‘WFPC2u5780205r_c0fx.fits’ IDL>Fits_Info, ‘WFPC2u5780205r_c0fx.fits’ WFPC2u5780205r_c0fx.fits has 1 extensions Primary header: 263 records Image -- Real*4 array ( 200 200 4 ) Extension 1 -- u5780205r_cvt.c0h.tab Header : 354 records ASCII Table ( 796 4 ) IDL>file=‘EUVEngc4151imgx.fits’ IDL>fits_help, file XTENSION EXTNAME EXTVER EXTLEVEL BITPIX GCOUNT PCOUNT NAXIS NAXIS* 0 8 0 0 0 1 IMAGE ds 16 1 0 2 512 x 512 2 IMAGE sw_night 16 1 0 2 2048 x 300 3 IMAGE mw 16 1 0 2 2048 x 300 4 IMAGE lw 16 1 0 2 2048 x 300 5 BINTABLE ds_limits 8 1 0 2 16 x 3 6 BINTABLE sw_night_limits 8 1 0 2 20 x 2 7 BINTABLE mw_limits 8 1 0 2 20 x 2 8 BINTABLE lw_limits 8 1 0 2 20 x 2

  38. Fits Example: Orion Nebula http://casa.colorado.edu/~bally/ HST/HST/master/ data=readfits('masterf673.fits.gz',head) loadct,10 & tvlct,r,g,b,/get r[0]=34 & g[0]=139 & b[0]=34 & tvlct,r,g,b plot_image,alog10(data>5e-1) loadct,10 plot_image,alog10(data[1000:2000,1200:2200]>0.5)

  39. Solar Map Software An IDL map is a structure that contains two-dimensional (2-d) image data with accompanying pixel coordinate and spatial scale information. The latter parameters are defined as properties of the map and are unique for each image source. Defined in this manner, an arbitrary image can be manipulated or transformed in a manner that is independent of the image source. http://orpheus.nascom.nasa.gov/~zarro/idl/maps.html http://hesperia.gsfc.nasa.gov/~ptg/trace-align/

  40. V. Image Processing A Short Introduction to Digital Image Processing http://web.uct.ac.za/depts/physics/laser/hanbury/intro_ip.html

  41. BYTSCL (Byte Scale) mr_brain BYTSCL READ_DICOM(FILEPATH('mr_brain.dcm', subdir = ['examples', 'data']))

  42. 01/23/2001 19:19:43 Dragon in EIT 304 http://umbra.nascom.nasa.gov/eit/eit-catalog.html

  43. Color Table Highlights, Contrast file = FILEPATH('mineral.png', subdir = ['examples', 'data']) image = READ_PNG(file, r, g, b) colorLevel = [[0, 0, 0], [255, 0, 0], [255, 255, 0], [0, 255, 0], $ [0, 255, 255], [0, 0, 255], [255, 0, 255], [255, 255, 255]] numberLevel = CEIL(!D.TABLE_SIZE/8.) level = INDGEN(!D.TABLE_SIZE)/numberLevel newR = colorLevel(0, level) & newR[!D.TABLE_SIZE - 1] = 255 newG = colorLevel(1, level) & newG[!D.TABLE_SIZE - 1] = 255 newB = colorLevel(2, level) & newB[!D.TABLE_SIZE - 1] = 255 TVLCT, r,g,b TVLCT, 13 TVLCT, newR, newG, newB

  44. Histogram Equalization Left:BYTSCL(HISTOGRAM(image)) Right: Mineral Left: BYTSCL(HISTOGRAM(hq_image)) Right: hq_image = HIST_EQUAL(image) also, H_EQ_CT, H_EQ_INT ADAPT_HIST_EQUAL Left: HISTOGRAM(equalizedImage) Right: ADAPT_HIST_EQUAL(image)

  45. Adjust Histogram http://web.uct.ac.za/depts/physics/laser/hanbury/intro_ip.html

  46. Fast Fourier Transform de-noise file= FILEPATH('abnorm.dat', subdir = ['examples', 'data']) powerSpec = ALOG10(SHIFT(FFT(image), zz[0]/2, zz[0]/2)) CONGRID(powerSpec, zz2[0], zz2[1]) mask = FLOAT(scaledSpec) GT 2.6 maskedSpec = scaledSpec * mask + MIN(powerSpec) inverseTran = ABS(FFT(SHIFT(10.^(maskedSpec), $ zz[0]/2,zz[1]/2), /inverse)) scaledSpec = powerSpec - Min(powerSpec)

  47. FFT: removing corrugated effect http://web.uct.ac.za/depts/physics/laser/hanbury/intro_ip.html

  48. Inverse Laplace Trans inverseTrans maskedSpec mask = HANNING(imagesize[0], imagesize[1]) maskedSpec = (scaledSpec * mask) + MIN(powerSpec) inverseTrans = ABS(FFT(SHIFT(10.^(maskedSpec), $ imagesize[0]/2, imagesize[1]/2), /inverse))

  49. Median Smoothing MEDIAN(rbcells, 5) rbcells file= FILEPATH('rbcells.jpg', subdir = ['examples', 'data'])

  50. Rim Enhancing croppedImage croppedSize = [96, 96] file= FILEPATH('nyny.dat', subdir = ['examples', 'data']) croppedImage = image[200 : croppedSize[0] - 1 + 200, 180 : croppedSize[0] - 1 + 180] ROBERTS SOBEL

More Related