%
%       M-File to normalize Paroscientific frequency output
%       press = pnorm(count, N, Fcl, temp, sensor)
%       press = normalized pressure in db
%       count = input data in counts
%       N = period count
%       Fcl = clock frequency, Hz
%       temp = temperature,
%               if a scalar, it is applied to all pressure points
%               if a vector, each pressure point is temperature compensated
%                       for each corresponding temperature
%		sensor=sensor serial number
%		NOTE:  9.75 db = approx. 14.7 psia at sea level
%
%       written by Marinna Martini & Fran Hotchkiss
%               U.S. Geological Survey, Woods Hole, MA

function press=pnorm(count, N, Fcl, temp, sensor);

%  set up pertinent constants
if (exist('N') ~= 1) | (exist('Fcl') ~= 1), % no period count info, so make some up
  N=6208;        % number of pressure pulses counted, typical
  Fcl=2457600;    % clock frequency, typical Tattletale
end

G=0.68946498;     % conversion constant, should never change
                  % psi to decibars (approx meters)

% matrix of calibration constants as of 7/8/94 from Fran for older sensors
% serial  AH  AC  BH  BC  T0H  T0C  T1H  T1C  TR  Date(y m d)
ser1=[660   8665.87 8697.1  6106.69 6117.89 25.78366 25.79047 29.06395 29.05716 21.00 94 7 8];
ser1=[568   3620.68 3624.92 2198.85 2191.82 25.9927  25.9944  29.5062  29.5016  21.00 94 7 8; ser1];
% #571 now has the new constants
%ser1=[571   3966.33 3971.77 2427.78 2422.04 25.4261  25.4274  28.5038  28.4992  21.00 94 7 8; ser1];
ser1=[679   2218.81 2221.33 1294.09 1295.86 25.17697 25.17864 27.83350 27.83246 21.00 94 7 8; ser1];
%ser1=[1274  2555.17 2555.36 1589.12 1566.53 24.13494 24.13572 26.30554 26.30428 21.00 94 7 8; ser1];
%ser1=[1703  2059.95 2060.74 1139.25 1140.38 24.91012 24.90976 27.76735 27.76469 24.00 94 7 8; ser1];
ser1=[6763  8124.33 8119.95 3971.34 3985.86 26.34139 26.34207 26.34139 26.34207 21.00 94 7 8; ser1];
ser1=[12980 9756.34 9753.68 5320.66 5301.23 25.12163 25.12196 25.12163 25.12196 21.00 94 7 8; ser1];
%C	SET PRESSURE CONSTANTS FOR PRESSURE SENSORS
%C	IN THE SENSOR ARRAY 4=679,5=688
%C	6=1045,7=1047,8=1151,9=1557,10=1558,11=1274,12=1703
%C       13=6763,14=12980
%        DATA (SARR(4,J),J=1,8) //
%        DATA (SARR(5,J),J=1,8) /2316.34,2318.71,1359.52,1363.01,
%     +  24.50147,24.50364,26.96107,26.96077/
%        DATA (SARR(6,J),J=1,8) /2313.03,2315.54,1368.50,1372.23,
%     +  24.44891,24.45069,26.90864,26.90775/
%        DATA (SARR(7,J),J=1,8) /2342.65,2345.08,1387.61,1394.77,
%     +  24.28864,24.28999,26.69644,26.69566/
%        DATA (SARR(8,J),J=1,8) /4532.163,4538.507,2943.719,2937.814,
%     +  24.44050,24.44361,26.97610,26.97469/

% Calibration constants for Serial number 1151, 6/13/79
% serial A B T0 Date(y m d)
ser2=[571 3967.827 2435.961 25.42178 99 2 26];
%ser2=[1045 2314.508 1370.476 24.45001 03 9 9; ser2];
ser2=[1045 2314.38215 1370.40123 24.45062 09 10 23; ser2];
ser2=[1151 4532.163 2943.719 24.44361 79 6 13; ser2];
%ser2=[1703 2060.139 1140.825 24.91017 92 12 2];
ser2=[1274 2554.656 1588.374 24.13613 92 12 2; ser2];
ser2=[1557 2332.39  1397.186 24.51533 99 2 26; ser2];
ser2=[1558 2152.52871 1278.06539 24.99152 87 7 28; ser2];

