pro read_vdf_2
;------------------------------------------------------------------------------
;
;   programme to read the Helios ion/electron velocity distribution
;   functions, vdfp, vdfe and the corresponding velocity components.
;   program for the 3-D proton (2-D electron) velocity distribution function
;   and determination of the velocity ranges and maximum vdf value
;
;-----------------------------------------------------------------------------


;------------------- variables and arrays ------------------------------------

comment=strarr(1)



vdfp=dblarr(12,12,32)      ; protons/alpha-particles  3-D cubes
vxp =dblarr(12,12,32)
vyp =  dblarr(12,12,32)
vzp = fltarr(12,12,32)

vbp_ll =dblarr(12,12,32)
vbp_perp =dblarr(12,12,32)

countp=intarr(12,12,32)

vdfe=dblarr(8,16)           ; electrons  2-D squares
vxe =dblarr(8,16)
vye =dblarr(8,16)
counte=intarr(8,16)

print,'*********'
print,123
print, 1

end







l = 12*12*32         ;= 4608  maximum proton arry size
print, 'lp=', l
l = 8*16             ;= 128   maximum electron array size
print, 'le=', l
ndimp = 4608
ndime = 128


end

zp =  fltarr(4608)    ; 1-D array protons
wxp = fltarr(4608)
wyp = fltarr(4608)
wzp = fltarr(4608)

wp_ll=fltarr(4608)
wp_perp=fltarr(4608)

ze =  fltarr(128)     ; 1-D array electrons
wxe = fltarr(128)
wye = fltarr(128)
print,2
;---------------------- definitions of variables ------------------------------
;
;       n_e                  = number of electron data
;       n_p                  = number of proton data
;	iearcarrot           = Carrington rotation number from Earth
;       ihelcarrot           = Carrington rotation number from spacecraft
;       i1a_3_on             = instrument i1a on (=0) or i3 on (=1)
;       index_shift
;       minus                = indicates a HDM file which contained bad
;                              data (frames), but could be handled as NDM file
;	mode                 = HDM or NMD
;       iperihelion_shift    = on or of (0 or 1)
;       idummy, iyesno       = dummy variable
;       indexi1.             = position of maximum value of instrument
;
;integer, ianzahl_e, ianzahl_p, idim, iearcarrot, ihelcarrot, i1a_3_on
;integer, index_shift, minus, mode
;integer, iperihelion_shift, idummy, ivesno
;integer, index_i1a(1), index_i1b(1)
;
;       earcarlat,earcarlng = Carrington longitude and latitude from Earth
;       earsdissun          = distance Earth-Sun (always zero, historical)
;       helcarlat,helcarlng = Carrington longitude and latitude from Helios
;       heldissun           = distance Helios S/C to Sun
;       vnorm, vrad    = S/C tangential (in ecliptic), radial velocity(AU/day)
;       angle               = spacecraft-Sun-Earth angle
;       np                  = proton density [cm^{-3}]
;       vxs,vys,vzs         = center of mass velocity (1-fluid velocity)
;       x,y,z,zfp,zfe       = dummy parameters used to read in the 3-D VDF
;       sinfi, cosfi        = sine and cosine of azimuth of proton velocity
;       sintheta, costheta  = sine and cosine of elevation of proton velocity
;       background_i1a,b    = background of instruments (here the first two
;                             and last two channels)
;
;real, earcarlat, earcarlng, earsdissun
;real, dissun, carlat, carlng, vnorm, vrad, angle
;real, np, vxs, vys, vzs, x, y, z, zfp, zfe
;real, sintheta, costheta, sinfi, cosfi, background_i1a, background_i1b

;       b                     = 3 components of magnetic field [nT]
;       bsig                  = stanadrad deviation of the 3 components
;                               of magnetic field [nT]
;       i1aint                = 32 onboard integrated enerry channels of
;                               protons insrtument i1a
;       i1bint                = 32 onboard integrated energy channels of
;                               protons insrtument i1b
;       i1av, iabv            = corresponding velocities to i1aint, i1bint
;       i1abv, i1int          = sum and ratio of i1aint, i1bint (must be
;                               improved)
;       i1aalp                = 3 dimenional array: density [cm^{-3}],
;                                velocity [km/s] and temperature [K] of alphas
;	i1apro                = 5 dimenional array: density [cm^{-3}],
;                                velocity [km/s], temperature [K] of protons
;	i1bpro                = 3 dimenional array: density [cm^{-3}],
;                                velocity [km/s], temperature [K] of protons
;


