Understanding lat/lon to pixel conversion and vice-versa in a geotiff image

Question

I'm relatively new to GIS and trying to understand some basic concepts from first principles.

I'm writing a small javascript tool to display a geotiff image on a page, and show the lat/lon of the cursor as I move it over the image. The maps I'm using are FAA georeferenced charts. I ran gdalinfo on the chart I'm interested in (apologies for the size, it's about 70mb), and this is the relevant output:

Size is 16658, 12340
Coordinate System is:
PROJCS["Lambert Conformal Conic",
    GEOGCS["NAD83",
        DATUM["North_American_Datum_1983",
            SPHEROID["GRS 1980",6378137,298.2572221010042,
                AUTHORITY["EPSG","7019"]],
            AUTHORITY["EPSG","6269"]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433],
        AUTHORITY["EPSG","4269"]],
    PROJECTION["Lambert_Conformal_Conic_2SP"],
    PARAMETER["standard_parallel_1",38.66666666666666],
    PARAMETER["standard_parallel_2",33.33333333333334],
    PARAMETER["latitude_of_origin",38.33333333333334],
    PARAMETER["central_meridian",-121.4166666666667],
    PARAMETER["false_easting",0],
    PARAMETER["false_northing",0],
    UNIT["metre",1,
    AUTHORITY["EPSG","9001"]]]
Origin = (-384321.710823537898250,256254.594788954302203)
Pixel Size = (42.336600677923037,-42.337168148839147)
Corner Coordinates:
Upper Left  ( -384321.711,  256254.595) (125d56'47.40"W, 40d33'33.03"N)
Lower Left  ( -384321.711, -266186.060) (125d40'38.73"W, 35d51'25.41"N)
Upper Right (  320921.383,  256254.595) (117d37'59.81"W, 40d35' 0.49"N)
Lower Right (  320921.383, -266186.060) (117d51'29.18"W, 35d52'48.01"N)
Center      (  -31700.164,   -4965.733) (121d46'44.80"W, 38d17'17.00"N)

From this, I take it that the lat/lons on the map are relative to NAD83, and the map uses an LCC projection (all the parameters are there).

Now, I've tried to do some research on how this works, and my understanding so far is that there are essentially two steps to go from pixel x/y coordinate to lat/lon:

1) apply the affine GeoTransform (in this chart it's defined by the values in Origin and Pixel Size, since it's a north up map)

2) apply the inverse LCC Projection (since I'm going from pixel to lat/lon) to those values.

I ran some tests for step 1) and the affine projection in my code, and indeed, for pixel (0, 0) I got (-384321.711, 256254.595), and for pixel (16658, 12340) I got (320921.383, -266186.060), which suggests I'm on the right track, as those values match the Upper Left and Lower Right data from gdalinfo, respectively.

Now here's where I get lost. I'm using d3's geo module to setup an LCC projection and its inverse, like so:

const projection = (d3
    .geoConicConformal()
    .parallels([stdParallel1, stdParallel2])
    .rotate([-centralMeridian, -latitudeOfOrigin])
    .scale(1)
    .translate([0, 0]);

If I try to use this projection's inverse on the values obtained by my affine transformation (e.g. [-384321.711, 256254.595] for the upper left corner), I get nonsensical results that are very different from what I would have expected. There's a gap in my understanding here, as I think I'm missing some conversion factor that would take these seemingly large coordinates to values that the projection would then map correctly to the expected lat/lon. What am I doing wrong/what concept am I missing and/or misunderstanding?

I would also love an explanation of what these intermediate coordinates (the result of applying the affine transformation) represent.

Answer 1

After searching further, I found better documents describing how to compute the Lambert Conformal Conic projection (see page 19 of that document). In particular, these include parameters for the ellipsoid from the georeferenced coordinate system (GRS 1980 in my case). These computations are more involved than what d3 supports, or the definition of the projection on Wikipedia (the 'a spherical datum' part should have clued me in). And sure enough, writing some test code to implement those calculations generates the expected results.

The fact that these 'correct' computations involve the ellipsoid parameters also answers the second part of my question - what do the large numbers from the affine transformation (e.g. [-384321.711, 256254.595]) mean: they're easting/northing meters from the origin, on the ellipsoid.