(* ::Package:: *)

BeginPackage[ "CORDS`",{"ErrorBarPlots`"}]

Cords::usage = 
	"Code for Originslab Routine Double-Spike (CORDS). "

Begin["Private`"]

Cords[inib_:1.1,itnb_:3]:=Module[{afrisotope,anmisotope,blist,blklist,blknb,colnb,coltitle,corcup,cupconfig,defparameter,exppath,frisotope,frmass,frsp,frstd,frtitle,i,intcoef,intcup,intnb,inttitle,isocup,isocuporder,isomass,isoname,isonb,isotitle,j,l,maincolnb,maincoltitle,massnb,massratio,mratio,mscoef,nmisotope,nmmass,nmnb,nmsp,nmstd,nmtitle,norisotope,offcup,optol,optprm,prmtable,rationb,ratiotitle,settings,simblist,simnb,smplist,smpnb,spcomp,spratio,stdcomp,stdratio,subcolnb,subcupnb,summary,sumtitle1,sumtitle2,ver,workparameter,extparameter},
ver="V 1.0.1";
(* ZoomGraphics function is from Szabolcs Horv\[AAcute]t <szhorvat@gmail.com> *)
(* Http://forums.wolfram.com/mathgroup/archive/2008/Jan/msg00009.html *)
ZoomGraphics[graph_Graphics]:=With[{gr=First[graph],opt=DeleteCases[Options[graph],PlotRange->_],plr=PlotRange/.Options[graph,PlotRange],rectangle={Dashing[Small],Line[{#1,{First[#2],Last[#1]},#2,{First[#1],Last[#2]},#1}]}&},DynamicModule[{dragging=False,first,second,range=plr},Panel@EventHandler[Dynamic@Graphics[If[dragging,{gr,rectangle[first,second]},gr],PlotRange->range,Sequence@@opt],{{"MouseDown",1}:>(first=MousePosition["Graphics"]),{"MouseDragged",1}:>(dragging=True;
second=MousePosition["Graphics"]),{"MouseUp",1}:>If[dragging,dragging=False;
range=Transpose@{first,second},range=plr]}]]];
(************************************************************** Default Parameters ************************************************************)
defparameter={{"Mass Number","*","Element","","","","","","","","*","Atomic mass","","","","","","","","*","Natural abundance %","","","","","","","","*"},{"Pos. ","*",1.`,2.`,3.`,4.`,5.`,6.`,7.`,8.`,"*",1.`,2.`,3.`,4.`,5.`,6.`,7.`,8.`,"*",1.`,2.`,3.`,4.`,5.`,6.`,7.`,8.`,"*"},{1.`,"*","H","","","","","","","","*",1.007825035`,"","","","","","","","*",99.9885`,"","","","","","","","*"},{2.`,"*","H","","","","","","","","*",2.014101779`,"","","","","","","","*",0.0115`,"","","","","","","","*"},{3.`,"*","","He","","","","","","","*","",3.01602931`,"","","","","","","*","",0.000137`,"","","","","","","*"},{4.`,"*","","He","","","","","","","*","",4.00260324`,"","","","","","","*","",99.999863`,"","","","","","","*"},{5.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{6.`,"*","","","Li","","","","","","*","","",6.0151214`,"","","","","","*","","",7.59`,"","","","","","*"},{7.`,"*","","","Li","","","","","","*","","",7.016003`,"","","","","","*","","",92.41`,"","","","","","*"},{8.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{9.`,"*","","","","Be","","","","","*","","","",9.0121822`,"","","","","*","","","",100.`,"","","","","*"},{10.`,"*","","","","","B","","","","*","","","","",10.0129369`,"","","","*","","","","",19.9`,"","","","*"},{11.`,"*","","","","","B","","","","*","","","","",11.0093054`,"","","","*","","","","",80.1`,"","","","*"},{12.`,"*","","","","","","C","","","*","","","","","",12.`,"","","*","","","","","",98.93`,"","","*"},{13.`,"*","","","","","","C","","","*","","","","","",13.003354826`,"","","*","","","","","",1.07`,"","","*"},{14.`,"*","","","","","","","N","","*","","","","","","",14.003074002`,"","*","","","","","","",99.632`,"","*"},{15.`,"*","","","","","","","N","","*","","","","","","",15.00010897`,"","*","","","","","","",0.368`,"","*"},{16.`,"*","","","","","","","","O","*","","","","","","","",15.99491463`,"*","","","","","","","",99.757`,"*"},{17.`,"*","","","","","","","","O","*","","","","","","","",16.9991312`,"*","","","","","","","",0.038`,"*"},{18.`,"*","","","","","","","","O","*","","","","","","","",17.9991603`,"*","","","","","","","",0.205`,"*"},{19.`,"*","F","","","","","","","","*",18.99840322`,"","","","","","","","*",100.`,"","","","","","","","*"},{20.`,"*","","Ne","","","","","","","*","",19.9924356`,"","","","","","","*","",90.48`,"","","","","","","*"},{21.`,"*","","Ne","","","","","","","*","",20.9938428`,"","","","","","","*","",0.27`,"","","","","","","*"},{22.`,"*","","Ne","","","","","","","*","",21.9913831`,"","","","","","","*","",9.25`,"","","","","","","*"},{23.`,"*","","","Na","","","","","","*","","",22.9897677`,"","","","","","*","","",100.`,"","","","","","*"},{24.`,"*","","","","Mg","","","","","*","","","",23.9850423`,"","","","","*","","","","","","","","","*"},{25.`,"*","","","","Mg","","","","","*","","","",24.9858374`,"","","","","*","","","","","","","","","*"},{26.`,"*","","","","Mg","","","","","*","","","",25.9825937`,"","","","","*","","","",78.99`,"","","","","*"},{27.`,"*","","","","","Al","","","","*","","","","",26.9815386`,"","","","*","","","",10.`,100.`,"","","","*"},{28.`,"*","","","","","","Si","","","*","","","","","",27.9769271`,"","","*","","","",11.01`,"",92.2297`,"","","*"},{29.`,"*","","","","","","Si","","","*","","","","","",28.9764949`,"","","*","","","","","",4.6832`,"","","*"},{30.`,"*","","","","","","Si","","","*","","","","","",29.9737707`,"","","*","","","","","",3.0872`,"","","*"},{31.`,"*","","","","","","","P","","*","","","","","","",30.973762`,"","*","","","","","","",100.`,"","*"},{32.`,"*","","","","","","","","S","*","","","","","","","",31.9720707`,"*","","","","","","","",94.93`,"*"},{33.`,"*","","","","","","","","S","*","","","","","","","",32.97145843`,"*","","","","","","","",0.76`,"*"},{34.`,"*","","","","","","","","S","*","","","","","","","",33.96786665`,"*","","","","","","","",4.29`,"*"},{35.`,"*","Cl","","","","","","","","*",34.968852721`,"","","","","","","","*",75.78`,"","","","","","","","*"},{36.`,"*","","Ar","","","","","","S","*","",35.96754552`,"","","","","",35.96708062`,"*","",0.3365`,"","","","","",0.02`,"*"},{37.`,"*","Cl","","","","","","","","*",36.96590262`,"","","","","","","","*",24.22`,"","","","","","","","*"},{38.`,"*","","Ar","","","","","","","*","",37.9627325`,"","","","","","","*","",0.0632`,93.2581`,"","","","","","*"},{39.`,"*","","","K","","","","","","*","","",38.9637074`,"","","","","","*","","",0.0117`,"","","","","","*"},{40.`,"*","","Ar","K","Ca","","","","","*","",39.9623837`,39.9639992`,39.9625906`,"","","","","*","",99.6003`,6.7302`,96.941`,"","","","","*"},{41.`,"*","","","K","","","","","","*","","",40.9618254`,"","","","","","*","","","","","","","","","*"},{42.`,"*","","","","Ca","","","","","*","","","",41.9586176`,"","","","","*","","","",0.647`,"","","","","*"},{43.`,"*","","","","Ca","","","","","*","","","",42.9587662`,"","","","","*","","","",0.135`,"","","","","*"},{44.`,"*","","","","Ca","","","","","*","","","",43.9554806`,"","","","","*","","","",2.086`,"","","","","*"},{45.`,"*","","","","","Sc","","","","*","","","","",44.95591`,"","","","*","","","","",100.`,"","","","*"},{46.`,"*","","","","Ca","","Ti","","","*","","","",45.953689`,"",45.9526294`,"","","*","","","",0.004`,"",8.25`,"","","*"},{47.`,"*","","","","","","Ti","","","*","","","","","",46.951764`,"","","*","","","","","",7.44`,"","","*"},{48.`,"*","","","","Ca","","Ti","","","*","","","",47.952533`,"",47.9479473`,"","","*","","","",0.187`,"",73.72`,"","","*"},{49.`,"*","","","","","","Ti","","","*","","","","","",48.9478711`,"","","*","","","","","",5.41`,"","","*"},{50.`,"*","","","","","","Ti","V","Cr","*","","","","","",49.9447921`,49.9471609`,49.9460464`,"*","","","","","",5.18`,0.25`,4.345`,"*"},{51.`,"*","","","","","","","V","","*","","","","","","",50.9439617`,"","*","","","","","","",99.75`,"","*"},{52.`,"*","","","","","","","","Cr","*","","","","","","","",51.9405098`,"*","","","","","","","",83.789`,"*"},{53.`,"*","","","","","","","","Cr","*","","","","","","","",52.9406513`,"*","","","","","","","",9.501`,"*"},{54.`,"*","","Fe","","","","","","Cr","*","",53.9396127`,"","","","","",53.9388825`,"*","",5.845`,"","","","","",2.365`,"*"},{55.`,"*","Mn","","","","","","","","*",54.9380471`,"","","","","","","","*",100.`,"","","","","","","","*"},{56.`,"*","","Fe","","","","","","","*","",55.9349393`,"","","","","","","*","",91.754`,"","","","","","","*"},{57.`,"*","","Fe","","","","","","","*","",56.9353958`,"","","","","","","*","",2.119`,"","","","","","","*"},{58.`,"*","","Fe","","Ni","","","","","*","",57.9332773`,"",57.9353462`,"","","","","*","",0.282`,"",68.0769`,"","","","","*"},{59.`,"*","","","Co","","","","","","*","","",58.9331976`,"","","","","","*","","",100.`,"","","","","","*"},{60.`,"*","","","","Ni","","","","","*","","","",59.9307884`,"","","","","*","","","",26.2231`,"","","","","*"},{61.`,"*","","","","Ni","","","","","*","","","",60.9310579`,"","","","","*","","","",1.1399`,"","","","","*"},{62.`,"*","","","","Ni","","","","","*","","","",61.9283461`,"","","","","*","","","",3.6345`,"","","","","*"},{63.`,"*","","","","","Cu","","","","*","","","","",62.9295989`,"","","","*","","","","",69.17`,"","","","*"},{64.`,"*","","","","Ni","","Zn","","","*","","","",63.9279679`,"",63.9291448`,"","","*","","","",0.9256`,"",48.63`,"","","*"},{65.`,"*","","","","","Cu","","","","*","","","","",64.9277929`,"","","","*","","","","",30.83`,"","","","*"},{66.`,"*","","","","","","Zn","","","*","","","","","",65.9260347`,"","","*","","","","","",27.9`,"","","*"},{67.`,"*","","","","","","Zn","","","*","","","","","",66.9271291`,"","","*","","","","","",4.1`,"","","*"},{68.`,"*","","","","","","Zn","","","*","","","","","",67.9248459`,"","","*","","","","","",18.75`,"","","*"},{69.`,"*","","","","","","","Ga","","*","","","","","","",68.92558`,"","*","","","","","","",60.108`,"","*"},{70.`,"*","","","","","","Zn","","Ge","*","","","","","",69.925325`,"",69.9242497`,"*","","","","","",0.62`,"",20.84`,"*"},{71.`,"*","","","","","","","Ga","","*","","","","","","",70.9247005`,"","*","","","","","","",39.892`,"","*"},{72.`,"*","","","","","","","","Ge","*","","","","","","","",71.9220789`,"*","","","","","","","",27.54`,"*"},{73.`,"*","","","","","","","","Ge","*","","","","","","","",72.9234626`,"*","","","","","","","",7.73`,"*"},{74.`,"*","","Se","","","","","","Ge","*","",73.9224746`,"","","","","",73.9211774`,"*","",0.89`,"","","","","",36.28`,"*"},{75.`,"*","As","","","","","","","","*",74.9215942`,"","","","","","","","*",100.`,"","","","","","","","*"},{76.`,"*","","Se","","","","","","Ge","*","",75.919212`,"","","","","",75.9214016`,"*","",9.37`,"","","","","",7.61`,"*"},{77.`,"*","","Se","","","","","","","*","",76.9199125`,"","","","","","","*","",7.63`,"","","","","","","*"},{78.`,"*","","Se","","Kr","","","","","*","",77.9173076`,"",77.920396`,"","","","","*","",23.77`,"",0.35`,"","","","","*"},{79.`,"*","","","Br","","","","","","*","","",78.9183361`,"","","","","","*","","",50.69`,"","","","","","*"},{80.`,"*","","Se","","Kr","","","","","*","",79.9165196`,"",79.91638`,"","","","","*","",49.61`,"",2.28`,"","","","","*"},{81.`,"*","","","Br","","","","","","*","","",80.916289`,"","","","","","*","","",49.31`,"","","","","","*"},{82.`,"*","","Se","","Kr","","","","","*","",81.9166978`,"",81.913482`,"","","","","*","",8.73`,"",11.58`,"","","","","*"},{83.`,"*","","","","Kr","","","","","*","","","",82.914135`,"","","","","*","","","",11.49`,"","","","","*"},{84.`,"*","","","","Kr","","Sr","","","*","","","",83.911507`,"",83.91343`,"","","*","","","",57.`,"",0.56`,"","","*"},{85.`,"*","","","","","Rb","","","","*","","","","",84.911794`,"","","","*","","","","",72.17`,"","","","*"},{86.`,"*","","","","Kr","","Sr","","","*","","","",85.910616`,"",85.9092672`,"","","*","","","",17.3`,"",9.86`,"","","*"},{87.`,"*","","","","","Rb","Sr","","","*","","","","",86.909187`,86.9088841`,"","","*","","","","",27.83`,7.`,"","","*"},{88.`,"*","","","","","","Sr","","","*","","","","","",87.9056188`,"","","*","","","","","",82.58`,"","","*"},{89.`,"*","","","","","","","Y","","*","","","","","","",89.905849`,"","*","","","","","","",100.`,"","*"},{90.`,"*","","","","","","","","Zr","*","","","","","","","",89.9047026`,"*","","","","","","","",51.45`,"*"},{91.`,"*","","","","","","","","Zr","*","","","","","","","",90.9056439`,"*","","","","","","","",11.22`,"*"},{92.`,"*","","Mo","","","","","","Zr","*","",91.906809`,"","","","","",91.9050386`,"*","",14.84`,"","","","","",17.15`,"*"},{93.`,"*","Nb","","","","","","","","*",92.9063772`,"","","","","","","","*",100.`,"","","","","","","","*"},{94.`,"*","","Mo","","","","","","Zr","*","",93.9050853`,"","","","","",93.9063148`,"*","",9.25`,"","","","","",17.38`,"*"},{95.`,"*","","Mo","","","","","","","*","",94.9058411`,"","","","","","","*","",15.92`,"","","","","","","*"},{96.`,"*","","Mo","Ru","","","","","Zr","*","",95.9046785`,95.907599`,"","","","",95.908275`,"*","",16.68`,5.54`,"","","","",2.8`,"*"},{97.`,"*","","Mo","","","","","","","*","",96.9060205`,"","","","","","","*","",9.55`,"","","","","","","*"},{98.`,"*","","Mo","Ru","","","","","","*","",97.9054073`,97.905287`,"","","","","","*","",24.13`,1.87`,"","","","","","*"},{99.`,"*","","","Ru","","","","","","*","","",98.9059389`,"","","","","","*","","",12.76`,"","","","","","*"},{100.`,"*","","Mo","Ru","","","","","","*","",99.907477`,99.9042192`,"","","","","","*","",9.63`,12.6`,"","","","","","*"},{101.`,"*","","","Ru","","","","","","*","","",100.9055819`,"","","","","","*","","",17.06`,"","","","","","*"},{102.`,"*","","","Ru","","Pd","","","","*","","",101.9043485`,"",101.905634`,"","","","*","","",31.55`,"",1.02`,"","","","*"},{103.`,"*","","","","Rh","","","","","*","","","",102.9055`,"","","","","*","","","",100.`,"","","","","*"},{104.`,"*","","","Ru","","Pd","","","","*","","",103.905424`,"",103.904029`,"","","","*","","",18.62`,"",11.14`,"","","","*"},{105.`,"*","","","","","Pd","","","","*","","","","",104.905079`,"","","","*","","","","",22.33`,"","","","*"},{106.`,"*","","","","","Pd","","Cd","","*","","","","",105.903478`,"",105.906461`,"","*","","","","",27.33`,"",1.25`,"","*"},{107.`,"*","","","","","","Ag","","","*","","","","","",106.905092`,"","","*","","","","","",51.839`,"","","*"},{108.`,"*","","","","","Pd","","Cd","","*","","","","",107.903895`,"",107.904176`,"","*","","","","",26.46`,"",0.89`,"","*"},{109.`,"*","","","","","","Ag","","","*","","","","","",108.904756`,"","","*","","","","","",48.161`,"","","*"},{110.`,"*","","","","","Pd","","Cd","","*","","","","",109.905167`,"",109.903005`,"","*","","","","",11.72`,"",12.49`,"","*"},{111.`,"*","","","","","","","Cd","","*","","","","","","",110.904182`,"","*","","","","","","",12.8`,"","*"},{112.`,"*","Sn","","","","","","Cd","","*",111.904826`,"","","","","",111.902757`,"","*",0.97`,"","","","","",24.13`,"","*"},{113.`,"*","","","","","","","Cd","In","*","","","","","","",112.9044`,112.904061`,"*","","","","","","",12.22`,4.29`,"*"},{114.`,"*","Sn","","","","","","Cd","","*",113.902784`,"","","","","",113.903357`,"","*",0.66`,"","","","","",28.73`,"","*"},{115.`,"*","Sn","","","","","","","In","*",114.903348`,"","","","","","",114.903882`,"*",0.34`,"","","","","","",95.71`,"*"},{116.`,"*","Sn","","","","","","Cd","","*",115.901747`,"","","","","",115.904755`,"","*",14.54`,"","","","","",7.49`,"","*"},{117.`,"*","Sn","","","","","","","","*",116.902956`,"","","","","","","","*",7.68`,"","","","","","","","*"},{118.`,"*","Sn","","","","","","","","*",117.901609`,"","","","","","","","*",24.22`,"","","","","","","","*"},{119.`,"*","Sn","","","","","","","","*",118.903311`,"","","","","","","","*",8.59`,"","","","","","","","*"},{120.`,"*","Sn","","Te","","","","","","*",119.9021991`,"",119.904048`,"","","","","","*",32.58`,"",0.09`,"","","","","","*"},{121.`,"*","","Sb","","","","","","","*","",120.9038212`,"","","","","","","*","",57.21`,"","","","","","","*"},{122.`,"*","Sn","","Te","","","","","","*",121.9034404`,"",121.90305`,"","","","","","*",4.63`,"",2.55`,"","","","","","*"},{123.`,"*","","Sb","Te","","","","","","*","",122.904216`,122.904271`,"","","","","","*","",42.79`,0.89`,"","","","","","*"},{124.`,"*","Sn","","Te","","Xe","","","","*",123.9052743`,"",123.902818`,"",123.9058942`,"","","","*",5.79`,"",4.74`,"",0.09`,"","","","*"},{125.`,"*","","","Te","","","","","","*","","",124.9044285`,"","","","","","*","","",7.07`,"","","","","","*"},{126.`,"*","","","Te","","Xe","","","","*","","",125.9033095`,"",125.904281`,"","","","*","","",18.84`,"",0.09`,"","","","*"},{127.`,"*","","","","I","","","","","*","","","",126.904473`,"","","","","*","","","",100.`,"","","","","*"},{128.`,"*","","","Te","","Xe","","","","*","","",127.904463`,"",127.9035312`,"","","","*","","",31.74`,"",1.92`,"","","","*"},{129.`,"*","","","","","Xe","","","","*","","","","",128.9047801`,"","","","*","","","","",26.44`,"","","","*"},{130.`,"*","","","Te","","Xe","","Ba","","*","","",129.906229`,"",129.9035094`,"",129.906282`,"","*","","",34.08`,"",4.08`,"",0.106`,"","*"},{131.`,"*","","","","","Xe","","","","*","","","","",130.905072`,"","","","*","","","","",21.18`,"","","","*"},{132.`,"*","","","","","Xe","","Ba","","*","","","","",131.904144`,"",131.905042`,"","*","","","","",26.89`,"",0.101`,"","*"},{133.`,"*","","","","","","Cs","","","*","","","","","",132.905429`,"","","*","","","","","",100.`,"","","*"},{134.`,"*","","","","","Xe","","Ba","","*","","","","",133.905395`,"",133.904486`,"","*","","","","",10.44`,"",2.417`,"","*"},{135.`,"*","","","","","","","Ba","","*","","","","","","",134.905665`,"","*","","","","","","",6.592`,"","*"},{136.`,"*","Ce","","","","Xe","","Ba","","*",135.90714`,"","","",135.907214`,"",135.904553`,"","*",0.185`,"","","",8.87`,"",7.854`,"","*"},{137.`,"*","","","","","","","Ba","","*","","","","","","",136.905812`,"","*","","","","","","",11.232`,"","*"},{138.`,"*","Ce","","","","","","Ba","La","*",137.905985`,"","","","","",137.905232`,137.907105`,"*",0.251`,"","","","","",71.698`,0.09`,"*"},{139.`,"*","","","","","","","","La","*","","","","","","","",138.906347`,"*","","","","","","","",99.91`,"*"},{140.`,"*","Ce","","","","","","","","*",139.905433`,"","","","","","","","*",88.45`,"","","","","","","","*"},{141.`,"*","","Pr","","","","","","","*","",140.907647`,"","","","","","","*","",100.`,"","","","","","","*"},{142.`,"*","Ce","","Nd","","","","","","*",141.909241`,"",141.907719`,"","","","","","*",11.114`,"",27.151999999999997`,"","","","","","*"},{143.`,"*","","","Nd","","","","","","*","","",142.90981`,"","","","","","*","","",12.174`,"","","","","","*"},{144.`,"*","","","Nd","Sm","","","","","*","","",143.910083`,143.911998`,"","","","","*","","",23.798`,3.07`,"","","","","*"},{145.`,"*","","","Nd","","","","","","*","","",144.912569`,"","","","","","*","","",8.293000000000001`,"","","","","","*"},{146.`,"*","","","Nd","","","","","","*","","",145.913112`,"","","","","","*","","",17.189`,"","","","","","*"},{147.`,"*","","","","Sm","","","","","*","","","",146.914894`,"","","","","*","","","",14.99`,"","","","","*"},{148.`,"*","","","Nd","Sm","","","","","*","","",147.916889`,147.914819`,"","","","","*","","",5.756`,11.24`,"","","","","*"},{149.`,"*","","","","Sm","","","","","*","","","",148.91718`,"","","","","*","","","",13.82`,"","","","","*"},{150.`,"*","","","Nd","Sm","","","","","*","","",149.920887`,149.917273`,"","","","","*","","",5.638`,7.38`,"","","","","*"},{151.`,"*","","","","","Eu","","","","*","","","","",150.919702`,"","","","*","","","","",47.81`,"","","","*"},{152.`,"*","","","","Sm","","Gd","","","*","","","",151.919728`,"",151.919786`,"","","*","","","",26.75`,"",0.2`,"","","*"},{153.`,"*","","","","","Eu","","","","*","","","","",152.921225`,"","","","*","","","","",52.19`,"","","","*"},{154.`,"*","","","","Sm","","Gd","","","*","","","",153.922205`,"",153.920861`,"","","*","","","",22.75`,"",2.18`,"","","*"},{155.`,"*","","","","","","Gd","","","*","","","","","",154.922618`,"","","*","","","","","",14.8`,"","","*"},{156.`,"*","","","","","","Gd","","Dy","*","","","","","",155.922118`,"",155.924277`,"*","","","","","",20.47`,"",0.06`,"*"},{157.`,"*","","","","","","Gd","","","*","","","","","",156.923956`,"","","*","","","","","",15.65`,"","","*"},{158.`,"*","","","","","","Gd","","Dy","*","","","","","",157.924019`,"",157.924403`,"*","","","","","",24.84`,"",0.1`,"*"},{159.`,"*","","","","","","","Tb","","*","","","","","","",158.925342`,"","*","","","","","","",100.`,"","*"},{160.`,"*","","","","","","Gd","","Dy","*","","","","","",159.927049`,"",159.925193`,"*","","","","","",21.86`,"",2.34`,"*"},{161.`,"*","","","","","","","","Dy","*","","","","","","","",160.92693`,"*","","","","","","","",18.91`,"*"},{162.`,"*","","Er","","","","","","Dy","*","",161.928775`,"","","","","",161.926795`,"*","",0.14`,"","","","","",25.51`,"*"},{163.`,"*","","","","","","","","Dy","*","","","","","","","",162.928728`,"*","","","","","","","",24.9`,"*"},{164.`,"*","","Er","","","","","","Dy","*","",163.929198`,"","","","","",163.929171`,"*","",1.61`,"","","","","",28.18`,"*"},{165.`,"*","Ho","","","","","","","","*",164.930319`,"","","","","","","","*",100.`,"","","","","","","","*"},{166.`,"*","","Er","","","","","","","*","",165.93029`,"","","","","","","*","",33.61`,"","","","","","","*"},{167.`,"*","","Er","","","","","","","*","",166.932046`,"","","","","","","*","",22.93`,"","","","","","","*"},{168.`,"*","","Er","","Yb","","","","","*","",167.932368`,"",167.933894`,"","","","","*","",26.78`,"",0.13`,"","","","","*"},{169.`,"*","","","Tm","","","","","","*","","",168.934212`,"","","","","","*","","",100.`,"","","","","","*"},{170.`,"*","","Er","","Yb","","","","","*","",169.935461`,"",169.934759`,"","","","","*","",14.93`,"",3.04`,"","","","","*"},{171.`,"*","","","","Yb","","","","","*","","","",170.936323`,"","","","","*","","","",14.28`,"","","","","*"},{172.`,"*","","","","Yb","","","","","*","","","",171.936378`,"","","","","*","","","",21.83`,"","","","","*"},{173.`,"*","","","","Yb","","","","","*","","","",172.938208`,"","","","","*","","","",16.13`,"","","","","*"},{174.`,"*","","","","Yb","","Hf","","","*","","","",173.938859`,"",173.940044`,"","","*","","","",31.83`,"",0.16`,"","","*"},{175.`,"*","","","","","Lu","","","","*","","","","",174.94077`,"","","","*","","","","",97.41`,"","","","*"},{176.`,"*","","","","Yb","Lu","Hf","","","*","","","",175.942564`,175.942679`,175.941406`,"","","*","","","",12.76`,2.59`,5.26`,"","","*"},{177.`,"*","","","","","","Hf","","","*","","","","","",176.943217`,"","","*","","","","","",18.6`,"","","*"},{178.`,"*","","","","","","Hf","","","*","","","","","",177.943696`,"","","*","","","","","",27.28`,"","","*"},{179.`,"*","","","","","","Hf","","","*","","","","","",178.9458122`,"","","*","","","","","",13.62`,"","","*"},{180.`,"*","","","","","","Hf","Ta","W","*","","","","","",179.9465457`,179.947462`,179.946701`,"*","","","","","",35.08`,0.012`,0.12`,"*"},{181.`,"*","","","","","","","Ta","","*","","","","","","",180.947992`,"","*","","","","","","",99.988`,"","*"},{182.`,"*","","","","","","","","W","*","","","","","","","",181.948202`,"*","","","","","","","",26.5`,"*"},{183.`,"*","","","","","","","","W","*","","","","","","","",182.95022`,"*","","","","","","","",14.31`,"*"},{184.`,"*","","Os","","","","","","W","*","",183.952488`,"","","","","",183.950928`,"*","",0.02`,"","","","","",30.64`,"*"},{185.`,"*","Re","","","","","","","","*",184.952951`,"","","","","","","","*",37.4`,"","","","","","","","*"},{186.`,"*","","Os","","","","","","W","*","",185.95383`,"","","","","",185.954357`,"*","",1.59`,"","","","","",28.43`,"*"},{187.`,"*","Re","Os","","","","","","","*",186.955744`,186.955741`,"","","","","","","*",62.6`,1.96`,"","","","","","","*"},{188.`,"*","","Os","","","","","","","*","",187.95583`,"","","","","","","*","",13.24`,"","","","","","","*"},{189.`,"*","","Os","","","","","","","*","",188.958137`,"","","","","","","*","",16.15`,"","","","","","","*"},{190.`,"*","","Os","","Pt","","","","","*","",189.958436`,"",189.959917`,"","","","","*","",26.26`,"",0.014`,"","","","","*"},{191.`,"*","","","Ir","","","","","","*","","",190.960584`,"","","","","","*","","",37.3`,"","","","","","*"},{192.`,"*","","Os","","Pt","","","","","*","",191.961467`,"",191.961019`,"","","","","*","",40.78`,"",0.782`,"","","","","*"},{193.`,"*","","","Ir","","","","","","*","","",192.962917`,"","","","","","*","","",62.7`,"","","","","","*"},{194.`,"*","","","","Pt","","","","","*","","","",193.962655`,"","","","","*","","","",32.967`,"","","","","*"},{195.`,"*","","","","Pt","","","","","*","","","",194.964766`,"","","","","*","","","",33.832`,"","","","","*"},{196.`,"*","","","","Pt","","Hg","","","*","","","",195.964926`,"",195.965807`,"","","*","","","",25.242`,"",0.15`,"","","*"},{197.`,"*","","","","","Au","","","","*","","","","",196.966543`,"","","","*","","","","",100.`,"","","","*"},{198.`,"*","","","","Pt","","Hg","","","*","","","",197.967869`,"",197.966743`,"","","*","","","",7.163`,"",9.97`,"","","*"},{199.`,"*","","","","","","Hg","","","*","","","","","",198.968254`,"","","*","","","","","",16.87`,"","","*"},{200.`,"*","","","","","","Hg","","","*","","","","","",199.9683`,"","","*","","","","","",23.1`,"","","*"},{201.`,"*","","","","","","Hg","","","*","","","","","",200.970277`,"","","*","","","","","",13.18`,"","","*"},{202.`,"*","","","","","","Hg","","","*","","","","","",201.970617`,"","","*","","","","","",29.86`,"","","*"},{203.`,"*","","","","","","","Tl","","*","","","","","","",202.97232`,"","*","","","","","","",29.524`,"","*"},{204.`,"*","","","","","","Hg","","Pb","*","","","","","",203.973467`,"",203.97302`,"*","","","","","",6.87`,"",1.4`,"*"},{205.`,"*","","","","","","","Tl","","*","","","","","","",204.974401`,"","*","","","","","","",70.476`,"","*"},{206.`,"*","","","","","","","","Pb","*","","","","","","","",205.97444`,"*","","","","","","","",24.1`,"*"},{207.`,"*","","","","","","","","Pb","*","","","","","","","",206.975872`,"*","","","","","","","",22.1`,"*"},{208.`,"*","","","","","","","","Pb","*","","","","","","","",207.976627`,"*","","","","","","","",52.4`,"*"},{209.`,"*","Bi","","","","","","","","*",208.980374`,"","","","","","","","*",100.`,"","","","","","","","*"},{210.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{211.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{212.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{213.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{214.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{215.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{216.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{217.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{218.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{219.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{220.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{221.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{222.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{223.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{224.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{225.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{226.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{227.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{228.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{229.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{230.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{231.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{232.`,"*","","Th","","","","","","","*","",232.0380508`,"","","","","","","*","",100.`,"","","","","","","*"},{233.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{234.`,"*","","","U","","","","","","*","","",234.0409468`,"","","","","","*","","",0.0055`,"","","","","","*"},{235.`,"*","","","U","","","","","","*","","",235.0439242`,"","","","","","*","","",0.72`,"","","","","","*"},{236.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{237.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{238.`,"*","","","U","","","","","","*","","",238.0507847`,"","","","","","*","","",99.2745`,"","","","","","*"},{239.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"},{240.`,"*","","","","","","","","","*","","","","","","","","","*","","","","","","","","","*"}};
(********************************************************** Load settings *****************************************************************)
loadconstant[constant_]:=Module[{cuparr,extpath},
(* Read Isotope Information: isotope name, mass number, isotopic number, istopic ratio number *)
isoname=constant[[3,3]];
massnb=Select[#,NumberQ]&@constant[[5]];
isonb=Length@massnb;
rationb=isonb-1;
(* Read spike and standard composition *)
spcomp=Take[#,{3,isonb+2}]&@constant[[6]];
stdcomp=Take[#,{3,isonb+2}]&@constant[[7]];
(* Read user preferences: optional parameter, outlier test, MonteCarlo simulation times, normalized isotope, isotopes for findroot, isotope for minimization, number of minimization isotopes *)
optprm="y";
If[constant[[12,3]]=="n"||constant[[12,3]]=="N",optprm="n"];
simnb=500;
If[Positive@constant[[14,3]],simnb=Ceiling@constant[[14,3]]];
optol="y";
If[constant[[13,3]]=="n"||constant[[13,3]]=="N",optol="n"];
norisotope=(Position[#,"x"]&@constant[[9]]-2)[[1,1]];
frisotope=Flatten@Position[#/."X"->"x","x"]&@Drop[#,{norisotope}]&@Drop[#,2]&@constant[[10]];
nmisotope=Flatten@Position[#/."X"->"x","x"]&@Drop[#,{norisotope}]&@Drop[#,2]&@constant[[11]];
nmnb=Length@nmisotope;
(* Read cup configuration: cup configuration, column title for all, isotopes for calculation, isotope for interference correction, isotope for offset correction and their numbers *)
cupconfig=Flatten@Position[#,"Sub"]&@Drop[#,2]&@constant[[17]];
coltitle=Drop[#,2]&@constant[[18]];
subcupnb=Length@cupconfig;
colnb=Length@coltitle;
(* Calculate number of columns in main and sub configuration *)
{maincolnb,subcolnb}={First[#],Differences[#]}&@Append[#, colnb] &@(cupconfig - 1);
cuparr=Take[#,{3,colnb+2}]&@constant[[19]];
isocup=Flatten@Position[#/."X"->"x","x"]&@cuparr;
intcup=Flatten@Position[#/."I"->"i","i"]&@cuparr;
offcup=Flatten@Position[#/."O"->"o","o"]&@cuparr;
intnb=Length@intcup;
(* Optional: all user-defined parameter *)
If[optprm=="n",
extpath=SystemDialogInput["FileOpen", {#, {"xlsx, xls"->{"*.xlsx", "*.xls"}}}, WindowTitle->"Select file for external parameter"]&@Directory[]; 
extparameter=Import@extpath;
workparameter=Flatten[#,1]&@extparameter,
workparameter=defparameter];
];
(* Outlier test *)
oltest[table_]:=Module[{row0,row,avg,stdev,dmax},
row0={};
row=Range@Length@table;
While[row!=row0,
avg=Mean@table[[row]];
stdev=StandardDeviation@table[[row]];
dmax=(Abs[#-avg]/stdev)&/@table;
row0=row;
row=Complement[row,#]&@#[[All,1]]&@Position[dmax,x_/;x>Abs[InverseCDF[NormalDistribution[0,1],1/(4 Length@row)]]];
];
Return@row;
];
(*************************************************** Initialization of data table ************************************************************)
newblk[page_]:=Module[{table,il},
table=DeleteCases[#,Table["",colnb]]&@Drop[#,11]&@page[[All,Range[2,colnb+1]]];
(* Lines that pass the outlier test *)
If[optol=="y",If[MemberQ[StandardDeviation@table,0.],il=oltest@table[[All,Range@maincolnb]],il=oltest@table],il=Range@table];
Return@<|
"name"->page[[1,2]],
"type"->page[[2,2]],
"cynb"->Length@table,
"irow"->il,
"avg"->Mean@table[[il]],
"2sd"->2*StandardDeviation@table[[il]],
"raw"->table
|>;
];
newsmp[page_]:=Module[{table,nb},
table=DeleteCases[#,Table["",colnb]]&@Drop[#,11]&@page[[All,Range[2,colnb+1]]];
nb=Length@table;
Return@<|
"name"->page[[1,2]],
"type"->page[[2,2]],
"smpms"->page[[3,2]],
"smpmserr"->page[[3,4]],
"spcc"->page[[4,2]],
"spccerr"->page[[4,4]],
"spms"->page[[5,2]],
"spmserr"->page[[5,4]],
"blk"->DeleteCases[#,""]&@Select[#,StringQ]&@Rest@page[[6]],
"std"->DeleteCases[#,""]&@Select[#,StringQ]&@Rest@page[[7]],
"blkcor"->{},
"anomaly"->page[[9,Range[2,rationb+1]]]/.""->0.,
"anomalyerr"->page[[10,Range[2,rationb+1]]]/.""->0.,
"franomaly"->{},
"nmanomaly"->{},
"raw"->table,
"bctable"->{},
"otable"->{},
"isotable"->{},
"rtable"->{},
"frdata"->{},
"nmdata"->{},
"frcc"->{},
"frsol"->{},
"frtable"->{},
"nmcc"->{},
"nmsol"->{},
"nmtable"->{},
"frrow"->Range@nb,
"nmrow"->Range@nb,
"simanomaly"->{},
"simfranomaly"->{},
"simnmanomaly"->{},
"simsmpms"->{},
"simspms"->{},
"simspcc"->{},
"simfrbc"->{},
"simnmbc"->{},
"simfrdata"->{},
"simnmdata"->{},
"simfrcc"->{},
"simfrsol"->{},
"simfrratio"->{},
"simfrtable"->{},
"simnmcc"->{},
"simnmsol"->{},
"simnmratio"->{},
"simnmtable"->{},
"frmean"->{},
"nmmean"->{},
"frci"->{},
"nmci"->{},
"frd"->{},
"fr2se"->{},
"nmd"->{},
"nm2se"->{}
|>;
];
(*************************************************** Load sample and blank data ************************************************************)
load[data_]:=Module[{datanb,typelist,blkseq,smpseq},
typelist=data[[All,2,2]];
datanb=Length@typelist;
blkseq=Flatten@Position[#,"b"]&@typelist;
smpseq=Complement[Range@datanb,#]&@blkseq;
blknb=Length@blkseq;
smpnb=Length@smpseq;
(* Create sample list and blank list *)
blklist=newblk/@data[[blkseq]];
smplist=newsmp/@data[[smpseq]];
];
(********************************************** Set working directory and load input data **************************************************)
readfile[]:=Module[{imppath,input,data},
(* Load file *)
imppath=SystemDialogInput["FileOpen", {#, {"xlsx, xls"->{"*.xlsx", "*.xls"}}}, WindowTitle->"Select the input file"]&@NotebookDirectory[]; 
If[imppath==$Canceled,Quit[]];
(* Set working directory *)
SetDirectory@DirectoryName@imppath;
exppath=SystemDialogInput["FileSave",{Directory[],{"xlsx"->{"*.xlsx"}}},WindowTitle->"Save Simulated Data "];
If[exppath==$Canceled,Quit[]];
If[exppath==imppath,Print["WARNING: You are overwriting the input file! Reduction canceled"];Abort[]];
input=Import@imppath; 
settings=First@input;
loadconstant@settings;
data=Rest@input;
load[data];
(* Setting display *)
Print["Code for Originslab Routine Double-Spike"];
Print[ver];
Print[TableForm@(Rest/@settings[[Range@19]]),"\n\n"];
];

(*********************************************************** Settings preparation **************************************************************)
setprep[]:=Module[{isocupmnb,intcupmnb,intelement,wt,wisomass,wabundance,isoprmtable,wisoname,wintelement},
(* Specify main cup configuration *)
maincoltitle=Take[#,maincolnb]&@coltitle;
(* Cups correponding to off cup *)
corcup=Position[maincoltitle,coltitle[[#]]][[1,1]]&/@offcup;
(* Create titles for isotope calculation and for interference correction *)
isotitle=coltitle[[isocup]];
inttitle=coltitle[[intcup]];
(* Separate isotope mass number and element in isotitle and inttitle *)
isocupmnb=ToExpression@StringDrop[#,-2]&/@isotitle;
intcupmnb=ToExpression@StringDrop[#,-2]&/@inttitle;
intelement=StringTake[#,-2]&/@inttitle;
(* Check the sequence of isotitle *)
isocuporder=Ordering@isocupmnb;
(* Rearrange the sequence of isotitle and isotable from small to large based on mass number *)
isotitle=isotitle[[isocuporder]];
(* Row position: wt stands for work parameter title *)
wt=2;
(* Column position: mass number and natural abundance *)
wisomass=11;
wabundance=20;
(* Call work parameters for interference correction: *)
prmtable=workparameter[[Union[isocupmnb,intcupmnb]+wt]];
isoprmtable=workparameter[[massnb+wt]];
(* Determine column position of isotopes and elements; minus 2 for the first two column titles *)
wisoname=(Position[#,isoname]&@isoprmtable[[1]])[[1,1]]-2;
wintelement=Table[(Position[#,intelement[[i]]]&@workparameter[[intcupmnb[[i]]+wt]])[[1,1]],{i,1,intnb}]-2;
(* Read standard composition and isotope mass from work parameter *)
isomass=isoprmtable[[All,wisomass+wisoname]];
(* Transfer values to ratios: mass, spike and standard *)
massratio=Drop[#,{norisotope}]&@isomass/isomass[[norisotope]];
spratio=Drop[#,{norisotope}]&@spcomp/spcomp[[norisotope]];
stdratio=Drop[#,{norisotope}]&@stdcomp/stdcomp[[norisotope]];
(* Create isotopic ratio title *)
ratiotitle=#<>"/"<>isotitle[[norisotope]]&/@(Drop[#,{norisotope}]&@isotitle);
(* Create ratios for exact solving and minimization *)
frmass=massratio[[frisotope]];
frsp=spratio[[frisotope]];
frstd=stdratio[[frisotope]];
nmmass=massratio[[nmisotope]];
nmsp=spratio[[nmisotope]];
nmstd=stdratio[[nmisotope]];
(* Create titles for exact solving and minimization *)
frtitle=ratiotitle[[frisotope]];
nmtitle=ratiotitle[[nmisotope]];
(* Determine coefficient for interference correction *)
intcoef=Table[(isoprmtable[[Range@isonb,wabundance+wintelement[[i]]]]/.""->0)/workparameter[[intcupmnb[[i]]+wt,wabundance+wintelement[[i]]]],{i,intnb}];
mscoef=Table[(isoprmtable[[Range@isonb,wisomass+wintelement[[i]]]]/.""->0)/workparameter[[intcupmnb[[i]]+wt,wisomass+wintelement[[i]]]],{i,intnb}];
(* Principle: {2,4,0}=({1,2,""}/.""\[Rule]0)/0.5; *)
(* Calculate weighted mass ratio of spike and standard *)
mratio=Total[isomass*spcomp]/Total[isomass*stdcomp];
(* Calculate isotopes not used in exact solving *)
afrisotope=Complement[Range@rationb,frisotope];
(* Calculate isotopes not used in minimization, t stands for transposed *)
anmisotope=Complement[Range@rationb,nmisotope];
(* Parameter Display *)
Print[TableForm[Prepend[{ratiotitle,massratio,spratio,stdratio}[[#]],{"ISOTOPE","MASS","SPIKE","STANDARD"}[[#]]]&/@Range@4],"\n\n"];
Print[TableForm[Prepend[Prepend[intcoef,isotitle][[#]],Prepend[inttitle,"INTERFERENCE"][[#]]]&/@Range[intnb+1]],"\n\n"];
];
(*********************************************************** Blank correction ***************************************************************)
blkcor[]:=Module[{blkcorlist,namelist,avglist,exception,rule},
namelist=#["name"]&/@blklist;
avglist=#["avg"]&/@blklist;
exception=Complement[namelist,Flatten[#[["blk"]]&/@smplist]];
If[exception!={},Print[TableForm[exception,TableDirections->Row],"not found! "];Quit[]];
rule=Table[namelist[[i]]->avglist[[i]],{i,blknb}];
For[i=1,i<=smpnb,i++,
smplist[[i]][["blkcor"]]=Mean[smplist[[i]]["blk"]/.rule]/.Mean[{}]->Table[0.,{colnb}];
smplist[[i]][["bctable"]]=#-smplist[[i]]["blkcor"]&/@smplist[[i]]["raw"];
];
];
(********************************************************** Offset correction *************************************************************)
offset[smpinfo_]:=Module[{otable,offcoef,subcuparr},
(* Create table for data after offset correction *)
otable=smpinfo;
If[subcupnb==0,,
For[i=1,i<=subcupnb,i++,
(* Calculated offset coefficient *)
offcoef=smpinfo[[All,corcup[[i]]]]/smpinfo[[All,offcup[[i]]]];
(* Calculate columns to be corrected *)
subcuparr=Table[k,{k,cupconfig[[i]],cupconfig[[i]]+subcolnb[[i]]-1}];
(* Do offset correction *)
otable[[All,subcuparr]]=smpinfo[[All,subcuparr]]*offcoef;
];
];
Return@otable
];
(**************************************************** Interference correction **************************************************************)
intcor[smpinfo_,beta_]:=Module[{isotable,inttable,b},
(* Create tables for isotope calculation and for interference correction *)
inttable=smpinfo[[All,intcup]];
isotable=Part[#,All,isocuporder]&@smpinfo[[All,isocup]];
If[Length@beta==0,b=Table[beta,{Length@inttable}],b=beta];
(* Correct interferences *)
b=Abs@b;
For[i=1,i<=intnb,i++,
isotable-=Table[inttable[[All,i]][[k]]*intcoef[[i]]*(mscoef[[i]]^b[[k]]),{k,Length@inttable}]
];
Return@isotable
];
(*************************************************** Convert intensities to ratios ********************************************************)
ctor[smpinfo_]:=Drop[#,{norisotope}]/#[[norisotope]]&/@smpinfo;
(********************************************* FindRoot and NMinimization for reduction ****************************************************) 
(* Use FindRoot function to do the reduction *)
fr[smp_,sp_,std_,ratio_,anom_]:=Quiet[FindRoot[{smp[[1]]==(p sp[[1]]+(1-p) std[[1]] Exp[anom[[1]]/10000] ratio[[1]]^\[Alpha]) ratio[[1]]^\[Beta],smp[[2]]==(p sp[[2]]+(1-p) std[[2]] Exp[anom[[2]]/10000] ratio[[2]]^\[Alpha]) ratio[[2]]^\[Beta],smp[[3]]==(p sp[[3]]+(1-p) std[[3]] Exp[anom[[3]]/10000] ratio[[3]]^\[Alpha]) ratio[[3]]^\[Beta]},{p,0.1`},{\[Alpha],0.01`},{\[Beta],2`},WorkingPrecision->36]];
(* Use NMinimize function to do the reduction *)
nmin[smp_,sp_,std_,err_,ratio_,nminnb_,anom_]:=Quiet[NMinimize[{Sum[(smp[[i]]-(p*sp[[i]]+(1-p)*std[[i]]*Exp[anom[[i]]/10000]*ratio[[i]]^\[Alpha])*ratio[[i]]^\[Beta])^2/(err[[i]]^2),{i,nminnb}]},{p,\[Alpha],\[Beta]},WorkingPrecision->36,AccuracyGoal->Infinity]];
(********************************************************* Loop for FindRoot ***************************************************************)
frloop[rtable_,frdata_,anomaly_,smpmass_,spmass_,spconc_]:=Module[{frsol,tfrtable1,tfrtable2,frtable,frconc},
(* Calculate p, \[Alpha] and \[Beta] in every sample *)
frsol=If[Depth@anomaly==2,{p,\[Alpha],\[Beta]}//.fr[#,frsp,frstd,frmass,anomaly]&/@frdata,Table[{p,\[Alpha],\[Beta]}//.fr[frdata[[i]],frsp,frstd,frmass,anomaly[[i]]],{i,Length@anomaly}]];
(* Calculate isotopes used in exact solving, t stands for transposed *)
tfrtable1=frstd[[#]]*frmass[[#]]^frsol[[All,2]]&/@Range[3];
tfrtable2=(rtable[[All,afrisotope[[#]]]]/massratio[[afrisotope[[#]]]]^frsol[[All,3]]-frsol[[All,1]]*spratio[[afrisotope[[#]]]])/(1-frsol[[All,1]])&/@Range[rationb-3];
(* Create table to store the final results of fr and rearrange the sequence of isotopes *)
frtable=Transpose@Join[tfrtable1,tfrtable2][[Ordering@Join[frisotope,afrisotope]]];
(* Calculated the concentration of samples *)
frconc=spconc*mratio*((frsol[[All,1]]/spcomp[[norisotope]])/((1-frsol[[All,1]])/(1/(Total/@frtable+1))))/(spmass/smpmass);
Return@{frconc,frsol,frtable}
];
(******************************************************* Loop for NMinimization ************************************************************)
nmloop[rtable_,nmdata_,anomaly_,smpmass_,spmass_,spconc_]:=Module[{datasd,nm,nmsol,tnmtable1,tnmtable2,nmtable,nmconc},
(* Calculate min result, p, \[Alpha] and \[Beta] in every sample *)
datasd=StandardDeviation@nmdata;
nm=If[Depth@anomaly==2,nmin[#,nmsp,nmstd,datasd,nmmass,nmnb,anomaly]&/@nmdata,Table[nmin[nmdata[[i]],nmsp,nmstd,datasd,nmmass,nmnb,anomaly[[i]]],{i,Length@anomaly}]];
nmsol=Flatten@{First[#],{p,\[Alpha],\[Beta]}//.Rest[#]}&/@nm;
(* Calculate isotopes used in minimization, t stands for transposed *)
tnmtable1=nmstd[[#]]*nmmass[[#]]^nmsol[[All,3]]&/@Range[nmnb];
tnmtable2=(rtable[[All,anmisotope[[#]]]]/massratio[[anmisotope[[#]]]]^nmsol[[All,4]]-nmsol[[All,2]]*spratio[[anmisotope[[#]]]])/(1-nmsol[[All,2]])&/@Range[rationb-nmnb];
(* Create table to store the final results of nm and rearrange the sequence of isotopes *)
nmtable=Transpose@Join[tnmtable1,tnmtable2][[Ordering@Join[nmisotope,anmisotope]]];
(* Calculated the concentration of samples *)
nmconc=spconc*mratio*((nmsol[[All,2]]/spcomp[[norisotope]])/((1-nmsol[[All,2]])/(1/(Total/@nmtable+1))))/(spmass/smpmass);
Return@{nmconc,nmsol,nmtable};
];
(********************************************** Prepare data for MonteCarlo simulation *****************************************************)
simnd[val_,err_,nb_]:=RandomReal[#,nb]&@NormalDistribution[val,err];
gensim[table_,row_,nb_]:=Module[{oltable,cynb,mean,cov,maincuparr,elsecuparr,simtable},
oltable=table[[row]];
cynb=Length@oltable;
mean=Mean@oltable;
cov=Covariance@oltable;
(* Generate data for MC simulation. If the subconfig data is the same in every cycle, just use the mean of subconfig data *)
If[subcupnb!=0&&cov[[1,cupconfig[[1]]]]<=10^-24,
maincuparr=Range@maincolnb;
elsecuparr=Range[maincolnb+1,colnb];
simtable=Join[#,mean[[elsecuparr]]]&/@RandomReal[#,nb]&@MultinormalDistribution[#1,#2/#3]&[mean[[maincuparr]],cov[[maincuparr,maincuparr]],cynb],
simtable=RandomReal[#,nb]&@MultinormalDistribution[#1,#2/#3]&[mean,cov,cynb];
];
Return@simtable;
];
(* Calculate confidence interval *)
getci[table_]:=(Quantile[table,0.975]-Quantile[table,0.025])/2;
(******************************************************** Double spike bracketing ****************************************************************)
errprop[num_,denom_,numerr_,denomlist_]:=Module[{denomerr},
denomerr=Sqrt[Total[denomlist^2]/Length@denomlist];
Return[(1000+num/denom)*Sqrt[(numerr/num)^2+(denomerr/denom)^2]];
];
bracket[]:=Module[{namelist,frmlist,nmmlist,frclist,nmclist,frmrule,nmmrule,frcrule,nmcrule,frval,frerrlist,nmval,nmerrlist},
namelist=#[["name"]]&/@smplist;
frmlist=#["frmean"]&/@smplist;
frclist=#["frci"]&/@smplist;
nmmlist=#["nmmean"]&/@smplist;
nmclist=#["nmci"]&/@smplist;

frmrule=Table[namelist[[i]]->frmlist[[i]],{i,smpnb}];
frcrule=Table[namelist[[i]]->frclist[[i]],{i,smpnb}];
nmmrule=Table[namelist[[i]]->nmmlist[[i]],{i,smpnb}];
nmcrule=Table[namelist[[i]]->nmclist[[i]],{i,smpnb}];

frval=Mean/@(#[["std"]]&/@smplist/.frmrule/.{}->{stdratio});
nmval=Mean/@(#[["std"]]&/@smplist/.nmmrule/.{}->{stdratio});

For[i=1,i<=smpnb,i++,
smplist[[i,"frd"]]=1000*(smplist[[i,"frmean"]]/frval[[i]]-1);
smplist[[i,"nmd"]]=1000*(smplist[[i,"nmmean"]]/nmval[[i]]-1);
];
frerrlist=Mean/@(#[["std"]]&/@smplist/.frcrule/.{}->{{Table[0,{rationb}]}});
nmerrlist=Mean/@(#[["std"]]&/@smplist/.nmcrule/.{}->{{Table[0,{rationb}]}});
For[i=1,i<=smpnb,i++,
smplist[[i,"fr2se"]]=errprop[smplist[[i,"frmean"]],frval[[i]],smplist[[i,"frci"]],frerrlist[[i]]];
smplist[[i,"nm2se"]]=errprop[smplist[[i,"nmmean"]],nmval[[i]],smplist[[i,"nmci"]],nmerrlist[[i]]];
];
];
(************************************************************ Display result ******************************************************************)
cmb[t1_,t2_]:=Table[Flatten@{t1[[i]],"",t2[[i]]},{i,Length@t1}];
space[k_]:=Table["",{k}];
cuttail[list_]:=Prepend[Drop[#,-2]&/@Rest@list,First@list];
plot[table_]:=Module[{xaxes},
xaxes=Transpose@{Range@smpnb,Rest@table[[All,1]]};
For[i=1,i<=rationb,i++,
Print@ZoomGraphics@ErrorListPlot[#,ImageSize->Large,PlotLabel->table[[1,i+1]],Frame->True,FrameTicks->{{Automatic,Automatic},{xaxes,Automatic}},PlotRange->{{0,smpnb+1},{Min@(#[[All,1]]-2*#[[All,2]]),Max@(#[[All,1]]+2*#[[All,2]])}}]&@Rest@table[[All,{i+1,i+rationb+2}]];
];
];
sumdata[]:=Module[{etitle,mtitle,blksumtitle,blksum,rsumtitle,frsum,nmsum,stdlist,frdelta,nmdelta},
etitle="EXACT SOLVING USING "<>StringJoin@@Riffle[ratiotitle[[frisotope]]," "];
mtitle="MINIMIZATION USING "<>StringJoin@@Riffle[ratiotitle[[nmisotope]]," "];
sumtitle1=Flatten@{"","RAW DATA",space[colnb],"DATA AFTER BLANK CORRECTION",space[colnb],"DATA AFTER BLANK AND OFFSET CORRECTION",space[colnb],"DATA AFTER BLANK, OFFSET AND INTERFERENCE CORRECTION",space[isonb],"ISOTOPIC RATIO",space[rationb],etitle,space[5+rationb],mtitle,space[5+rationb]};
sumtitle2=Flatten@{"Cycle",coltitle,"",coltitle,"",coltitle,"",isotitle,"",ratiotitle,"",{"Ol","Conc. (ppm)","p","\[Alpha]","\[Beta]"},ratiotitle,"",{"Ol","Conc. (ppm)","Min","p","\[Alpha]","\[Beta]"},ratiotitle};
blksumtitle=Flatten@{"Name",coltitle,"",coltitle};
blksum=Prepend[#,blksumtitle]&@Table[Flatten@{blklist[[k,"name"]],blklist[[k,"avg"]],"",blklist[[k,"2sd"]]},{k,blknb}];
rsumtitle=Flatten@{"Name",ratiotitle,"",ratiotitle};
frsum=Prepend[#,rsumtitle]&@Table[Flatten@{smplist[[k,"name"]],smplist[[k,"frmean"]],"",smplist[[k,"frci"]]},{k,smpnb}];
nmsum=Prepend[#,rsumtitle]&@Table[Flatten@{smplist[[k,"name"]],smplist[[k,"nmmean"]],"",smplist[[k,"nmci"]]},{k,smpnb}];
stdlist="Bracketed by "<>(StringJoin@@Riffle[smplist[[#,"std"]]," "]/.""->"Default")&/@Range@smpnb;
frdelta=Prepend[#,rsumtitle]&@Table[Flatten@{smplist[[k,"name"]],smplist[[k,"frd"]],"",smplist[[k,"fr2se"]],"",stdlist[[k]]},{k,smpnb}];
nmdelta=Prepend[#,rsumtitle]&@Table[Flatten@{smplist[[k,"name"]],smplist[[k,"nmd"]],"",smplist[[k,"nm2se"]],"",stdlist[[k]]},{k,smpnb}];
Print@"DATA SUMMARY\n";
Print["Blank\n",TableForm@blksum,"\n\n"];
Print[etitle,"\n"];
Print["Ratio\n",TableForm@frsum,"\n"];
Print["\[Delta]\n",TableForm@frdelta,"\n\n"];
Print[mtitle,"\n"];
Print["Ratio\n",TableForm@nmsum,"\n"];
Print["\[Delta]\n",TableForm@nmdelta,"\n\n"];
Print["\[Delta] reduced by exact solving using "<>StringJoin@@Riffle[ratiotitle[[frisotope]]," "]];
plot@cuttail@frdelta;
Print["\[Delta] reduced by minimization using "<>StringJoin@@Riffle[ratiotitle[[nmisotope]]," "]];
plot@cuttail@nmdelta;
Return@Join[{{"DATA SUMMARY"}},{{""}},{{""}},{{"Blank"}},blksum,{{""}},{{""}},{{etitle}},{{""}},{{"Ratio"}},frsum,{{""}},{{"\[Delta]"}},frdelta,{{""}},{{""}},{{mtitle}},{{""}},{{"Ratio"}},nmsum,{{""}},{{"\[Delta]"}},nmdelta];
];
expblk[blk_]:=Module[{nb,irow,head,body},
nb=blk[["cynb"]];
irow=Table["x",{nb}];
irow[[blk[["irow"]]]]="o";
head={{"Name",blk[["name"]]},{"Smp/Blk (s/b)",blk[["type"]]},{},{},{},{},{},Prepend[coltitle,"Isotope"],Prepend[blk[["avg"]],"Avg"],Prepend[blk[["2sd"]],"2sd"]};
body=Prepend[#,Flatten@{"Cycle",coltitle,"Oltest"}]&@Table[Flatten@{i,blk[["raw"]][[i]],irow[[i]]},{i,nb}];
Return@Join[head,body];
];
expsmp[smp_]:=Module[{nb,frol,nmol,frsum,nmsum,head,body},
nb=Length@smp[["raw"]];
frol=Table["x",{nb}];
nmol=Table["x",{nb}];
frol[[smp[["frrow"]]]]="o";
nmol[[smp[["nmrow"]]]]="o";
frsum=Table[Flatten@{frol[[i]],smp[["frcc",i]],smp[["frsol",i]],smp[["frtable",i]]},{i,nb}];
nmsum=Table[Flatten@{nmol[[i]],smp[["nmcc",i]],smp[["nmsol",i]],smp[["nmtable",i]]},{i,nb}];
head={{"Name",smp[["name"]]},{"Smp/Blk (s/b)",smp[["type"]]},{"Mass (mg)",smp[["smpms"]],"\[PlusMinus]",smp[["smpmserr"]]},
{"Sp. conc. (ppm)",smp[["spcc"]],"\[PlusMinus]",smp[["spccerr"]]},
{"Sp. mass (mg)",smp[["spms"]],"\[PlusMinus]",smp[["spmserr"]]},
Prepend[smp[["blk"]],"Blk name"],
Prepend[smp[["std"]],"Std name"],
Prepend["e"<>#&/@ratiotitle,"Anomaly"],
Prepend[smp[["anomaly"]],"Val"],
Prepend[smp[["anomalyerr"]],"2SE"]
};
body=Prepend[#,sumtitle1]&@Prepend[#,sumtitle2]&@cmb[cmb[cmb[cmb[cmb[cmb[Table[Prepend[smp[["raw"]][[i]],i],{i,nb}],smp[["bctable"]]],smp[["otable"]]],smp[["isotable"]]],smp[["rtable"]]],frsum],nmsum];
Return@Join[head,body];
];
(************************************************************ Main function *******************************************************************)
readfile[];
setprep[];
blkcor[];
For[j=1,j<=smpnb,j++,
smplist[[j,"otable"]]=offset[smplist[[j,"bctable"]]];
smplist[[j,"franomaly"]]=smplist[[j,"anomaly",frisotope]];
smplist[[j,"nmanomaly"]]=smplist[[j,"anomaly",nmisotope]];
];
blist=Table[inib,{smpnb}];
For[l=1,l<=itnb,l++,
For[j=1,j<=smpnb,j++,
smplist[[j,"isotable"]]=intcor[smplist[[j,"otable"]],blist[[j]]];
smplist[[j,"rtable"]]=ctor[smplist[[j,"isotable"]]];
smplist[[j,"frdata"]]=smplist[[j,"rtable",All,frisotope]];
blist[[j]]=#[[2,All,3]]&@frloop[smplist[[j,"rtable"]],smplist[[j,"frdata"]],smplist[[j,"franomaly"]],smplist[[j,"smpms"]],smplist[[j,"spms"]],smplist[[j,"spcc"]]];
];
];
For[j=1,j<=smpnb,j++,
smplist[[j,"isotable"]]=intcor[smplist[[j,"otable"]],blist[[j]]];
smplist[[j,"rtable"]]=ctor[smplist[[j,"isotable"]]];

smplist[[j,"frdata"]]=smplist[[j,"rtable",All,frisotope]];
{smplist[[j,"frcc"]],smplist[[j,"frsol"]],smplist[[j,"frtable"]]}=frloop[smplist[[j,"rtable"]],smplist[[j,"frdata"]],smplist[[j,"franomaly"]],smplist[[j,"smpms"]],smplist[[j,"spms"]],smplist[[j,"spcc"]]];

smplist[[j,"nmdata"]]=smplist[[j,"rtable",All,nmisotope]];
{smplist[[j,"nmcc"]],smplist[[j,"nmsol"]],smplist[[j,"nmtable"]]}=nmloop[smplist[[j,"rtable"]],smplist[[j,"nmdata"]],smplist[[j,"nmanomaly"]],smplist[[j,"smpms"]],smplist[[j,"spms"]],smplist[[j,"spcc"]]];
];
If[optol=="y",
For[j=1,j<=smpnb,j++,
smplist[[j,"frrow"]]=oltest@smplist[[j,"frtable",All,frisotope]];
smplist[[j,"nmrow"]]=oltest@smplist[[j,"nmtable",All,nmisotope]];
]
];
For[j=1,j<=smpnb,j++,
smplist[[j,"frmean"]]=Mean@smplist[[j,"frtable",smplist[[j,"frrow"]]]];
smplist[[j,"nmmean"]]=Mean@smplist[[j,"nmtable",smplist[[j,"nmrow"]]]];
];
For[j=1,j<=smpnb,j++,
smplist[[j,"simsmpms"]]=simnd[smplist[[j,"smpms"]],smplist[[j,"smpmserr"]],simnb];
smplist[[j,"simspms"]]=simnd[smplist[[j,"spms"]],smplist[[j,"spmserr"]],simnb];
smplist[[j,"simspcc"]]=simnd[smplist[[j,"spcc"]],smplist[[j,"spccerr"]],simnb];

smplist[[j,"simanomaly"]]=Table[If[smplist[[j,"anomalyerr"]][[k]]<10^-6,Table[smplist[[j,"anomaly"]][[k]],{simnb}],simnd[smplist[[j,"anomaly"]][[k]],smplist[[j,"anomalyerr"]][[k]],simnb]],{k,rationb}];
smplist[[j,"simfranomaly"]]=Transpose@smplist[[j,"simanomaly",frisotope]];
smplist[[j,"simnmanomaly"]]=Transpose@smplist[[j,"simanomaly",nmisotope]];

smplist[[j,"simfrbc"]]=gensim[smplist[[j,"bctable"]],smplist[[j,"frrow"]],simnb];
smplist[[j,"simnmbc"]]=gensim[smplist[[j,"bctable"]],smplist[[j,"nmrow"]],simnb];
];
simblist=Table[inib,{smpnb}];
For[l=1,l<=itnb,l++,
For[j=1,j<=smpnb,j++,
smplist[[j,"simfrratio"]]=ctor@intcor[offset@smplist[[j,"simfrbc"]],simblist[[j]]];
smplist[[j,"simfrdata"]]=#[[All,frisotope]]&@smplist[[j,"simfrratio"]];
simblist[[j]]=#[[2,All,3]]&@frloop[smplist[[j,"simfrratio"]],smplist[[j,"simfrdata"]],smplist[[j,"simfranomaly"]],smplist[[j,"simsmpms"]],smplist[[j,"simspms"]],smplist[[j,"simspcc"]]];
];
];
For[j=1,j<=smpnb,j++,
smplist[[j,"simfrratio"]]=ctor@intcor[#,simblist[[j]]]&@offset@smplist[[j,"simfrbc"]];
smplist[[j,"simfrdata"]]=#[[All,frisotope]]&@smplist[[j,"simfrratio"]];

smplist[[j,"simnmratio"]]=ctor@intcor[#,simblist[[j]]]&@offset@smplist[[j,"simnmbc"]];
smplist[[j,"simnmdata"]]=#[[All,nmisotope]]&@smplist[[j,"simnmratio"]];
];
For[j=1,j<=smpnb,j++,
{smplist[[j,"simfrcc"]],smplist[[j,"simfrsol"]],smplist[[j,"simfrtable"]]}=frloop[smplist[[j,"simfrratio"]],smplist[[j,"simfrdata"]],smplist[[j,"simfranomaly"]],smplist[[j,"simsmpms"]],smplist[[j,"simspms"]],smplist[[j,"simspcc"]]];
smplist[[j,"frci"]]=getci@smplist[[j,"simfrtable"]];

{smplist[[j,"simnmcc"]],smplist[[j,"simnmsol"]],smplist[[j,"simnmtable"]]}=nmloop[smplist[[j,"simnmratio"]],smplist[[j,"simnmdata"]],smplist[[j,"nmanomaly"]],smplist[[j,"simsmpms"]],smplist[[j,"simspms"]],smplist[[j,"simspcc"]]];
smplist[[j,"nmci"]]=getci@smplist[[j,"simnmtable"]];
];
bracket[];
summary=sumdata[];
Export[exppath,Join[{settings},{summary},expblk/@blklist,expsmp/@smplist]];
Speak["Reduction complete. "];
];

End[]

EndPackage[]