require'plot' NB. Plot polarization over planet's projection on sky NB. Usage: (pensize;tilt mag;vect) pol_pict (long lat S phase) NB. vect= scale of line segments showing polarization NB. "tilt mag" > 1 amplifies deviation from 0 deg (vertical) polarization NB. S is the array of Stokes parameters from "PLANET.ijs" pol_pict=: dyad define 'll pp S phase'=. y 'I p a'=. FracPol S 'f avp avang'=. ,FracPol +/"1 +/"1 S lc=. ll- phase zz=. (lc<360)*/(sin DTR*pp) NB. projected latitudes xx=. -(sin DTR*lc)*/(cos DTR*pp) NB. projected longitudes a=. DTR*a+(a>90)*_180 NB. polarization angle (rad) 'sz tilt scale'=. x a=. tilt*a NB. magnify angle dx=. p* sin a dz=. p* cos a options=. 'title Polarization from Illuminated Scattering Atmosphere' options=. options,';axes 0 0;grids 0 0;aspect 0.5' size=. 'pensize ',":sz options=. options,';type point;color blue;',size options=. options,';frame 1;xrange _1 1;yrange 0 1' ggg=. ';xcaption Phase= ',(":phase) ggg=. ggg,', Pol. tilt scale= ',(":tilt),', Ave Pol= ',":avp options=. options,ggg dels=. 0.005*scale* _8+ i. 17 xvec=. , (((#dels),$xx)$,xx)+ dels*/dx zvec=. , (((#dels),$zz)$,zz)+ dels*/dz pd 'reset' pd options pd xvec;zvec roundup=. (>. 10*avp)%10 rule=. 0.005*scale*roundup*(_8+16*int01 10) xspot=: 11#0.93 zspot=: 0.96+rule kname=: (":100*roundup),'% polarization:' pd xspot;zspot pd 'textc 850 790 ',kname pd 'show' pd 'pdf' NB. write plot to .pdf file )