pro euv_dimming_cme1,dir,iev,vers,test ;+ ; Project : SDO/AIA, ; ; Name : EUV_DIMMING_CME ; ; Category : data analysis of EUV dimming and CME mass, speed, and energy ; ; Explanation : Reads AIA images and performs EUV dimming analysis ; ; Syntax : IDL>euv_dimming_disp,dir,iev,vers,test ; ; Inputs : dir,iev,vers,test ; ; Outputs ; IDL save file: *.sav ; ; History : 24-Feb-2016, Version 1 written by Markus J. Aschwanden ; ; Contact : aschwanden@lmsal.com ;- common dimming,model,time,a_t,v_t,x_t,q_t,q_obs,q_sig,h0,ind_fit print,'__________________EUV_DIMMING_CME________________________ string0,3,iev,iev_str catfile ='dimming_'+iev_str+vers+'.txt' readcol,dir+catfile,event_,dateobs_,noaa_,helpos_,instr_,$ w1_,w2_,w3_,w4_,w5_,w6_,w7_,w8_,fov_,$ format='(I,A,A,A,A,A,A,A,A,A,A,A,A,A)',/silent nt =n_elements(event_) ;________________SAVE PARAMETERS______________________________ savefile='dimming_'+iev_str+vers+'.sav' restore,dir+savefile ;-->input,para,exptime_,flux_grid,x_grid,$ ;y_grid,flux_prof,emtot_prof,telog_prof,tsig_prof,$ ;chi_prof,wave_,time,dem_grid,dem_t,te_mk,telog,em_xy print,'read parameter file = ',savefile dim =size(flux_grid) nx =dim(1) ny =dim(2) nwave =dim(3) helpos =input.helpos nbin =input.nbin lon =long(strmid(helpos,4,2)) lat =long(strmid(helpos,1,2)) l =lon*!pi/180. b =lat*!pi/180. angle =sqrt(l^2+b^2) r0 =6.96e10 ;cm cadence =input.cadence ;_________________SELECTED LOCATION OF MAXIMUM DIMMING________________ chi_best=9999. q_peak =0.5 chi_mask=fltarr(nx,ny)+9999. dem_max =max(dem_grid) n =nx*ny for k=0,n do begin i =(k mod nx) j =k/nx if (k le n-1) then begin if (max(dem_grid(i,j,*)) lt dem_max*q_peak) then goto,skip_fit dem_xy=reform(dem_grid(i,j,*)) ;single macropixel profile endif if (k eq n) then begin dem_xy=fltarr(nt) chi_min=min(chi_mask) for ii=0,nx-1 do begin for jj=0,ny-1 do begin dem_ij=reform(dem_grid(ii,jj,*)) ;single macropixel profile weight=1./(chi_mask(ii,jj)/chi_min)^2 dem_xy=dem_xy+dem_ij*weight ;weighted emission measure profile endfor endfor endif ;__________________INTERPOLATE LOW_EM OUTLIERS__________________ q_out =0.5 for it=1,nt-2 do begin em_avg=(dem_xy(it-1)+dem_xy(it+1))/2. if (dem_xy(it) lt em_avg/2.) then dem_xy(it)=em_avg endfor ;___________________SELECT FITTING RANGE________________________ em_max=max(dem_xy,it_max) &it_max=(it_max > 1) < (nt-2) em_min=min(dem_xy(it_max:nt-1),di) it_min=it_max+di tfit =900. ;min 0.25 hours nfit =long(it_min-it_max+1) > (tfit/cadence) ind_fit=long(it_max+findgen(nfit)) q_obs =(dem_xy-em_min)/(em_max-em_min) em_obs =em_min+(em_max-em_min)*q_obs sig =q_obs(it_max:nt-2)-(q_obs(it_max-1:nt-3)+q_obs(it_max+1:nt-1))/2. q_sig =fltarr(nt)+median(abs(sig)) ;_________________CME VELOCITY AND ACCELERATION___________ model =4 dt =time(1)-time(0) t0 =time(it_max) ;s it =findgen(nt) tau =1000. ;s h0 =4.7e9*te_mk ;[cm] density scale height d0 =(r0+h0)*sin(angle) ;distance from disk center v0 =3.0e7 ;cm/s a0 =0.1*v0/dt qmin =0.01 coeff =[a0/1.e5,tau/1.e3,t0/1.e3,qmin] ncoeff =n_elements(coeff) xi =fltarr(ncoeff,ncoeff) for m=0,ncoeff-1 do xi(m,m)=1.0 ftol =1.e-4 powell,coeff,xi,ftol,fmin,'euv_dimming_powell' chi =euv_dimming_powell(coeff) ;-->a_model,v_model,x_model,q_t if (k le n-1) then chi_mask(i,j)=chi SKIP_FIT: endfor ;____________________BEST SOLUTION____________________________ chi =euv_dimming_powell(coeff) ;-->a_model,v_model,x_model,q_t acc =coeff(0)*1.e5 tau =coeff(1)*1.e3 t0 =coeff(2)*1.e3 ;s d0 =(r0+h0)*sin(angle) ;distance from disk center n_model =1 ;4 model functions a_model =fltarr(nt,n_model) v_model =fltarr(nt,n_model) x_model =fltarr(nt,n_model) em_model =fltarr(nt,n_model) a0_model =fltarr(n_model) t0_model =fltarr(n_model) ta_model =fltarr(n_model) qe_model =fltarr(n_model) ekin_model=fltarr(n_model) chi_model =fltarr(n_model) x_lasco =fltarr(n_model) v_lasco =fltarr(n_model) t_lasco =fltarr(n_model) m =0 a_model(*,m)=a_t v_model(*,m)=v_t x_model(*,m)=x_t+d0 em_model(*,m)=em_min+(em_max-em_min)*q_t a0_model(m)=acc t0_model(m)=t0 ta_model(m)=(t0+tau) ;s chi_model(m)=chi ;__________________LASCO/C2 TIME AND VELOCITY_____________________ lasco_file='~/work_global/ecme_lasco.dat' readcol,lasco_file,flare,dateobs_lasco,start_goes,start_lasco,dt_lasco,h_lasco,$ v1_lasco,v2_lasco,v3_lasco,m_lasco,e_lasco,skipline=3,format='(I,A,A,A,F,F,I,I,I,F,F)' ind0 =where(flare eq iev,n) i =ind0(0) d_lasco =h_lasco(i) x_lasco =r0*h_lasco(i) ;cm for m=0,n_model-1 do begin ind1 =where(x_model(*,m) ge x_lasco,n_lasco) if (max(x_model(*,m)) lt x_lasco) then ind1=[nt-1] v_lasco(m)=v_model(ind1(0),m) ;cm/s t_lasco(m)=time(ind1(0)) ;s print,m,v_lasco(m)/1.e5,' km/s' endfor ;_________________DIMMING AREA_________________________________ em_tot =em_max-em_min e1_xy =em_xy(*,*,it_max) e2_xy =em_xy(*,*,nt-1) de_med =median(e1_xy-e2_xy) e3_xy =e1_xy-de_med*2. ind_dimm=where(e2_xy lt e3_xy,npix_dimm) q_dimm =float(npix_dimm)/float(nx*ny) wid_dimm=sqrt(float(npix_dimm)) ;width in units of grid bins dx_grid_cm=r0*(x_grid(1)-x_grid(0)) ;grid scale in units of cm l_cm =wid_dimm*dx_grid_cm ;cm area =l_cm^2 ;equivalent area [cm>2] vol =l_cm^3 ;volume [cm3] angle =90.*(l_cm/(r0*!pi/2.)) ;angular width [deg] ;_________________DIMMING MASS AND KINETIC ENERGY________________ nel =sqrt(em_tot/l_cm) mp =1.67e-24 ;proton mass mass =mp*nel*vol ekin_model=(1./2.)*mass*v_lasco^2 ;_________________GRAVITATIONAL ENERGY________________________ g =6.672e-8 ;dyne cm&2 g^-2 gravitational constant m_sun =1.99e33 ;g solar mass e_grav =g*(m_sun/r0)*mass ;gravitational energy v_esc =sqrt(2.*g*(m_sun/r0)) ;escape velocity [cm/s] ;________________SAVE PARAMETERS______________________________ save,filename=dir+savefile,input,para,exptime_,flux_grid,x_grid,y_grid,$ flux_prof,emtot_prof,telog_prof,tsig_prof,chi_prof,wave_,time,$ dem_grid,dem_t,te_mk,telog,em_xy,$ l_cm,angle,area,vol,em_tot,dem_t,nel,te_mk,telog,mass,em_model,$ d0,x_lasco,v_lasco,t_lasco,a_model,v_model,x_model,ekin_model,$ q_t,q_obs,em_obs,a0_model,t0_model,ta_model,chi_model,ind_fit,h0,$ e_grav,v_esc print,'parameters saved in file = ',savefile end