W_sigma not being updated

Share icon

I have a procedure that involves multiple curve fittings. I call the fitting function with the command

FuncFit/H="01"/NTHR=0 cloudExpansion W_coefX averagesX /X=averagesXTimes /W=averagesXSigmas

and the sigma for each coefficient associated with the fitting is saved in W_sigma. However, when I call the function again on another wave

FuncFit/NTHR=0 cloudExpansion W_coefY averagesY /Y=averagesYTimes /W=averagesYSigmas

W_sigma still contains exactly the same value as it did after the first fit (even though different fit values are produced, and the data stored in averagesXTimes, averagesYTimes, averagesYSigmas, and averagesXSigmas are all different). The coefficients saved in W_coefY and W_coefX are different, but how can W_sigma not change? If I kill the wave W_sigma before calling the second FuncFit, it is not recreated.
Is the current data folder by any chance different when you perform the two fits? W_Sigma will be created in the current data folder. If you view W_sigma in a table, then change the current data folder, the table won't change because that's not the same W_sigma.

If that's not the explanation, please post a copy of your experiment file so I can investigate.

John Weeks
WaveMetrics, Inc.
support@wavemetrics.com

January 20, 2014 at 09:32 am - Permalink
I don't believe the data folder is ever being changed. Everything is being handled in the root folder.