% Calibration constants for temperature compensated calibrations
%	using the CDT0 equation
% ser3=[sensor C1 C2 C3 D1 D2 T1 T2 T3 T4 T5 C D T0 Date(y m d)]
%ser3=[47037 175.1086 -5.70706E-3 -1.87896E-5 0.045638 0.0 24.37787 -1.35891E-4 1.19770E-6 1.23228E-9 0.0 174.9805 0.045638 24.37555 1991 10 28];
ser3=[47037 175.1351 -5.70706e-3 -1.87896e-5 0.045638 0.0 24.37616 -1.35891e-4 1.19770e-6 1.23228e-9 0.0 174.9805 0.045638 24.37555 2003 06 27];
%ser3=[52103 981.2742 -1.76254E-3 -9.66989E-5 0.026698 0.0 27.68403 -6.97677E-5 9.30850E-7 1.10183E-9 0.0 981.1946 0.026698 27.68299 1992 10 26; ser3];
ser3=[52103 981.456	-1.76254E-3 -9.66989E-5 0.026698 0.0 27.68158 -6.97677E-5 9.30850E-7 1.10183E-9	0.0 NaN      NaN      NaN      2009 08 19; ser3];
ser3=[60005 590.4284 -4.05709E-4 -6.89342E-5 0.020903 0.0 27.83571 -1.12340E-4 1.02844E-6 1.30711E-9 0.0 NaN      NaN      NaN      1995 01 16; ser3];
ser3=[60006 587.4469 -1.77288E-3 -7.31431E-5 0.021239 0.0 27.92124 -7.23719E-5 9.83991E-7 1.45109E-9 0.0 NaN      NaN      NaN      1995 01 16; ser3];
ser3=[62493 624.4604 -2.27217e-3 -6.80591e-5 0.029671 0.0 27.70230 -6.21501e-5 1.01351e-6 1.32810e-9 0.0 NaN      NaN      NaN      1995 10 02; ser3];
ser3=[62495 589.1866  4.24177E-4 -6.25844E-5 0.026538 0.0 27.78682 -4.32422E-5 9.63356E-7 1.41121E-9 0.0 NaN      NaN      NaN      1995 10 02; ser3];
%ser3=[68022 1024.946 -3.76508E-3 -1.12975E-4 0.040873 0.0 27.74432 -5.64250E-5 9.08633E-7 1.68019E-9 0.0 NaN      NaN      NaN      1997 10 13; ser3];
ser3=[68022 1024.995 -3.76508E-3 -1.12975E-4 0.040873 0.0 27.74370 -5.64250E-5 9.08633E-7 1.68019E-9 0.0 NaN      NaN      NaN      2004 08 10; ser3];
ser3=[68093 2046.515 -4.52377E-3 -2.06479E-4 0.021522 0.0 27.79637 -6.83120E-5 9.86269E-7 1.26877E-9 0.0 NaN      NaN      NaN      2004 08 10; ser3];
ser3=[69367 230.5027 -4.20275e-3 -7.58202e-6 0.036831 0.0 27.65536 -8.77608e-5 7.73802e-7 2.25430e-9 0.0 NaN      NaN      NaN      2003 09 05; ser3];
%ser3=[69506 991.7378 -8.68840E-3 -1.23004E-4 0.036005 0.0 27.77617 -6.44352E-5 9.71216E-7 1.34311E-9 0.0 NaN      NaN      NaN      1998 02 06; ser3];
%ser3=[69506 991.7185 -8.68840E-3 -1.23004E-4 0.036005 0.0 27.77643 -6.44352E-5 9.71216E-7 1.34311E-9 0.0 NaN      NaN      NaN      2003 09 05; ser3];
ser3=[69506 991.7335 -8.68840E-03 -1.23004E-04 0.036005 0.0 27.77623 -6.44352E-05 9.71216E-07 1.34311E-09	0 NaN      NaN      NaN      2009 8 19; ser3];
ser3=[70390 606.2350  1.58605E-3 -5.55501E-5 0.024865 0.0 27.75276 -7.33351E-5 1.06256E-6 9.14717E-10 0.0 NaN     NaN      NaN      1998 02 06; ser3];
ser3=[38196 2185.636 -1.53026e-2 -3.39564e-4 -0.034337 0.0 24.12266 -5.52091e-5 1.11985e-6  0.0 0.0 NaN     NaN      NaN      0000 00 00; ser3];
ser3=[45417 167.6087 -5.70595e-3 -1.57493e-5 0.052763 0.0 24.48838  -1.12994e-4 1.12587e-6  9.29819e-10  0.0 NaN     NaN      NaN    2004 08 09; ser3];
% Calibration constants for non-temperature compensated calibrations
%	using the CDT0 equation
% ser4=[sensor C D T0 Date(y m d)]
ser4=[1703  1030.427  0.031500  24.91021  1992 12 2];
%ser4=[1274  1277.960  0.069064  24.13617  92 12 2];
%ser4=[45417 167.4283 0.052763 24.49106 1996 11 4; ser4];
%ser4=[454171 167.4819 0.052763 24.48651 2004 8 9; ser4];

% determine if sensor contains valid data
if exist('sensor') ~= 1,
	sensor = 52103;
end
if isstr(sensor),
   sensor=str2num(sensor);
elseif iscell(sensor),
   sensor=str2num(char(sensor));
