PRO elplot_setsym, index, dots, fill
  CASE index OF
     -1: BEGIN
        psympoly,SIZE=dots,FILL=fill,N=3,ALPHA=0.
     END
     -2: BEGIN
        psympoly,SIZE=dots,FILL=fill,N=3,ALPHA=60.
     END
     -3: BEGIN
        psympoly,SIZE=dots*0.8,FILL=fill,N=30,ALPHA=0.
     END
     -4: BEGIN
        psympoly,SIZE=dots,FILL=fill,N=4,ALPHA=0.
     END
     -5: BEGIN
        psympoly,SIZE=dots,FILL=fill,N=4,ALPHA=45.
     END
     -6: BEGIN
        psymdot,SIZE=dots,FILL=fill
     END
  ENDCASE
END

;-----------------------------------------------------------------------

PRO elplot, $
   yield0, elz0, $
   COLOR=color, $
   XRANGE=xrange, $
   STYLE=style, $
   LEGEND=legend, $
   YRANGE=yrange, $
   OBSDATA=obsdata0, $
   NORMALIZE=normalize, $
   LEG_CHARSIZE=leg_charsize, $
   LEG_LENGTH=leg_length, $
   LEG_XPOS=leg_xpos, $
   LEG_YPOS=leg_ypos, $
   XMARGIN=xmargin, $
   YMARGIN=ymargin, $
   YTITLE=ytitle, $
   XTITLE=xtitle, $
   COMBINE=combine, $
   EPSILON=epsilon, $
   NEP=nep, $
   FITIONS=fitions, $
   WEIGHTS=weights, $
   LIMITDOTS=limitdots, $
   DOTSIZE=dotsize0, $
   LABELSIZE=labelsize, $
   LABELTHICK=labelthick, $
   HOLLOW=hollow, $
   LINESTYLE=linestyle, $
   LINECOLOR=linecolor, $
   LINETHICK=linethick, $
   YMADD=ymarginadd, $
   LEG_STYLE=legstyle, $ ; new
   LEG_COLUMNS=leg_columns ; new
  

@compile_opt

; make a copy as to not overwrite original data
yield=yield0

elements_up=1
@elements
el=elements[elz0]

;combine=[{z:[6,7,8]}]
;combine=[{z:[6,7]}]
;combine=[{el:['C','N','O'],name:'C+N+O'}]

IF N_ELEMENTS(obsdata0) NE 0 THEN BEGIN
;  make local copy
    obsdata=obsdata0
END

IF N_ELEMENTS(epsilon) EQ 0 THEN epsilon=0
ctype='pf'
IF epsilon EQ 1 THEN BEGIN
    ctype='Y'
    IF N_ELEMENTS(obsdata) NE 0 THEN BEGIN
        obsdata=pf2eps(obsdata,LOG=1)
    ENDIF
ENDIF

IF N_ELEMENTS(obsdata) NE 0 THEN BEGIN
    eldatacombine, $
      DATA=obsdata, $
      CDATA=cdata, $
      COMBINE=combine, $
      CTYPE=ctype, $
      LOG=1
    obsdata=cdata
ENDIF

x=SIZE(yield,/DIMENSIONS)
IF N_ELEMENTS(x) EQ 1 THEN nyields=1 ELSE nyields=x[1]

data=REPLICATE({z:0L,sigma:0.0D0,abu:0.0D0,el:''},N_ELEMENTS(elz0),nyields)
FOR i=0,nyields-1 DO BEGIN
    data[*,i].z=elz0[*]
    data[*,i].el=el[*]
    data[*,i].abu=yield[*,i]
    
    IF epsilon EQ 1 THEN BEGIN
        data[*,i]=pf2eps(data[*,i],LOG=0)
    ENDIF
    eldatacombine, $
      DATA=data[*,i], $
      CDATA=cdata, $
      COMBINE=combine, $
      CTYPE=ctype, $
      LOG=0
    ncdata=N_ELEMENTS(cdata)
    yield[0:ncdata-1,i]=cdata.abu
