The Origin Forum
File Exchange
Try Origin for Free
The Origin Forum
Home | Profile | Register | Active Topics | Members | Search | FAQ | Send File to Tech support
 All Forums
 Origin Forum for Programming
 Forum for Origin C
 Using a convolution for fitting

Note: You must be registered in order to post a reply.
To register, click here. Registration is FREE!

Screensize:
UserName:
Password:
Anti-Spam Code:
Format Mode:
Format: BoldItalicizedUnderlineStrikethrough Align LeftCenteredAlign Right Horizontal Rule Insert HyperlinkUpload FileInsert Image Insert CodeInsert QuoteInsert List
   
Message:

* HTML is OFF
* Forum Code is ON
Smilies
Smile [:)] Big Smile [:D] Cool [8D] Blush [:I]
Tongue [:P] Evil [):] Wink [;)] Clown [:o)]
Black Eye [B)] Eight Ball [8] Frown [:(] Shy [8)]
Shocked [:0] Angry [:(!] Dead [xx(] Sleepy [|)]
Kisses [:X] Approve [^] Disapprove [V] Question [?]

 
Check here to subscribe to this topic.
   

T O P I C    R E V I E W
nmr2013 Posted - 06/18/2013 : 06:49:10 AM
Origin Ver. and Service Release (Select Help-->About Origin): 8.5.0G SR1
Operating System: Windows XP


Hi everybody,

I'm trying to use a convolution of a rectangular and a gaussian as fitting function for some experimental data. I used the codes given in the online tutorials (link1 and link2) as a starting point, but unfortunately I can't get the code to work properly, even the one given in the tutorials doesn't produce a proper fit when I follow it. I printed the functions and their convolution results to a datasheet and this part seems to work properly. Unfortunately fitting only calculates the first x/y pair and then stops saying that the fit doesn't converge. The online tutorial files produce the same error (only the first x/y pair is calculated) and I can't identify the reason. Any ideas?


Thanks a lot

Sven

________________________

#pragma warning(error : 15618)
#include <origin.h>
#include <ONLSF.h>
#include <fft_utils.h>
NLFitContext *pCtxt = Project.GetNLFitContext();
if (pCtxt )
{
Worksheet wks = Project.ActiveLayer();
vector<double> vNu;
Dataset Nu(wks, 0);
Dataset test1(wks, 2);
Dataset test2(wks, 3);
Dataset test3(wks, 4);
Dataset test4(wks, 5);
Dataset test5(wks, 6);
Dataset test6(wks, 7);
vNu=Nu;
int iSize=vNu.GetSize();
vector<double> vGauss(iSize), vA(iSize), vFalt(2*iSize-1), vFalt2(2*iSize-1), vRet(iSize);
vGauss=A*exp(-1.0*(vNu-mu)^2/(2.0*(M2*10^3)^2));
double lower=mu-d;
double upper=mu+d;
for(int ii=0; ii<iSize; ii++) { vA[ii]=A;}
for(ii=0; ii<iSize; ii++) {
if(vNu[ii]<lower) vA[ii]=0;
if(vNu[ii]>upper) vA[ii]=0;
continue;
}
BOOL bIsNewParamValues = pCtxt->IsNewParamValues();
if ( bIsNewParamValues )
{
vFalt2.SetSubVector(vA);
vFalt.SetSubVector(vGauss);
test1=vNu; test2=vGauss; test3=vA; test4=vFalt; test5=vFalt2;
int iRet = fft_fft_convolution(2*iSize-1, vFalt, vFalt2);
vFalt.GetSubVector(vRet, floor(iSize/2), iSize + floor(iSize/2)-1 );
test6=vRet;
}
NLSFCURRINFO stCurrInfo;
pCtxt->GetFitCurrInfo(&stCurrInfo);
int nCurrentIndex = stCurrInfo.nCurrDataIndex;
y=vRet[nCurrentIndex];
x;



}
// End of editable part
}
4   L A T E S T    R E P L I E S    (Newest First)
nmr2013 Posted - 06/24/2013 : 03:39:27 AM
Thank you very much, the code you posted worked just fine!
Penn Posted - 06/20/2013 : 03:52:29 AM
Hi Sven,

