In what format should I pass data when calculating PostGIS distance?

Question

satellite navigation

Hello.
Help me understand what the problem is - I spent a lot of time - in the end there is a solution, but still not what is needed.
There are two points in nmea format that came from a GPS device:

// 1.
  [x] 		=> 	4644.15310
  [x_direct] 	=> 	N
  [y] 		=> 	3322.23860
  [y_direct] 	=> 	E
// 2.	
  [x] 		=> 	4644.15830
  [x_direct] 	=> 	N
  [y] 		=> 	3322.29570                    
  [y_direct] 	=> 	E

By simple calculations, I got the representation in a normal format - the calculations are below:
Point1:
X: 46 + ( 44.1531 / 60 ) = 46.735885
Y: 33 + ( 22.2386 / 60 ) = 33.370643333
Point2:
X: 46 + ( 44, 1583 / 60 ) = 46.735971667
Y: 33 + ( 22.2957 / 60 ) = 33.371595 There are
two points in total, which I can display on Google and see the distance between them - 73 m.
Excellent. The desktop implementation of my project uses the following function:

CREATE OR REPLACE FUNCTION calcdistance(x1 double precision, x2 double precision, y1 double precision, y2 double precision)
    RETURNS double precision AS
  $BODY$
  DECLARE
    fPhimean DOUBLE PRECISION;
    fdLambda DOUBLE PRECISION;
    fdPhi DOUBLE PRECISION;
    fAlpha DOUBLE PRECISION;
    fRho DOUBLE PRECISION;
    fNu DOUBLE PRECISION;
    fR DOUBLE PRECISION;
    fz DOUBLE PRECISION;
    fTemp DOUBLE PRECISION;
    Distance DOUBLE PRECISION;
    d2r DOUBLE PRECISION;
    a DOUBLE PRECISION;
    e2 DOUBLE PRECISION;
    
  BEGIN
  d2r = 0.017453;
  a = 6378137;
  e2 = 0.006739496742337;
    
  IF ((x1=x2) and (y1=y2)) then return(0);
  END IF;

  fdLambda = (y1 - y2) * d2r;
  fdPhi = (x1 - x2) * d2r;
  fPhimean = ((x1 + x2) / 2) * d2r;
  fTemp = 1 - e2 * (Pow(Sin(fPhimean), 2));
  fRho = (a * (1 - e2)) / Pow(fTemp, 1.5);

  fNu = a / (Sqrt(1 - e2 * (Sin(fPhimean) * Sin(fPhimean))));
  fz = Sqrt(Pow(Sin(fdPhi / 2), 2) + Cos(x2 * d2r) * Cos(x1 * d2r) *Pow(Sin(fdLambda / 2), 2));
  fz = 2 * Asin(fz);
  fAlpha = Cos(x2 * d2r) * Sin(fdLambda) * 1 / Sin(fz);
  fAlpha =Asin(fAlpha);
  fR = (fRho * fNu) / ((fRho * Pow(Sin(fAlpha), 2)) + (fNu * Pow(Cos(fAlpha), 2)));
  Distance = (fz * fR);
  return(Distance);

  END;

She counts the distance between two points. The result of her work is correct, according to Google - 73 m.
But I also have PostGIS configured and I would like to solve distance issues, etc. with his help.
I tried different options, for example, like this:

SELECT ROUND(
  CAST(
    ST_Distance_Sphere(
      ST_GeomFromText('POINT(46.735885 33.370643333)',4326), 
      ST_GeomFromText('POINT(46.735971667 33.371595)',4326)
    ) As numeric
  ), 2) as distance;

But the result is not 73 meters, but 106 meters.
Tell me, what's the matter, where to dig - I really want to use PostGIS calculations - because there will be many.

Answer 1

In the example above, X and Y are mixed up. Format difference. (((
Example of correct calculation:

SELECT ROUND(
  CAST(
    ST_Length_Spheroid(
      ST_GeomFromText('LINESTRING(33.370643333 46.735885, 33.371595 46.735971667)',4326),
      'SPHEROID["GRS_1980",6378137,298.257222101]' 
    ) As numeric
  ), 2) as distance;