b=dblarr(3)
bsig=dblarr(3)

i1aint=dblarr(32)
i1av=dblarr(32)
i1bint=dblarr(32)
i1bv=dblarr(32)
i1abv=dblarr(32)
i1int=dblarr(32)
i1apro=fltarr(5)
i1aalp=fltarr(3)
i1bpro=fltarr(3)

;---------------------- definitions ----------------------------------------
;
;       wxp, wyp, wzp      = velocity components corresponding
;	                        to respective distribution element
;                               1-dimensional arrays (protons) [km/s]
;       zp                 = value of distribution function
;                               1-dimensional arrays (protons)
;                       	[cm^[-6] s^[-3])]
;       wxe, wye, ze       = 2-D electron parameters (see above)
;
;
;-------------------------  open the file ----------------------------------
;read, 'filename'
;openr,1, 'filename'

;openr,1, 'P:\data-vdf\ions_Marsch_02\HELIOS\copy.dat'

;openr,1, 'E:\data-helios\2002-new\idl-vdf-02\h2y76d018h00m00s37_hdm.0
 openr,1, '/home/peters/helios/idl-vdf-02/h2y76d105h00m00s31_hdm.1

;       now read in the data. The first data set, the correction from
;       Pizzo et al. and for the finite solar diameter will not be used
;       (the first three records).
;
;	All other parameters are read in. For a definition see top.
;
;   A 'read' without a parameter reads of a 'comment' in the data file
;
print, 'enter reading'
;
;---------------------------------------------------------------------------

readf,1, comment
print, comment
readf,1, comment
print, comment
readf,1, comment
print, comment
readf,1, imode, ishift, iperihelion_shift, minus, i1a_3_on
print, imode, ishift, iperihelion_shift, minus, i1a_3_on
readf,1, heldissun, helcarlng, helcarlat, helcarrot
print, heldissun, helcarlng, helcarlat, helcarrot
readf,1, eardissun, earcarlng, earcarlat, angle, earcarrot
print, eardissun, earcarlng, earcarlat, angle, earcarrot
readf,1, vrad, vnorm
print, vrad, vnorm
readf,1, i1apro
print, i1apro
readf,1, i1aalp
print, i1aalp
readf,1, i1bpro
print, i1bpro
readf,1, b, bsig
print, b, bsig
readf,1, comment
print, comment
readf,1, i1aint
print, i1aint
readf,1, comment
print, comment
readf,1, i1av
print, i1av
readf,1, comment
print, comment
readf,1, i1bint
print, i1bint
readf,1, comment
print, comment
readf,1, i1bv
print, i1bv
readf,1, comment
print, comment
readf,1, comment
print, comment
readf,1, comment
print, comment


;---------------------------------------------------------------------------
;       Define velocity with respect to center of mass (1-D velocity).
;       Transform from degree to radian. Calculate center of mass
;       velocity 'v.s'
;
pi180 = 3.141592653589793238462643d0/180.d0
;print, pi180

        i1av = i1av - i1apro(1)
	i1bv = i1bv - i1apro(1)
	sinfi = sin(i1apro(3) * pi180)
	cosfi = cos(i1apro(3) * pi180)
	sintheta = sin(i1apro(4) * pi180)
	costheta = cos(i1apro(4) * pi180)
	vxs = i1apro(1) * cosfi * costheta
	vys = i1apro(1) * sinfi * costheta
	vzs = i1apro(1) * sintheta

print, 'solar wind (proton) velocity components'
print, vxs,vys,vzs


;---------------------------- Read VDFs --------------------------------------
;
;       now loop through the distribution fucntion. First read in the
;       the indices i,j,k of distribution function and the velocity
;       components. The dummy variables z,x,y,zfp,zfe contain the value of
;       the distribution function, and the respective velocity components.
;       The variable 'dummy' contains the count rate, which can be
;       used for statistics. The dummy variables are then
;       transformed in the center of mass system and written to the
;       arrays z,wvx,wvy,wvz. If an error occurs during read, the program
;       will continue at label electrons. The next data will be the
;       electrons, which are handled later.
;

;========================= protons ===========================================

print,'****************** protons **********************'

for mi=0,11 do begin
for mj=0,11 do begin
for mk=0,31 do begin
vdfp(mi,mj,mk) = 0.d0
vxp(mi,mj,mk) = 0.d0
vyp(mi,mj,mk) = 0.d0
vzp(mi,mj,mk) = 0.d0