end
% determine which sensor type we're dealing with
if find(ser1(:,1)==sensor), 
	stype=1;
	sptr=find(ser1(:,1)==sensor);
elseif find(ser2(:,1)==sensor), 
	stype=2;
	sptr=find(ser2(:,1)==sensor);
elseif find(ser3(:,1)==sensor), 
	stype=3;
	sptr=find(ser3(:,1)==sensor);
elseif find(ser4(:,1)==sensor),
	stype=4;
	sptr=find(ser4(:,1)==sensor);
else disp(['Unrecognized Sensor: ' int2str(sensor) ' using default constants.']),
	stype=2;
	sptr=1;
end

if exist('temp') ~= 1,
        temp = 10.5;      % if not provided, make it up!
end

freq=(2*N*Fcl)./count; % freq = sensor output frequency, hz

if stype==1,
	AH=ser1(sptr,2);  AC=ser1(sptr,3);
	BH=ser1(sptr,4);  BC=ser1(sptr,5);
	T0H=ser1(sptr,6);	T0C=ser1(sptr,7);  T1H=ser1(sptr,8);  T1C=ser1(sptr,9);
	TR=ser1(sptr,10);
	% adjust calibration coefficients for temperature dependence
	if length(temp) == 1,  % if temperature is scalar
	  AD=((AC-AH)/TR)*temp;
    BD=((BH-BC)/TR)*temp;
    T0D=((T0H-T0C)/TR)*temp;
    A=AC+AD;
    B=BC+BD;
    T0=T0C+T0D;
    FO=(1/T0)*(10^6);
    % Calculate pressure in psi
    C = 1.0-(freq/FO);
    press = (A.*C) - (B.*(C.^2));
	elseif (length(temp) > 1) & (length(temp) == length(count)),
    A=AC+((AC-AH)./TR).*temp;
    B=BC+((BH-BC)./TR).*temp;
    FO=(1./(T0C+((T0H-T0C)./TR).*temp)).*(10^6);
    % Calculate pressure in psi
    press = (A.*(1.0-freq./FO)) - (B.*((1.0-freq./FO).^2));
    clear A B FO
	end
elseif stype==2,
	A=ser2(sptr,2); B=ser2(sptr,3); T0=ser2(sptr,4);
	K = 1.0-(T0.*freq./1000000);
    % Calculate pressure in psi
	press = A.*K - B.*(K.^2);
elseif stype==3,
    C1=ser3(sptr,2); C2=ser3(sptr,3); C3=ser3(sptr,4); 
    D1=ser3(sptr,5); D2=ser3(sptr,6); 
    T1=ser3(sptr,7); T2=ser3(sptr,8); T3=ser3(sptr,9); T4=ser3(sptr,10); T5=ser3(sptr,11); 
    C=ser3(sptr,12); D=ser3(sptr,13); T0=ser3(sptr,14); 
    % adjust calibration coefficients for temperature dependence
    if length(temp) == 1,  % if temperature is scalar
        C=C1+(C2.*temp)+(C3.*(temp.^2));	%	psia
        D=D1+(D2.*temp);	% no units
        T0=T1+(T2.*temp)+(T3.*(temp.^2))+(T4.*(temp.^3))+(T5.*(temp.^4));	% usec
        % Calculate pressure in psi
        %T0=T0./(1000000);	% convert to sec
        K = 1.0-((T0.*freq./1000000).^2);
        %period=(1./freq).*1000000;
        %K = 1.0-((T0.^2)./(period.^2));
        press = (C.*K).*(1.0 - (D.*K));
    elseif (length(temp) > 1) & (length(temp) == length(count)),
        C=C1+C2.*temp+C3.*(temp.^2);	%	psia
        D=D1+D2.*temp;	% no units
        T0=T1+T2.*temp+T3.*(temp.^2)+T4.*(temp.^3)+T5.*(temp.^4);	% usec
        % Calculate pressure in psi
        %T0=T0./(1000000);	% convert to sec
        K = 1.0-((T0.*freq./1000000).^2);
        press = C.*K.*(1-D.*K);
    else
        % Calculate pressure in psi
        %T0=T0./(1000000);	% convert to sec
        K = 1.0-((T0.*freq./1000000).^2);
        press = C.*K.*(1-D.*K);
    end
elseif stype==4,
	C=ser4(sptr,2); D=ser4(sptr,3); T0=ser4(sptr,4); 
    K = 1.0-((T0.*freq./1000000).^2);
    press = C.*K.*(1-D.*K);
end
%end

% convert psi to decibars
press=press.*G;

% replace trash with NaN so it won't be plotted
flags=find(isfinite(press)==0);
if length(flags) == 1,
    press(flags)=NaN;
elseif isempty(flags) == 0,
    press(flags)=ones(size(flags)).*NaN;
end