ENDFOR
el=cdata.el
elz=cdata.z
yield=yield[0:ncdata-1,0:nyields-1]

obsdataplot=0
IF N_ELEMENTS(obsdata) NE 0 THEN obsdataplot=1

IF N_ELEMENTS(color) EQ 0 THEN color=1
cf=!P.COLOR
IF color EQ 1 THEN BEGIN
   IF !D.N_COLORS LE 256 THEN BEGIN
      cl=1
      cd=2
      ca=3
      TVLCT,255*[0,1,0],255*[0,0,0.5],255*[1,0,0],1
   ENDIF ELSE BEGIN
      cl=RGB(0.,0.,1.)
      cd=RGB(1.,0.,0.)
      ca=RGB(0.,0.,0.)
   ENDELSE
ENDIF ELSE BEGIN
   cl=!P.COLOR
   ca=!P.COLOR
   cd=!P.COLOR
ENDELSE
IF N_ELEMENTS(linestyle) EQ 0 THEN BEGIN
   l=indgen(nyields)
ENDIF ELSE BEGIN
   IF N_ELEMENTS(linestyle) LT nyields THEN BEGIN
      nl=N_ELEMENTS(linestyle)
      l=[linestyle,REPLICATE(linestyle[nl-1],nyields-nl)]
   ENDIF ELSE BEGIN
      l=linestyle
   ENDELSE
ENDELSE
IF N_ELEMENTS(linethick) EQ 0 THEN BEGIN
   tl=REPLICATE(!P.THICK,nyields)
ENDIF ELSE BEGIN
   IF N_ELEMENTS(linethick) LT nyields THEN BEGIN
      nl=N_ELEMENTS(linethick)
      tl=[linethick,REPLICATE(linethick[nl-1],nyields-nl)]
   ENDIF ELSE BEGIN
      tl=linethick
   ENDELSE
ENDELSE
IF N_ELEMENTS(linecolor) EQ 0 THEN BEGIN
   cl=REPLICATE(cl,nyields)
ENDIF ELSE BEGIN
   IF N_ELEMENTS(linecolor) LT nyields THEN BEGIN
      nl=N_ELEMENTS(linecolor)
      cl=[linecolor,REPLICATE(linecolor[nl-1],nyields-nl)]
   ENDIF ELSE BEGIN
      cl=linecolor
   ENDELSE
ENDELSE

gg=0.7
IF !D.N_COLORS LE 256 THEN BEGIN
    cg=4
    TVLCT,255*gg,255*gg,255*gg,4
ENDIF ELSE BEGIN
    cg=RGB(gg)
ENDELSE
t=!P.THICK
IF t EQ 0 THEN t=1

nmax=MAX(elz,MIN=nmin)
IF N_ELEMENTS(xrange) NE 2 THEN xrange=[nmin-1,nmax+1]

nel=N_ELEMENTS(elz)

; normalize data if requested
IF N_ELEMENTS(normalize) EQ 0 THEN normalize=0
type=SIZE(normalize,/TNAME)
IF type EQ 'STRING' THEN BEGIN
    normalize=el2z(normalize)
ENDIF
; ???
; IF (obsdataplot NE 0) AND nyields GT 1 THEN normalize=12

ytitle_norm=''
IF normalize NE 0 THEN BEGIN
    img=WHERE(elz EQ normalize)
    img=img[0]
    FOR iy=0,nyields-1 DO BEGIN
        yield[*,iy]/=yield[img,iy]
    ENDFOR
    img=WHERE(obsdata.z EQ normalize)
    img=img[0]
    obsdata.abu-=obsdata[img].abu
    ytitle_norm=elements[normalize]+' '
ENDIF