vbp_ll(mi,mj,mk) = 0.d0
vbp_perp(mi,mj,mk) = 0.d0

countp(mi,mj,mk) = 0
endfor
endfor
endfor

for n_p=0, ndimp do begin
        on_ioerror, electrons      ; jump to electron data
        readf,1, i, j, k, zfp, dummy, x, y, z
        i=fix(i)
        j=fix(j)
        k=fix(k)
        dummy = fix(dummy)
        print, i, j, k, zfp, dummy, x, y, z
        ii= i-1
        jj= j-1
        kk= k-1
       print, ii, jj, kk
	zp(n_p) = zfp
	wxp(n_p) = x - vxs
	wyp(n_p) = y - vys
	wzp(n_p) = z - vzs

;bx=dblarr(1)
;by=dblarr(2)
;bz=dblarr(3)
;babs=(bx^2+by^2+bz^2)^0.5
;angx=bx/babs
;angy=by/babs
;angz=bz/babs
 ; vbl=wxp*angx+wyp*angy+wzp*angz
  ;  v2=wxp^2 +wyp^2 +wzp^2
   ; vbp=( v2-vbl^2)^0.5

;wp_ll(n_p)=vbl
;wp_perp(n_p)=vbp

       print, zp(n_p), wxp(n_p), wyp(n_p), wzp(n_p)
        vdfp(ii,jj,kk) = zfp
        vxp(ii,jj,kk) = x
        vyp(ii,jj,kk) = y
        vzp(ii,jj,kk) = z

    ;  b=dblarr(3)


;vbp_ll(ii,jj,kk) =vbl
;vbp_perp(ii,jj,kk) =vbp

        countp(ii,jj,kk)= dummy
endfor

;============================ electrons ===========================
;
;       Now read electron data. Procedure similar to that of protons.
;       First print one coment line, then loop: Structure: i,j indices
;       of the distribution function, zfe value of distribution function,
;       dummy is the count rate, x,y are the electron velocities

print,'************* electrons *************'
electrons: readf,1, comment
           print, comment
;           readf,1, comment
;           print, comment

for mi=0,7 do begin
for mj=0,15 do begin
vdfe(mi,mj) = 0.d0
vxe(mi,mj) = 0.d0
vye(mi,mj) = 0.d0
counte(mi,mj) = 0
endfor
endfor

for n_e=0, ndime do begin
    while not EOF(1) do begin
        readf,1, i, j, zfe, dummy, x, y
        i = fix(i)
        j = fix(j)
        dummy = fix(dummy)
        print, i, j, zfe, dummy, x, y
	ze(n_e) = zfe
        wxe(n_e) = x
        wye(n_e) = y
        ii = i - 1
        jj = j - 1
        vdfe(ii,jj) = zfe
        vxe(ii,jj) = x
        vye(ii,jj) = y
        counte(ii,jj) = dummy
    endwhile
endfor

close,1

;------------- determine the velocity ranges and maxima -----------------------
print, 'velocity ranges and maxima'

; protons

vdfpmax = max(vdfp)  ;imaxp,jmaxp,kmaxp
vxpmin = min(vxp)
vxpmax = max(vxp)
vypmin = min(vyp)
vypmax = max(vyp)
vzpmin = min(vzp)
vzpmax = max(vzp)
print, vdfpmax
;print, imaxp,jamxp,kmaxp
print, vxpmin,vxpmax,vypmin,vypmax,vzpmin,vzpmax

; electrons

vdfemax = max(vdfe)   ;imaxe,jmaxe
vxemin = min(vxe)
vxemax = max(vxe)
vyemin = min(vye)
vyemax = max(vye)
print, vdfemax
;print, imaxe,jamxe
print, vxemin,vxemax,vyemin,vyemax
;b=dblarr(3)

;set_plot,'win'
retain=2
;vdfp(ii,jj,kk), vxp(ii,jj,kk), vyp(ii,jj,kk), vzp(ii,jj,kk)
;contour,vdfp(*,*,15)
;xrange=[0.0,500.]
;!xmax=1000

;wp_ll
;wp_perp

plot,wyp,zp,psym=1
;plot,wxp,zp,psym=1

stop:
;---------------------------------------------------------------------------

end

;-------------------------- end programme -----------------------------------