Here is my code: 
#pragma rtGlobals=3		// Use modern global access method and strict wave access.
Function openTemperatureData(pathname)	
	string pathname
	variable binariesIndex = 0
	variable stringIndex
	string currentFilename
	
	
	variable uselessNumber, timeInMs, iterationNumber
	// Get a semicolon-separated list of all files in the folder
	String fileList = IndexedFile($pathName, -1, ".dat")
	

	make/N=500 times//A wave storing the times associated with data points in data waves
	make/N=500 iterations//A wave storing the iteration number associated with data points in data waves

	
	make/N=500 timesToUse//A wave storing the iteration number associated with data points in data waves
	make/N=500 FWHMsXToUse//A wave storing the FWHMs associated with each data wave
	make/N=500 FWHMsYToUse
	make/N=500 FWHMsXSigmaToUse
	make/N=500 FWHMsYSigmaToUse	
	Make/D/O/N=7 waveSizeCoefficients
	make averagesX//The average FWHM at each timestep
	make averagesY
	// Sort using combined alpha and numeric sort
	fileList = SortList(fileList, ";", 16)
	variable filledPoint = 0
	variable filledPointToUse = 0
	make W_sigma
	variable cropSizeX
	variable cropSizeY
	Make/N=(200,200)/W cropped
	variable lastTime = 0////The last time that holds data that can be fit
	variable firstTime = 9999 //The first time that holds data that can be fit	
	variable V_FitError//Holds error status of fit routines.  
	make temp
	make temp3
	
	
	for(binariesIndex = 0; binariesIndex < itemsInList(fileList)-1; binariesIndex += 1)//itemsInList(fileList)-1
		currentFilename = StringFromList(binariesIndex,fileList)
		stringIndex=0
		//FindListItem(currentFilename, fileList)
		sscanf currentFilename, "%d_%dms_a_%d", uselessNumber, timeInMs, iterationNumber
		
		GBLoadWave/Q/B=0/A=binary/T={16,16}/W=1/P=$pathName currentFilename
		string wavename = "binary" + num2str(filledPoint)	
		times[filledPoint] = timeInMs
		Redimension/N=(640,480) $wavename
		wave temp=$wavename
	 	cropped = temp[x+350][y+100]
		
				
		iterations[filledPoint]= iterationNumber
		CurveFit/NTHR=0/Q/N Gauss2D kwCWave=waveSizeCoefficients,  cropped
		if(V_FitError > 0.0002)
			print "Error encountered when attempting to fit to " + currentFilename+".  "  
			print "This is usually caused by a picture that does not contain a coherent cloud.  This file will not be used in analysis.  "
			V_FitError = 0
		endif

			
		if((waveSizeCoefficients[3] < 100) && (waveSizeCoefficients[3] > 0.0001))
			FWHMsXToUse[filledPointToUse]=waveSizeCoefficients[3]*6.7*3*10^-6
			timesToUse[filledPointToUse]=times[filledPoint]
			if (times[filledPoint] > lastTime)
				lastTime = times[filledPoint]
			endif
			if (times[filledPoint] < firstTime)
				firstTime = times[filledPoint]
			endif					
			FWHMsYToUse[filledPointToUse]=waveSizeCoefficients[5]*6.7*3*10^-6
			FWHMsXSigmaToUse[filledPointToUse]=W_sigma[3]*6.7*3*10^-6
			FWHMsYSigmaToUse[filledPointToUse]=W_sigma[5]*6.7*3*10^-6
			filledPointToUse = filledPointToUse + 1
		endif
			
				//Determine time
			filledPoint = filledPoint + 1
		//KillWaves temp3
	endfor
	
	
	
	variable timeXIndex = firstTime
	make averagesXTimes
	make averagesXSigmas
	variable averagesXIndex = 0//Index for array of averaged X data
	variable populationXSize
	variable averageXIndex//The index for finding the original data
	
	
	do
		averageXIndex = 0
		populationXSize = 0
			do
				if( (timeXIndex < timesToUse[averageXIndex] + 0.001) && (timeXIndex > timesToUse[averageXIndex] - 0.001) )
					averagesX[averagesXIndex] = averagesX[averagesXIndex] + (FWHMsXToUse[averageXIndex] * (1/(FWHMsXSigmaToUse[averageXIndex])^2))
					populationXSize = populationXSize + 1/((FWHMsXSigmaToUse[averageXIndex])^2)
				endif
				averageXIndex = averageXIndex + 1
			while(averageXIndex < filledPointToUse - 0.001)//Now we have sum of members for a particular time
			if(populationXSize > 0.0001)//To handle float rounding errors (Igor doesn't have int variables)
				averagesX[averagesXIndex] = averagesX[averagesXIndex]/populationXSize
				averagesXTimes[averagesXIndex] = timeXIndex//Will place the time into the space with the same index number as the x and y data corresponding to that time.
				averagesXTimes = averagesXTimes - averagesXTimes[0]
				averagesXSigmas[averagesXIndex] = sqrt(1/populationXSize)
				if(averagesXIndex > 4)
					averagesXSigmas[averagesXIndex] = averagesXSigmas[averagesXIndex]*2
				endif
				averagesXIndex = averagesXIndex + 1
			endif
		timeXIndex = timeXIndex + 1
	while(timeXIndex <= lastTime)//handles out of order times

	
	Make/D/N=2/O W_coefX
	Duplicate/D timesToUse,velocity_fit
	W_coefX = {2,0.000322167}
	FuncFit/H="01"/NTHR=0 cloudExpansion W_coefX  averagesX /X=averagesXTimes /W=averagesXSigmas /I=1
	

	
	
	
	variable velocityX = W_coefX[0] * 1000
	print "velocityX: " + num2str(velocityX)
	variable temperatureX = 2*(85*1.67*10^-27*0.5*velocityX^2)/(1.38*10^-23)//2 * (85*1.67*10^-27*0.5*(W_coefX[0]*1000)^2) / (1.38 * 10^-23) 
	print "Temperature X: " + num2str(temperatureX) + " +/- " + num2str(W_sigma[0]) + " Kelvins"
	

	variable timeYIndex = firstTime
	make averagesYTimes
	make averagesYSigmas
	variable averagesYIndex = 0//Index for array of averaged Y data
	variable populationYSize
	variable averageYIndex//The index for finding the original data
	
	
	do
		averageYIndex = 0
		populationYSize = 0
			do
				if( (timeYIndex < timesToUse[averageYIndex] + 0.001) && (timeYIndex > timesToUse[averageYIndex] - 0.001) )
					averagesY[averagesYIndex] = averagesY[averagesYIndex] + (FWHMsYToUse[averageYIndex] * (1/(FWHMsYSigmaToUse[averageYIndex])^2))
					populationYSize = populationYSize + 1/((FWHMsYSigmaToUse[averageYIndex])^2)
				endif
				averageYIndex = averageYIndex + 1
			while(averageYIndex < filledPointToUse - 0.001)//Now we have sum of members for a particular time
			if(populationYSize > 0.0001)//To handle float rounding errors (Igor doesn't have int variables)
				averagesY[averagesYIndex] = averagesY[averagesYIndex]/populationYSize
				averagesYTimes[averagesYIndex] = timeYIndex//Will place the time into the space with the same index number as the x and y data corresponding to that time.
				averagesYTimes = averagesYTimes - averagesYTimes[0]
				averagesYSigmas[averagesYIndex] = sqrt(1/populationYSize)
				averagesYIndex = averagesYIndex + 1
			endif
		timeYIndex = timeYIndex + 1
	while(timeYIndex <= lastTime)//handles out of order times
	Make/D/N=2/O W_coefY
	W_coefY = {2,0.000248354}
	//killwaves W_sigma
	FuncFit/NTHR=0 cloudExpansion W_coefY averagesY /Y=averagesYTimes /W=averagesYSigmas /I=1
	


	
	variable velocityY = W_coefX[0] * 1000
	variable temperatureY = 2*(85*1.67*10^-27*0.5*velocityY^2)/(1.38*10^-23)
	print "Temperature Y: " + num2str(temperatureY) + " +/- " + num2str(W_sigma[0]) + " Kelvins"

	print "Average: " + num2str((temperatureY + temperatureX) / 2)