; normalize plot
IF N_ELEMENTS(nep) EQ 0 THEN nep=0
nep_type=SIZE(nep,/TNAME)
IF nep_type EQ 'STRING' THEN BEGIN
    nep=el2z(nep)
ENDIF

IF nep GT 0 THEN BEGIN
    inep=WHERE(obsdata.z EQ nep)
    inep=inep[0]
    IF inep NE -1 THEN BEGIN
        yoffset=obsdata[inep].abu
        yield/=10.D0^yoffset
        obsdata.abu-=yoffset
        ytitle_norm='/ '+elements[nep]+' '
    ENDIF ELSE BEGIN
        PRINT,' [ELPOT] ERROR: NEP = '+STRTRIM(STRING(nep),2)+' NOT FOUND.'
        nep=0
    ENDELSE
ENDIF

IF N_ELEMENTS(xmargin) EQ 0 THEN xmargin=[7,0.1]
IF N_ELEMENTS(ymargin) EQ 0 THEN ymargin=[3,0.25]
IF N_ELEMENTS(ymarginadd) EQ 0 THEN ymarginadd=[0,0]
ymargin+=ymarginadd
IF N_ELEMENTS(title) EQ 0 THEN BEGIN
    title=''
ENDIF ELSE BEGIN
    ymargin[1]+=1.5
ENDELSE

IF N_ELEMENTS(yrange) EQ 0 THEN yrange=[-2,2]
IF yrange[0] GT yrange[1] THEN yrange=yrange[[1,0]]
IF epsilon EQ 0 THEN BEGIN
    IF nep GT 0 THEN BEGIN
        ytitle0='[ X '+ytitle_norm+']'
    ENDIF ELSE BEGIN
        ytitle0='log( element number fraction relative to solar '+ytitle_norm+')'
    ENDELSE
ENDIF ELSE BEGIN
    IF ytitle_norm NE '' THEN ytitle_norm='relative to '+ytitle_norm+' ' 
    ytitle0='log( epsilon '+ytitle_norm+')'
ENDELSE
IF N_ELEMENTS(ytitle) EQ 0 THEN ytitle=ytitle0
IF N_ELEMENTS(xtitle) EQ 0 THEN xtitle='atomic number'

plot,[0,1],/nodata, $
  yr=yrange,/ys, $
  xr=xrange,/xs, $
  xtitle=xtitle, $
  ytitle=ytitle, $
  TITLE=title, $
  COLOR=cf, $
  XMARGIN=xmargin, $
  YMARGIN=ymargin

IF N_ELEMENTS(dotsize0) NE 0 THEN dotsize=dotsize0
IF N_ELEMENTS(dotsize) EQ 0 THEN dotsize=1.2
IF !D.NAME EQ 'X' THEN dotsize=2
dots=dotsize

y=ALOG10(yield)
FOR iy=0,nyields-1 DO BEGIN
    oplot,elz,y[*,iy],COLOR=cl[iy],LINE=l[iy],THICK=tl[iy]
ENDFOR

IF N_ELEMENTS(style) EQ 0 THEN style=0
CASE style OF
    0:BEGIN
        psym=[-1,-2,-3,-4,-5,REPLICATE(0,99)]
    END
    1:BEGIN
        psym=[-6,REPLICATE(0,99)]
    END
    2:BEGIN
        psym=REPLICATE(0,99)
    END
END

IF N_ELEMENTS(leg_charsize) EQ 0 THEN leg_charsize=1.D0
IF N_ELEMENTS(leg_length) EQ 0 THEN leg_length=-1
IF leg_length EQ -1 AND N_ELEMENTS(legend) NE 0 THEN leg_length=MAX(STRLEN(LEGEND))

IF N_ELEMENTS(leg_xpos) EQ 0 THEN leg_xpos=1
IF N_ELEMENTS(leg_ypos) EQ 0 THEN leg_ypos=1

