Using Scherrer Equation For Calculating The Grain Size

January 15, 2024 Post a Comment

I'm trying to calculate the grain size by Scherrer equation but I have stuck in FWHM. import numpy as np #import math k = 0.94 wave_length = 1.5406e-10 data = np.genfromtxt('G3.

Solution 1:

Here is a working example, assuming you have a normal distribution. I run this in a Jupyter console, so in case you don't, you have to skip the "magic line" (%matplotlib notebook) and add plt.show() at the very end.

%matplotlib notebook
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit
import numpy as np

numb = 500# data size
fwhm_in = 3# set FWHM for the artificial data
sigma = fwhm_in/2/np.sqrt(2*np.log(2))    # calculate sigma
xval = np.linspace(-10, 10, numb)         # calculate x and y values using the formula from Wikipedia (see link in question)
yval = (sigma*np.sqrt(2*np.pi))**(-1)*np.exp(-(xval)**2/(2*sigma**2))+np.random.normal(0, 0.03, numb)

deffitFunc(x, x0, sigm):                 # this defines the fit-functionreturn (sigm*np.sqrt(2*np.pi))**(-1)*np.exp(-(x-x0)**2/(2*sigm**2))

guess = (0.5, 2)                          # tell the code with which values it should start the iteration. Close but not equal to the real values
fitParams, fitCovariance = curve_fit(fitFunc, xval, yval, guess) # do the actual fitprint(fitParams)

print('FWHM_calc = {:.3f}'.format(fwhm_in))
fwhm_fit = 2*fitParams[1]*np.sqrt(2*np.log(2))  # calculate the FWHM from the fitted sigma ( = fitParams[1], since fitParams[0] is the offset x0)print('FWHM_fit = {:.3f}'.format(fwhm_fit))

plt.plot(xval,yval, 'r.', label='data')
plt.plot(xval, fitFunc(xval, fitParams[0], fitParams[1]), 'k-', label='fit', linewidth = 3)

plt.grid(True)
plt.legend()
ax = plt.gca()
ax.axvline(fwhm_fit/2, color='b')
ax.axvline(-fwhm_fit/2, color='b')

Solution 2:

I don't know how to solve this problem in python (at least for this moment). So I made it in Matlab.

    clear all
    clcA= dlmread('YOUR DATAS'); %Firstly add to path
    plot(A(:,1),A(:,2)) %Plotting the graph
    hold onmin_peak= input('Just write a value that is higher than minimum peak values: ');
%This value must be between requested peaks and non-requested peaks(you can see this in graph)
    [yval, yval_i] = findpeaks(A(:,2),'MinPeakHeight',min_peak); %Finding peaks
    scatter(A(yval_i,1), yval); %Showing peaksBeta= []; 
    xval = [];
    fork=1:size(yval_i,1) %Finding x values corresponding to y(peak)valuesxval1= A(yval_i(k),1);
        xval = [xval xval1];
    endTheta= xval / 2; 

    fori=1:size(yval,1) %Finding half of max. peak valuesyval_i1= yval_i(i,1);
        while (yval(i,1))/2 < A(yval_i1+1,2)
        yval_i1 = yval_i1+1;
        endyval_i2= yval_i(i,1);
        while (yval(i,1))/2 < A(yval_i2-1,2)
        yval_i2 = yval_i2-1;
        end

        plot(A(yval_i2,1)*ones(size(A(:,2))), A(:,2));
        plot(A(yval_i1,1)*ones(size(A(:,2))), A(:,2));
    %     hold on
    %     scatter(A(yval_i1,1),A(yval_i1,2))
    %     scatter(A(yval_i2,1),A(yval_i2,2))
        B = abs(A(yval_i1,1)-A(yval_i2,1));
        Beta = [Beta B];
    end
    BetaK=0.94;
    Lambda = 1.5406e-10;

    To = [];
    forj=1:size(Beta,2)
    To1 = (K*Lambda)/(Beta(j)*cos(Theta(j)));
    To = [To To1];
    endTo= abs(To)

Solution 3:

I'm sorry for sharing this as an answer but comments does not provide me upload a picture to clarify the problem. So I've illustrated the problem with an image below (BTW, I can't edit or delete my comments which I wrote reluctantly). I got this graph when I run your code