End

Make/D/N=2/O W_coef
W_coef[0] = {50,0.000322167}
FuncFit/NTHR=0 cloudExpansion W_coef  averagesY /X=averagesYTimes /W=averagesXSigmas /I=1 /F={0.950000, 4}


--Now let's all stay positive, and do some science.

January 20, 2014 at 07:59 pm - Permalink
I think it might fix your problem if you remove 
Make W_sigma
(Igor will make that wave whether you want it or not) and ADD 
WAVE W_sigma
AFTER the call to FuncFit.

It is important that the WAVE statement be after FuncFit- it must be after the wave has been created by FuncFit. To learn why, read this in the Igor help:
DisplayHelpTopic "Wave References"
The explanation starts with the sentence "A Wave statement has both a compile-time and a runtime effect."

John Weeks
WaveMetrics, Inc.
support@wavemetrics.com

January 21, 2014 at 09:44 am - Permalink
Are you saying to add "WAVE W_sigma" after each call to FuncFit? I did this, and removed the make W_sigma command at the beginning of the function, and the problem still persists.

I read the topic "Wave References", which doesn't explain why I would need to use "WAVE W_sigma" after calls to FuncFit, since I already know the name of the wave that is going to be created. I can just reference it as W_sigma even if I haven't used "WAVE W_sigma", right?

--Now let's all stay positive, and do some science.

January 21, 2014 at 11:01 am - Permalink
DimensionalDifficulties wrote: Are you saying to add "WAVE W_sigma" after each call to FuncFit? I did this, and removed the make W_sigma command at the beginning of the function, and the problem still persists.

Yes, that's what I was recommending.
I read the topic "Wave References", which doesn't explain why I would need to use "WAVE W_sigma" after calls to FuncFit,

It does, actually.
since I already know the name of the wave that is going to be created. I can just reference it as W_sigma even if I haven't used "WAVE W_sigma", right?

No. The WAVE statement has two functions, one at compile time and one at runtime. That's what is described in that help topic.

At compile time, the compiler encounters a name like W_Sigma and it has to know what type of object it is (wave, variable, string) to compile the right kind of code.
At runtime, the WAVE statement actually causes Igor to go look up the wave and store a pointer to the wave in the wave reference. This is why it has to be after FuncFit- the wave doesn't exist (or may be a different wave) before FuncFit runs, so if the WAVE statement is before, then it can't look it up.

Just to confuse the whole issue, there are certain cases where the compiler knows what kind of object is required and will generate the right kind of code even without the WAVE statement. But let's not worry about that detail just yet.