LEGENDE,leg_length,/RIGHT,/TOP,CHARSIZE=LEGEND_CHARSIZE, $
  XPOS=leg_xpos, YPOS=leg_ypos

IF N_ELEMENTS(hollow) EQ 0 THEN hollow=''

FOR iy=0,nyields-1 DO BEGIN
   IF dots NE 0 THEN BEGIN      
      FOR iel=0,nel-1 DO BEGIN
         fill=1
         IF (WHERE(hollow EQ el[iel]))[0] NE -1 THEN BEGIN
            fill=0
         ENDIF
         IF psym[iy] LT 0 THEN BEGIN
            elplot_setsym,psym[iy],dots,fill
            plots,elz[iel],y[iel,iy],psym=8,COLOR=cd,NOCLIP=0
         ENDIF ELSE if psym[iy] GT 0 THEN BEGIN
            plots,elz[iel],y[iel,iy],psym=psym[iy],COLOR=cd,NOCLIP=0
         ENDIF
      ENDFOR
   ENDIF

; ???   
   IF N_ELEMENTS(legend) EQ nyields THEN BEGIN
      IF dots NE 0 THEN nonewline=1 ELSE nonewline=0
      legende_line,l[iy],legend[iy],COLOR=cl,NONEWLINE=nonewline
      IF dots NE 0 THEN BEGIN
         psymx=psym[iy]         
         IF psym[iy] LT 0 THEN BEGIN
            psymx=8
            fill=1
            elplot_setsym,psym[iy],dots,fill
         ENDIF
         legende_symbol,psymx,'',COLOR=cd
      ENDIF
   ENDIF
ENDFOR

; annotation setup
IF N_ELEMENTS(labelsize) EQ 0 THEN labelsize=1
IF !P.charsize NE 0 THEN BEGIN
    charsize=!P.charsize*labelsize
ENDIF ELSE BEGIN
    charsize=labelsize
ENDELSE
cx=!D.X_CH_SIZE*charsize
cy=!D.Y_CH_SIZE*charsize
nc0=CONVERT_COORD([0,0],/DEVICE,/TO_DATA)
nc1=CONVERT_COORD([cx,cy],/DEVICE,/TO_DATA)
nc=nc1-nc0
flayout=[0.2,-0.9]*nc[1]
yframe=yrange+flayout+0.2*[+1,-1]*nc[0]

fmin=y[*,0]
fmax=fmin
FOR i=0,nel-1 DO BEGIN
    fmin[i]=MIN(y[i,*])
    fmax[i]=MAX(y[i,*])
ENDFOR

datadot=dots*0.3*[1,-1]*nc[1]
fmax+=datadot[0]
fmin+=datadot[1]

fmin0=fmin
fmax0=fmax

fborder=0.1*gettextwidth('x',/data)
FOR i=1,nel-2 DO BEGIN
   ffrac0=((0.5D0*gettextwidth(el[i],/data)+fborder) $
          /[elz[i]-elz[i-1],elz[i+1]-elz[i]]) $
         <1>0
   ffrac1=(1.D0-ffrac0)<1>0
   fpos=ffrac1*fmin0[i]+ffrac0*fmin0[[-1,+1]+i]
   fmax[i]>=MAX(fpos)-datadot[0]
   fmin[i]<=MIN(fpos)-datadot[1]
ENDFOR

