| speed, size, contouring [message #55187] |
Wed, 01 August 2007 17:44 |
Hawaiianite
Messages: 2
Registered: August 2007
|
Junior Member |
|
|
Hello all,
Curious if you have some suggestions about a program I'm trying to
run. Basic scenerio is I have a large array (176851) of mineralogy (4
minerals) and their model spectra. I take one of them replicate it to
the size of the original array and then compare it to all of them at
once. I take the ones that are within a criteria I call excursion
(band with max difference between original spectrum and another model
spectrum) and then I see how close their mineralogy is to that
spectrum of interest. I want to do that for each composition/spectrum
I'm looking at and make a contour diagram of my results for each
mineral.
There are three areas in which I am struggling. Speed, size, and the
contouring feature. It takes about 190 seconds just to do 100 model
compositions (for some reason I've struggled to get around for loops
in this case). In regards to size some spectra often have many close
matches especially if I'm only trying to match 4 bands rather than
eight of the spectrum. I've found once I get about 65 matches the
array becomes too large for my computer to handle. I'm also
struggling with the putting the results into a contoured ternary
diagram (this I have done before, but with different data a long time
ago). For this I keep getting the error "CONTOUR: Contour levels must
be in increasing order." I've just never gotten this error before and
I'm having a little trouble deciding what it is referring to.
thanks,
Hawaiianite
pro minexcur,bands,n,m
;=========================================================== =====================
;WAVES--THE WAVELENGTHS CHOSEN TO BE MODELED (NEEDS TO BE DEFINED
BELOW)
waves=[0.75, 0.90, 0.95, 1.0,1.1, 1.25, 1.50, 2.0]
numwaves=(size(waves))(1)
;BANDS--NUMBER OF BANDS IN SPECTRUM WILL BE COMPARING
;N--REFERS TO A NUMBER OF MODELS TO COMPARE (TO GET PROGRAM UP AND
RUNNING AND THEN INCREASE SPEED)
;M--REFERS TO THE MINERAL YOU WILL BE USING
;MG--THE SPECIFIC MG-NUMBER ENTERED FOR THIS MODEL(NEEDS TO BE
DEFINED BELOW)
;VEC--IS A 12 ELEMENT VECTOR OF MINERALOGY AND CHEMISTRY.
t1=systime(1)
;=========================================================== =====================
;READ IN MODEL SPECTRA
;=========================================================== =====================
vec=12.
index=176851.
nummat=1.
mgs=1 ;HOW
MANY MG' MODELS TO LOAD
mginc=20 ;MG'
INCREMENT TO LOAD
mod_out=fltarr(vec+numwaves,nummat*index,mgs) ;OUT_VEC+ABS.
REFLECTANCE
for i=0,mgs-1 do begin
vecspec=fltarr(vec+numwaves,index*nummat)
filename2=string('outvs_min_07112007_',strtrim(55+i*mginc,2) )
openr,2,filename2,get_lun
readu,2,vecspec
mod_out[*,*,i]=vecspec[*,*]
close,/all
endfor
t2=systime(1)
print,systime(1)-t1,' Seconds'
;=========================================================== ====================
;=========================================================== ====================
;CALCULATIONS
;=========================================================== ====================
excur=fltarr(nummat*index)
a=fltarr(bands)
out_arr=fltarr(21,nummat*index,65)
av_arr=fltarr(1,nummat*index)
for j=0L,n-1 do begin
a[*]=mod_out[12:bands+11,j]
a[*]=a[*]/mean(a)
reparr=rebin(a,bands,nummat*index)
for i=0L,nummat*index-1 do begin
mod_out[12:bands+11,i]=mod_out[12:bands+11,i]/mean(mod_out[1 2:bands
+11,i])
excur[i]=max(abs(mod_out[12:bands+11,i]-
reparr[0:bands-1,i]))
endfor
d=where(excur gt .005 and excur lt .00501)
s=(size(d))(0)
p=(size(d))(1)
print,string('Model#= ',strtrim(j,2),' Min= ',strtrim(mod_out[0:3,j],
2),' # of matches = ',strtrim(p,2))
if (s gt 0) then begin
out_arr[0:19,j,0:p-1]=mod_out[*,d[*]]
out_arr[20,j,0:p-1]=excur[d[*]]
modminrep=mod_out[0:3,j]
modminrep=rebin(modminrep,4,p)
matches=fltarr(4,p)
for k=0,p-1 do matches[*,k]=out_arr[0:3,j,k]
mindif=fltarr(4,p)
for k=0,p-1 do mindif[*,k]=abs(modminrep[*,k]-matches[0:3,k])
av_arr[*,j]=mean(mindif[m,0:p-1])
endif
endfor
t3=systime(1)
print,systime(1)-t2,' Seconds'
;=========================================================== ====================
;=========================================================== ====================
;TRIANGULATE THE DATA POINTS FOR CONTOUR PROCEDURE.
;=========================================================== ====================
x=(mod_out[2,0:n-1]+mod_out[3,0:n-1])*100
y=(mod_out[3,0:n-1]*(sqrt(3)/2))*100
z=av_arr[*,0:n-1]*100
window,2,xsize=640,ysize=640
contour,z,x,y,LEVELS=fix(findgen(15) * ((max(z) - min(z))/15) +
min(z)) ,$
background=255,c_color=0,c_linestyle=2,c_labels=fix(findgen( 15) *
((max(z) - min(z))/15) + min(z)),$
c_charsize=1.5,ymargin=7,xmargin=7,c_thick=1,/irregular
image=tvrd()
contour,z,x,y,/fill,LEVELS=fix(findgen(15) * ((max(z) - min(z))/
15) + min(z)) ,$
background=255,c_labels=fix(findgen(15) * ((max(z) - min(z))/15) +
min(z)),c_charsize=1.25,ymargin=7,xmargin=7;,/irregular
image2=tvrd()
tvscl,image/2+image2/2
;DEFINE OUTLINE OF THE TERNARY.
;Define bottom of triangle.
x1=findgen (101.0)
y1=fltarr(101.0)
oplot, x1,y1, color=0,thick=2
;Define OPX to Pl side of triangle.
x2=findgen(51.0)
y2=1.732*x2
oplot, x2,y2, color=0,thick=2
;Define CPX to Pl side of triangle
x3=findgen (101.0)
y3=-1.732*x3+173.2
oplot, x3[50:100],y3[50:100], color =0,thick=2
;oplot, x,y,psym=5, color=0
;Define Rock Type Regions.
y4=fltarr(101.0)+77.9
oplot, x3[45:55],y4, color=0, thick=1,linestyle=0
y5=fltarr(101.0)+67.1
oplot, x3[39:61],y5, color=0, thick=1,linestyle=0
y6=fltarr(101.0)+52
oplot, x3[30:70],y6, color=0, thick=1,linestyle=0
y7=fltarr(101.0)+8.7
oplot, x3[5:95],y7, color=0, thick=1,linestyle=0
y8=findgen(101.0)
x4=fltarr(101.0)+50
oplot, x4, y8[9.0:78.0], color=0, thick=1, linestyle=0
y11=-1.732*x3+155.88
oplot, x3[60:85],y11[60:85], color =0,thick=1, linestyle=0
y12=1.732*x3-17.32
oplot, x3[15:40],y12[15:40], color =0,thick=1, linestyle=0
AXIS, YAXIS =0, color=255
AXIS, YAXIS =1, color=255
AXIS, XAXIS=0, color=0,charthick=1.5, charsize=1.8
AXIS, XAXIS=1, color=255
xyouts, 315, 550, '!8PL', /device, color=0,charsize=3,charthick=2
xyouts, 40, 10, '!8OPX', /device, color=0,charsize=3,charthick=2
xyouts, 560, 10, '!8CPX', /device, color=0,charsize=3,charthick=2
;=========================================================== ====================
stop
end
|
|
|
|
comp.lang.idl-pvwave archive
Members
Search
Help
Login
Home
