tsa_saxs/Minimization/FitTools.py

#!/usr/bin/python

import sys
import os.path
import re
import time
import math

sys.path.append("/mntdirect/_users/semeraro/python_tools") 
import FormFactor
#import PLUV_POPC_RecBuf
from PLUV import SDP_base_POPC_RecBuf, SDP_POPC_RecBuf 

###########################################################################################
###########################################################################################
def TimeCount(t):	
	'Time counting'
	minutes=int(t/60)
	seconds=int(t-minutes*60)
	cents=(t-seconds-minutes*60)*100
	return (minutes,seconds,cents)

###########################################################################################
###########################################################################################
#def FreePar (PAR,FIX,L_LIM,H_LIM):
#	'Redefine all the sub-array related to free parameters'
#	FREE = PAR.copy() #list(PAR)
#	shift=0
#
#	for i in range (len(PAR)):
#		if FIX[i]=="f":
#			del FREE[i-shift]
#			del L_LIM[i-shift]
#			del H_LIM[i-shift]
#			shift=shift+1
#
#	return ( FREE , L_LIM , H_LIM )
	
def FreePar (PAR,FIX,L_LIM_org,H_LIM_org):
	'Redefine all the sub-array related to free parameters'
	FREE = []
	L_LIM = []
	H_LIM = []

	for i in range (len(PAR)):
		if FIX[i]!="f":
			FREE.append(PAR[i])
			L_LIM.append(L_LIM_org[i])
			H_LIM.append(H_LIM_org[i])
			
	return ( FREE , L_LIM , H_LIM )		

###########################################################################################
###########################################################################################
def MergePar (PAR,FREE,FIX):
	'Define the merged parameter list to evaluate the function'
	CALC=[]
	TEMP_PAR=list(PAR)
	TEMP_FREE=list(FREE)
	for el in (FIX):
		if el=="f":
			CALC.append(TEMP_PAR.pop(0)) 
		else:
			CALC.append(TEMP_FREE.pop(0)) 
			TEMP_PAR.pop(0)
	return CALC

###########################################################################################
###########################################################################################
def PrintResults(NAME,PAR,FIX,FREE,X2):	
	'Print Results'

	if X2>1e+7 or X2<1e-4:
		print("Reduced X^2 ...... %.1e\n" %X2 )
	else:
		print("Reduced X^2 ...... %.7s\n" %X2 )
	
	shift=0
	TEMP_NAME=list(NAME)
	TEMP_PAR=list(PAR)
	for i in range (len(TEMP_PAR)):
		if FIX[i]=="f":
			del TEMP_NAME[i-shift]
			del TEMP_PAR[i-shift]
			shift=shift+1
	for i in range(len(FREE)):
		print(" %d.\t%4s ...... %.4e" %(i,TEMP_NAME[i],FREE[i]) )
	print ()
	return

###########################################################################################
###########################################################################################
def PrintResultsFile(title,folder,datafile,function,NAME,PAR_RES,FIX,X2):	
	'Print Results on a File'
	localtime = time.asctime( time.localtime(time.time()) )
	

#	res_file = open(resfile,"a")
#	print(resfile)
#	res_file.write(localtime+"\n")
#	res_file.write("Title: %s\n" %title)
#	res_file.write("Data File: %s\n" %datafile)
#	res_file.write("Model: %s\n" %function)
	
	resfile="./"+folder+"/Results.dat"	
	with open(resfile, 'a') as fl:
		fl.writelines(localtime+"\n")
		fl.writelines("Title: %s\n" %title)
		fl.writelines("Data File: %s\n" %datafile)
		fl.writelines("Model: %s\n\n" %function)
		for i in range(len(NAME)) :
			fl.writelines( "%d\t%4s\t%.7f\n" %(i,NAME[i],PAR_RES[i]) ) 
	

###########################################################################################
###########################################################################################
def Convolution(q0,function,PAR,Dq):
	" Convolution "

	N= 10 
	step= 6 * Dq / N 
	Qv=[ (q0-3*Dq + q*step) for q in range(N) ]
	PAR[-1]=0.0

	CI = 0
	CI = CI + FormFactor.Normal(q0,Qv[0],Dq)*function(Qv[0],PAR)/2.0
	for i in range(1,len(Qv)-2):
		CI = CI + FormFactor.Normal(q0,Qv[i],Dq)*function(Qv[i],PAR)
	CI = CI + FormFactor.Normal(q0,Qv[len(Qv)-1],Dq)*function(Qv[len(Qv)-1],PAR)/2.0
	return CI*step

###########################################################################################
###########################################################################################
def Convolution2(q0,function,PAR,Dql,Dqt):
	" Convolution "
	lam=0.0955
	Dq = lambda Dql,Dqt,q: math.sqrt( ( q*Dql )**2 + ( 2*math.pi*Dqt )**2 )/lam
	N= 10 
	step= 6 * Dq(Dql,Dqt,q0) / N 
	Qv=list()
	for i in range(0,N):
		Qv.append( q0-3*Dq(Dql,Dqt,q0) + i*step )
	TEMP_PAR=list(PAR)
	TEMP_PAR.pop(-1)
	TEMP_PAR.append(0.0)
	CI = 0
	CI = CI + FormFactor.Normal(q0,Qv[0],Dq(Dql,Dqt,Qv[0]))*function(Qv[0],TEMP_PAR)/2.0
	for i in range(1,len(Qv)-2):
		CI = CI + FormFactor.Normal(q0,Qv[i],Dq(Dql,Dqt,Qv[i]))*function(Qv[i],TEMP_PAR)
	CI = CI + FormFactor.Normal(q0,Qv[len(Qv)-1],Dq(Dql,Dqt,Qv[len(Qv)-1]))*function(Qv[len(Qv)-1],TEMP_PAR)/2.0
	return CI*step

###########################################################################################
###########################################################################################
def ChooseFunction (function ):
	" Choose Function "
	
	pltOptions = { 'Form': 0, 'strucCtr': 0, 'parDiff': 0, 'porodCtr': 0, 'debyeCtr': 0 }

########### Sphere
	
	if function=="Sphere":
		intensity	= FormFactor.Sphere
		
	elif function=="Sphere_Normal":
		intensity	= FormFactor.Sphere_Normal
		
	elif function=="SDP_POPC_RecBuf":
		#intensity	= PLUV.SDP_POPC_RecBuf	
		intensity	= SDP_POPC_RecBuf

	elif function=="SDP_base_POPC_RecBuf":
		intensity	= SDP_base_POPC_RecBuf
				
	elif function=="SDP_POPC_RecBuf_conv":
		intensity	= PLUV_POPC_RecBuf.SDP_POPC_RecBuf_conv	
		
###########	Sticky Hard Sphere Potential

#	elif function=="SHS_MSA_a0_Schulz":
#		pltOptions['strucCtr']=1
#		pltOptions['parDiff'] = 3
#		if conv == 0:
#			intensity			= StickyHardSphere.SHS_MSA_a0_Schulz
#			pltOptions['Form']	= FormFactor.Sphere_Schulz
#		else:
#			intensity			= StickyHardSphere.SHS_MSA_a0_Schulz_CONV
#			pltOptions['Form']	= FormFactor.Sphere_Schulz_CONV
#			pltOptions['parDiff']+= 1

###########

	else:
		sys.exit("--- This function name does not exist")
	
	return ( intensity , pltOptions )
	

###########################################################################################
###########################################################################################