IF obsdataplot NE 0 THEN BEGIN
; this scale needs to be the same as in abudataplot.pro
    scale=CONVERT_COORD([0,1],[0,1],/DEVICE,/TO_DATA)
    dx=ABS(scale[0,0]-scale[0,1])
    dy=ABS(scale[1,0]-scale[1,1])
    dec=0.25D0
    scale=dec/dx*dy
    alen=7.5*scale
    fh=(flayout[0]-flayout[1]+0.2*nc[1])
    ch=CONVERT_COORD([0,!D.X_CH_SIZE],[0,!D.X_CH_SIZE],/DEVICE,/TO_DATA)
    ddot=ch[*,1]-ch[*,0]
    ndata=N_ELEMENTS(obsdata)

    IF N_ELEMENTS(fitions) EQ 0 THEN fitions=z2el(obsdata.z)

    filledions=fitions

    IF N_ELEMENTS(weights) EQ 0 THEN weights=1
    IF (N_ELEMENTS(weights) EQ 1) AND (ndata GT 1) THEN BEGIN
        weights=REPLICATE(weights,N_ELEMENTS(fitions))
    ENDIF

    FOR i=0,N_ELEMENTS(fitions)-1 DO BEGIN
        ion=fitions[i]
        filledions[i]=STRJOIN(STRSPLIT(ion,'+',/EXTRACT))
    ENDFOR

    filled=INTARR(ndata)
    obsweight=REPLICATE(1.D0,ndata)
    FOR i=0,ndata-1 DO BEGIN
        ii=(WHERE(obsdata[i].el EQ filledions))[0]
        IF ii NE -1 THEN BEGIN
            filled[i]=1
            obsweight[i]*=weights[ii]
        ENDIF
    ENDFOR
    dotsize=SQRT(obsweight)

    FOR i=0,ndata-1 DO BEGIN
        ii=WHERE(obsdata[i].z EQ elz)
        ii=ii[0]
        IF obsdata[i].sigma GT 0 THEN BEGIN
            lo=obsdata[i].abu-obsdata[i].sigma
            hi=obsdata[i].abu+obsdata[i].sigma
        ENDIF ELSE BEGIN
            IF obsdata[i].sigma GT -1.D99 THEN BEGIN 
                veclen=(-obsdata[i].sigma);<alen?
            ENDIF ELSE BEGIN
                veclen=alen
            ENDELSE
            lo=obsdata[i].abu-veclen
            hi=obsdata[i].abu+0.4*ddot[1]*dotsize[i]
        ENDELSE
        IF (fmin[ii] GT lo) AND (fmin[ii] LT hi+fh) THEN fmin[ii]<=lo 
        IF (fmax[ii] LT hi) AND (fmax[ii] GT lo-fh) THEN fmax[ii]>=hi 
        IF (fmin[ii] LT yframe[0]) THEN fmax[ii]>=hi
        IF (fmax[ii] GT yframe[1]) THEN fmin[ii]<=lo
    ENDFOR
ENDIF

yav=0.5D0*(fmin0+fmax0)

; show off-scale elements - maybe not needed
; fmax>=yframe[0]
; fmin<=yframe[1]

FOR i=0,nel-1 DO BEGIN
    pos=-1
    IF (i EQ 0) THEN BEGIN
        IF yav[i+1] GT yav[i] THEN pos=1 ELSE pos=0
    ENDIF ELSE IF (i LT nel-1) THEN BEGIN
        IF (yav[i+1]*(elz[i]-elz[i-1])+yav[i-1]*(elz[i+1]-elz[i])) GT (elz[i+1]-elz[i-1])*yav[i] THEN pos=1 ELSE pos=0
    ENDIF ELSE BEGIN
        IF yav[i-1] GT yav[i] THEN pos=1 ELSE pos=0
    ENDELSE


; IF i EQ 6 then xxx

    yyx=flayout+[fmax[i],fmin[i]]

    IF pos EQ 1 AND yyx[1] LT yframe[0] THEN pos=0
    IF pos EQ 0 AND yyx[0] GT yframe[1] THEN pos=1

    XYouts,elz[i],yyx[pos],el[i],alignment=0.5,COLOR=ca,NOCLIP=0, $
      CHARSIZE=labelsize, $
      CHARTHICK=labelthick

ENDFOR

IF obsdataplot NE 0 THEN BEGIN 
    abudataplot, $
      obsdata, $
      FILLED=filled, $
      DOTSIZE=dotsize, $
      LIMITDOTS=limitdots
ENDIF

END