I tried both tutorials in the same version as you are using, and both work.

If just using your code directly, there is a speed problem, so I change the code as below to improve the speed.

NLFitContext *pCtxt = Project.GetNLFitContext();
if (pCtxt )
{
	static vector vRet;
	
	BOOL bIsNewParamValues = pCtxt->IsNewParamValues();
	if ( bIsNewParamValues )
	{
		Worksheet wks = Project.ActiveLayer();
		vector<double> vNu;
		Dataset Nu(wks, 0);
		Dataset test1(wks, 2);
		Dataset test2(wks, 3);
		Dataset test3(wks, 4);
		Dataset test4(wks, 5);
		Dataset test5(wks, 6);
		Dataset test6(wks, 7);
		vNu=Nu;
		int iSize=vNu.GetSize();
		//vector<double> vGauss(iSize), vA(iSize), vFalt(2*iSize-1), vFalt2(2*iSize-1), vRet(iSize);
		vector<double> vGauss(iSize), vA(iSize), vFalt(2*iSize-1), vFalt2(2*iSize-1);
		vRet.SetSize( iSize );
		
		vGauss=A*exp(-1.0*(vNu-mu)^2/(2.0*(M2*10^3)^2));
		double lower=mu-d;
		double upper=mu+d;
		for(int ii=0; ii<iSize; ii++) { vA[ii]=A;}
		for(ii=0; ii<iSize; ii++) {
		if(vNu[ii]<lower) vA[ii]=0;
		if(vNu[ii]>upper) vA[ii]=0;
		continue;
		}
		
		vFalt2.SetSubVector(vA);
		vFalt.SetSubVector(vGauss);
		test1=vNu; test2=vGauss; test3=vA; test4=vFalt; test5=vFalt2;
		int iRet = fft_fft_convolution(2*iSize-1, vFalt, vFalt2);
		vFalt.GetSubVector(vRet, floor(iSize/2), iSize + floor(iSize/2)-1 );
		test6=vRet;
	}
	NLSFCURRINFO stCurrInfo;
	pCtxt->GetFitCurrInfo(&stCurrInfo);
	int nCurrentIndex = stCurrInfo.nCurrDataIndex;
	y=vRet[nCurrentIndex];
	x;
}

There are two main changes (do like the two tutorials) in the code to improve the speed.
1. Move most of the code into the if(bIsNewParamValues) branch.
2. vRet is defined as static.

With the following initial values for parameters, It is able to get the final results in 8.5.0 SR1, so you can have a try.

A = 4
mu = 777
d = 5249
M2 = 10


Also, there are some issues to be payed attention to.
1. To do such convolution as fitting function, the X data is better with the ascending order and even-space.
2. Please check whether you need to consider the dx, as the second tutorial does.

Penn
nmr2013 Posted - 06/19/2013 : 08:00:12 AM
Hi,

thanks for your help. Do you have the same Origin Version?

The desired fitting function is simply the convolution of a gaussian and a rectangular, as given in the code.

The parameters I used are to be interpreted as follows:

- A: Intensity of the functions (good starting value depends on dataset, the maximum intensity value should be ok)

- mu: Position of the gaussian/Rectangular (good starting value is around the position of the maximum of the experimental data)

- d: half width of the rectangular (s. lower and upper in the code, should be comparable to the full width at half height of the peak)
- M2: second moment of the gaussian (-> width parameter, 10 should be a good value)

The ASCII of the dataset can be downloaded here: http://rapidshare.com/files/61950759/spect.txt

Thanks again, your help is very much appreciated!

Sven
Penn Posted - 06/19/2013 : 02:26:55 AM
Hi Sven,

I have tried both tutorials you mentioned, and they work. Now, you are going to merge these two tutorials. Though the fitting function is able to define and compile successfully, without your data, I am not sure whether it works or not. If you can provide your data, definition of your fitting function, and the proper initial values for each parameters, and how you want to fit your data, then we can have a try.

Penn

The Origin Forum © 2020 Originlab Corporation Go To Top Of Page
Snitz Forums 2000