Project 1 - Parsing NASA thermodynamic data

%Program Matlab for Parsing Nasa Data

clear all
close all
clc


f_data = fopen('THERMO.dat','r');
fgetl(f_data);
TemperaturesString = fgetl(f_data);
TemperaturesStringArray = strsplit(TemperaturesString,' ');
lowTemperature = str2double(TemperaturesStringArray(2));
middleTemperature = str2double(TemperaturesStringArray(3));
highTemperature = str2double(TemperaturesStringArray(4));

T = linspace(lowTemperature,highTemperature,1000);
R = 8.314;


countMolecolars = 0;

fgetl(f_data);
fgetl(f_data);
fgetl(f_data);
molecularList = '';
index_coeff_in_row = 1;
index_coefficients = 0;

while(1)
	stringCurrent = fgetl(f_data);
	endOfFile = strcmp(stringCurrent,'END');
	if (endOfFile == 1)
		fprintf('End of parsing\n');
		break;
	end;
	if (stringCurrent(1) == ' ' || stringCurrent(1) == '-')
		fprintf('Parsing numbers: %s\n',stringCurrent);
		E = strfind(stringCurrent,'E');
        alfa = E(end);
        sizeE = find(E == alfa);
		indexE = 1;
		
		for c = 1: sizeE
            endNumber = E(c)+3
			value(index_coeff_in_row) = str2double(stringCurrent(indexE:endNumber));
			indexE = E(c)+4;
            index_coeff_in_row = index_coeff_in_row + 1;
		end
	else
        countMolecolars = countMolecolars + 1;
		%find element
		molecularString = strsplit(stringCurrent,' ');
        current = molecularString{1};
        space = ' ';
        molecularList = strcat(molecularList,current);
		molecularList = strcat(molecularList,{' '});
	end
end

%calculate molecular weights
elementsList = ['O','H','C','N','A','S'];
elementsWeightList = [15.99, 1.00, 12.01, 14, 39.94, 32.06];


molecularListChar = char(molecularList);
molecularList_Array = strsplit(molecularListChar,' ');

for i = 1:countMolecolars
	weightCurrent = 0;
	currentMolecular = molecularList_Array(i);
    speciesName = char(currentMolecular);
	temp = size(speciesName);
    size_currentMolecular = temp(2);
	for j = 1:size_currentMolecular
		single_element = speciesName(j);
		temp = size(elementsList);
        size_elementsList = temp(2);
		for k = 1:size_elementsList
            cmp = strcmp(elementsList(k),single_element);
			if (cmp == 1)
                if (j < size_currentMolecular)
                    currentMolCharNum = char(speciesName(j+1));
                    c = str2num(currentMolCharNum);
                    if (isnan(c) == 0) 
                        weightCurrent = weightCurrent + elementsWeightList(k)*c;
                    else
                        weightCurrent = weightCurrent + elementsWeightList(k);
                    end
                else
                    weightCurrent = weightCurrent + elementsWeightList(k);
                end
            end
		end
		
	end
	fprintf('Weight of species %s is: %d\n',speciesName,weightCurrent);
	
	%getting coefficients
	
	a1 = value(index_coefficients*14 + 8);
	a2 = value(index_coefficients*14 + 9);
	a3 = value(index_coefficients*14 + 10);
	a4 = value(index_coefficients*14 + 11);
	a5 = value(index_coefficients*14 + 12);
	a6 = value(index_coefficients*14 + 13);
	a7 = value(index_coefficients*14 + 14);
	a8 = value(index_coefficients*14 + 1);
	a9 = value(index_coefficients*14 + 2);
	a10 = value(index_coefficients*14 + 3);
	a11 = value(index_coefficients*14 + 4);
	a12 = value(index_coefficients*14 + 5);
	a13 = value(index_coefficients*14 + 6);
	a14 = value(index_coefficients*14 + 7);
	
    index_coefficients = index_coefficients + 1;
    
	for l = 1:1000
        if(T(l) < middleTemperature)
           	Cp(l) = R*(a8 + (a9*T(l)) + (a10*T(l).^2) + (a11*T(l).^3) + (a12*T(l).^4));
            S(l) = R*(a8.*log(T(l)) + (a9*T(l)) + ((a10*T(l).^2)/2) + ((a11*T(l).^3)/3) + ((a12*T(l).^4)/4)+ a14);
            H(l) = R*T(l).*(a8 + ((a9*T(l))/2) + ((a10*T(l).^2)/3) + ((a11*T(l).^3)/4) + ((a12*T(l).^4)/5) + (a13/T(l))); 
        else
            Cp(l) = R*(a1 + (a2*T(l)) + (a3*T(l).^2) + (a4*T(l).^3) + (a5*T(l).^4));
            S(l) = R*(a1.*log(T(l)) + (a2*T(l)) + ((a3*T(l).^2)/2) + ((a4*T(l).^3)/3) + ((a5*T(l).^4)/4) + a7);
            H(l) = R*T(l).*(a1 + ((a2*T(l))/2) + ((a3*T(l).^2)/3) + ((a4*T(l).^3)/4) + ((a5*T(l).^4)/5) + (a6./T(l)));
        end
    end
    cd('C:\Maurizio\MASTER POWER CONVERTERS\1 Matlab\challenges\Parsing_Nasa_Data');
    mkdir(speciesName);
    cd(speciesName);
    
    
	%graphics plot
	
	%plot Cp
	figure(1);
	plot(T,Cp,'linewidth',3,'color','b');
	xlabel('Temperature(K°)');
	ylabel('Specific Heat(KJ)');
	title(sprintf('Specific Heat vs Temperature of species %s',speciesName));
	grid on
	saveas(figure(1),'Specific Heat.jpg');
    
	%plot S
	figure(2);
	plot(T,S,'linewidth',3,'color','b');
	xlabel('Temperature(K°)');
	ylabel('Entropy(KJ/mol-K°)');
	title(sprintf('Entropy vs Temperature of species %s',speciesName));
	grid on
	saveas(figure(2),'Entropy.jpg');
    
	%plot H
	figure(3);
	plot(T,H,'linewidth',3,'color','b');
	xlabel('Temperature(K°)');
	ylabel('Enthalpy(KJ/mol)');
	title(sprintf('Enthalpy vs Temperature of species %s',speciesName));
	grid on
	saveas(figure(3),'Henthalpy.jpg');
end