Try to boil down the code to something I can run easily (I don't want to have to load the data set from an external file), make sure the code is copied into the experiment file, and post the file here. I'll try to take a look at it. Guidance is always appreciated in terms of instructions on how to invoke it, initial guesses if it's required, etc.

John Weeks
WaveMetrics, Inc.
support@wavemetrics.com

January 21, 2014 at 12:51 pm - Permalink
Alright, I've found a syntax error in my code that seems to have been preventing the FuncFit command from working properly. Fixing that, and utilizing your advice concerning the make command, I was able to get the function to work properly and update the W_sigma wave as expected. Thank you for your help. Here is my function in a form that works:


#pragma rtGlobals=3		// Use modern global access method and strict wave access.
Function openTemperatureData(pathname)	
	string pathname
	variable binariesIndex = 0
	variable stringIndex
	string currentFilename
	
	
	variable uselessNumber, timeInMs, iterationNumber
	// Get a semicolon-separated list of all files in the folder
	String fileList = IndexedFile($pathName, -1, ".dat")
	

	make/N=500 times//A wave storing the times associated with data points in data waves
	make/N=500 iterations//A wave storing the iteration number associated with data points in data waves

	
	make/N=500 timesToUse//A wave storing the iteration number associated with data points in data waves
	make/N=500 FWHMsXToUse//A wave storing the FWHMs associated with each data wave
	make/N=500 FWHMsYToUse
	make/N=500 FWHMsXSigmaToUse
	make/N=500 FWHMsYSigmaToUse	
	Make/D/O/N=7 waveSizeCoefficients
	make averagesX//The average FWHM at each timestep
	make averagesY
	// Sort using combined alpha and numeric sort
	fileList = SortList(fileList, ";", 16)
	variable filledPoint = 0
	variable filledPointToUse = 0
	
	variable cropSizeX
	variable cropSizeY
	Make/N=(200,200)/W cropped
	variable lastTime = 0////The last time that holds data that can be fit
	variable firstTime = 9999 //The first time that holds data that can be fit	
	variable V_FitError//Holds error status of fit routines.  
	make temp
	make temp3
	
	
	for(binariesIndex = 0; binariesIndex < itemsInList(fileList)-1; binariesIndex += 1)//itemsInList(fileList)-1
		currentFilename = StringFromList(binariesIndex,fileList)
		stringIndex=0
		//FindListItem(currentFilename, fileList)
		sscanf currentFilename, "%d_%dms_a_%d", uselessNumber, timeInMs, iterationNumber
		
		GBLoadWave/Q/B=0/A=binary/T={16,16}/W=1/P=$pathName currentFilename
		string wavename = "binary" + num2str(filledPoint)	
		times[filledPoint] = timeInMs
		Redimension/N=(640,480) $wavename
		wave temp=$wavename
	 	cropped = temp[x+350][y+100]
		
				
		iterations[filledPoint]= iterationNumber
		CurveFit/NTHR=0/Q/N Gauss2D kwCWave=waveSizeCoefficients,  cropped
		WAVE W_sigma
		if(V_FitError > 0.0002)
			print "Error encountered when attempting to fit to " + currentFilename+".  "  
			print "This is usually caused by a picture that does not contain a coherent cloud.  This file will not be used in analysis.  "
			V_FitError = 0
		endif

			
		if((waveSizeCoefficients[3] < 100) && (waveSizeCoefficients[3] > 0.0001))
			FWHMsXToUse[filledPointToUse]=waveSizeCoefficients[3]*6.7*3*10^-6
			timesToUse[filledPointToUse]=times[filledPoint]
			if (times[filledPoint] > lastTime)
				lastTime = times[filledPoint]
			endif
			if (times[filledPoint] < firstTime)
				firstTime = times[filledPoint]
			endif					
			FWHMsYToUse[filledPointToUse]=waveSizeCoefficients[5]*6.7*3*10^-6
			FWHMsXSigmaToUse[filledPointToUse]=W_sigma[3]*6.7*3*10^-6
			FWHMsYSigmaToUse[filledPointToUse]=W_sigma[5]*6.7*3*10^-6
			filledPointToUse = filledPointToUse + 1
		endif
			
				//Determine time
			filledPoint = filledPoint + 1
		//KillWaves temp3
	endfor
	
	
	
	variable timeXIndex = firstTime
	make averagesXTimes
	make averagesXSigmas
	variable averagesXIndex = 0//Index for array of averaged X data
	variable populationXSize
	variable averageXIndex//The index for finding the original data
	
	
	do
		averageXIndex = 0
		populationXSize = 0
			do
				if( (timeXIndex < timesToUse[averageXIndex] + 0.001) && (timeXIndex > timesToUse[averageXIndex] - 0.001) )
					averagesX[averagesXIndex] = averagesX[averagesXIndex] + (FWHMsXToUse[averageXIndex] * (1/(FWHMsXSigmaToUse[averageXIndex])^2))
					populationXSize = populationXSize + 1/((FWHMsXSigmaToUse[averageXIndex])^2)
				endif
				averageXIndex = averageXIndex + 1
			while(averageXIndex < filledPointToUse - 0.001)//Now we have sum of members for a particular time
			if(populationXSize > 0.0001)//To handle float rounding errors (Igor doesn't have int variables)
				averagesX[averagesXIndex] = averagesX[averagesXIndex]/populationXSize
				averagesXTimes[averagesXIndex] = timeXIndex//Will place the time into the space with the same index number as the x and y data corresponding to that time.
				averagesXTimes = averagesXTimes - averagesXTimes[0]
				averagesXSigmas[averagesXIndex] = sqrt(1/populationXSize)
				if(averagesXIndex > 4)
					averagesXSigmas[averagesXIndex] = averagesXSigmas[averagesXIndex]*2
				endif
				averagesXIndex = averagesXIndex + 1
			endif
		timeXIndex = timeXIndex + 1
	while(timeXIndex <= lastTime)//handles out of order times

	
	Make/D/N=2/O W_coefX
	Duplicate/D timesToUse,velocity_fit
	W_coefX = {2,0.000322167}
	

	KillWaves W_sigma
	FuncFit/Q/H="01"/NTHR=0 cloudExpansion W_coefX  averagesX[0,9] /X=averagesXTimes[0,9] /W=averagesXSigmas[0,9] /I=1
	print "X"
	print W_sigma
	WAVE W_sigma
	print W_sigma
	
	variable velocityX = W_coefX[0] * 1000
	print "velocityX: " + num2str(velocityX)
	variable temperatureX = 2*(85*1.67*10^-27*0.5*velocityX^2)/(1.38*10^-23)//2 * (85*1.67*10^-27*0.5*(W_coefX[0]*1000)^2) / (1.38 * 10^-23) 
	print "Temperature X: " + num2str(temperatureX) + " +/- " + num2str(W_sigma[0]) + " Kelvins"
	

	variable timeYIndex = firstTime
	make averagesYTimes
	make averagesYSigmas
	variable averagesYIndex = 0//Index for array of averaged Y data
	variable populationYSize
	variable averageYIndex//The index for finding the original data
	
	
	do
		averageYIndex = 0
		populationYSize = 0
			do
				if( (timeYIndex < timesToUse[averageYIndex] + 0.001) && (timeYIndex > timesToUse[averageYIndex] - 0.001) )
					averagesY[averagesYIndex] = averagesY[averagesYIndex] + (FWHMsYToUse[averageYIndex] * (1/(FWHMsYSigmaToUse[averageYIndex])^2))
					populationYSize = populationYSize + 1/((FWHMsYSigmaToUse[averageYIndex])^2)
				endif
				averageYIndex = averageYIndex + 1
			while(averageYIndex < filledPointToUse - 0.001)//Now we have sum of members for a particular time
			if(populationYSize > 0.0001)//To handle float rounding errors (Igor doesn't have int variables)
				averagesY[averagesYIndex] = averagesY[averagesYIndex]/populationYSize
				averagesYTimes[averagesYIndex] = timeYIndex//Will place the time into the space with the same index number as the x and y data corresponding to that time.
				averagesYTimes = averagesYTimes - averagesYTimes[0]
				averagesYSigmas[averagesYIndex] = sqrt(1/populationYSize)
				averagesYIndex = averagesYIndex + 1
			endif
		timeYIndex = timeYIndex + 1
	while(timeYIndex <= lastTime)//handles out of order times
	Make/D/N=2/O W_coefY
	

	
	W_coefY = {2,0.000248354} 
	KillWaves W_sigma
	FuncFit/H="01"/Q/NTHR=0 cloudExpansion W_coefY averagesY[0,9] /X=averagesYTimes[0,9] /W=averagesYSigmas[0,9] /I=1


	print "Y"
	print W_sigma
	WAVE W_sigma
	print W_sigma
	

	
	variable velocityY = W_coefX[0] * 1000
	variable temperatureY = 2*(85*1.67*10^-27*0.5*velocityY^2)/(1.38*10^-23)
	print "Temperature Y: " + num2str(temperatureY) + " +/- " + num2str(W_sigma[0]) + " Kelvins"

	print "Average: " + num2str((temperatureY + temperatureX) / 2)
End


--Now let's all stay positive, and do some science.

January 21, 2014 at 04:32 pm - Permalink