Tile translations and trouble with PointMatches

[martinschorb]

Hi,

we have a few datasets that require rotating individual tiles. To achieve this while maintaining their absolute position, I am applying the following transformation strategy:

• translate the tile by minus half its size to have the center at the origin.
• apply rotation matrix
• translate the tile by plus half its size to have the center back at the correct coordinates.
• translate tile to absolute position

This seems to screw up how the PointMatch Client finds its features. They seem only to be only placed in part of the image (indeed looks like the upper left quarter (consistent with the moves) and there are no useful matches found at all. For visualization purposes, I ran a MatchTrial without filtering for consensus sets.

with translations:

same tiles without:

It also screws up a stack where no rotation is applied, so basically just a translation back and forth. Why should this even be considered by the PointMatch Client? Is it considering the first transformation? Doesn't it just look at the pixel values of the images?

I have no idea what goes on, all other stack and tile parameters look sane to me.

Thanks for checking!

MT_translations.zip

[martinschorb]

Is it considering the first transformation?

That does not seem to be the case. If I add a unitary transform first, it still behaves the same...

"transforms" : {
      "type" : "list",
      "specList" : [ {
        "type" : "leaf",
        "className" : "mpicbg.trakem2.transform.AffineModel2D",
        "dataString" : "1.0000000000 0.0000000000 0.0000000000 1.0000000000 0.0000000000 0.0000000000"
      }, {
        "type" : "leaf",
        "className" : "mpicbg.trakem2.transform.AffineModel2D",
        "dataString" : "1.0000000000 0.0000000000 0.0000000000 1.0000000000 -2048.0000000000 -1536.0000000000"
      }, {
        "type" : "leaf",
        "className" : "mpicbg.trakem2.transform.AffineModel2D",
        "dataString" : "1.0000000000 0.0000000000 -0.0000000000 1.0000000000 0.0000000000 0.0000000000"
      }, {
        "type" : "leaf",
        "className" : "mpicbg.trakem2.transform.AffineModel2D",
        "dataString" : "1.0000000000 0.0000000000 0.0000000000 1.0000000000 2048.0000000000 1536.0000000000"
      }, {
        "type" : "leaf",
        "className" : "mpicbg.trakem2.transform.AffineModel2D",
        "dataString" : "1.0000000000 0.0000000000 0.0000000000 1.0000000000 4891.2000000000 -38491.3000000000"
      } ]
    }

[trautmane]

Hi @martinschorb,

Why should this even be considered by the PointMatch Client? Is it considering the first transformation? Doesn't it just look at the pixel values of the images?

Render was originally built to handle tiles that had two lens correction transforms in addition to stitching/alignment transforms. We wanted the point matches to be derived from lens corrected versions of the tiles instead of from the raw source data, so our default approach considered the first 2 transforms as lens correction and applied them before deriving matches. Later we had stacks with one or no lens correction transforms, so for those we simply removed the last transform (if there was one) before rendering for point matches. Finally, I introduced transform labels that allow us to explicitly identify which transforms are lens correction. I still needed to support the old hack-y approaches, so it's a bit of a confusing mess - sorry about that.

I'm guessing you don't have any lens correction. If that is the case, I think the easiest option for you is going to be to flatten your other transforms into a single one - but I'm not sure. Can you post here two tile specs as you have them now?
I'll take a look and let you know the best way to handle things once I see what you have.

Thanks,
Eric

[martinschorb]

Hi Eric,

there are two tilespecs of an identical tile in the zip file attached to my initial post. One has the multiple transforms one only a single translation placing it in space.

Let me know if you need more information.

Thanks,
Martin

[trautmane]

I only see translations in the tile specs in your zip file (no rotations). It's also not clear to me whether you simply want to rotate the aligned result or if you want to rotate some (or all) tiles before you derive matches. Can you post two tile specs along with a screen shot of the source images in their original orientation?

[martinschorb]

In our case, some of the tiles are simply rotated by 90 degrees to match the FOV with the cell shapes. There hardly is any change of a tile's rotation within a stack, so it does no matter whether the matches are derived before or after the rotation.
Since in this sample the tiles are rather sparsely placed but we still like to preserve their absolute positioning and orientation (it's a CLEM experiment where we have Fluorescence data on the whole block), I decided to pre-define the rotations during import.

What is your strategy to make Render handle a more general scenario with tiles changing rotation from slice to slice (as in ssTEM)? Are the tiles loaded with no rotation information and the SIFT PM client catches the rotations?

Here are TS of these 4 tiles, 3 of them are rotated:

[
  {
    "tileId": "0000.0000.00012",
    "layout": {
      "sectionId": "12",
      "temca": "3View",
      "camera": "3View",
      "stageX": 70465.2,
      "stageY": 61726.5,
      "rotation": 90,
      "pixelsize": 10
    },
    "z": 12,
    "minX": 70977,
    "minY": 61214,
    "maxX": 74049,
    "maxY": 65310,
    "width": 4096,
    "height": 3072,
    "minIntensity": 0,
    "maxIntensity": 255,
    "mipmapLevels": {
      "0": {
        "imageUrl": "file:///g/schwab/Beckwith_MSB/Microscopy_data/SBEM_data/230317_MSB30_4_CLEM/tiles/g0000/t0000/230317_MSB30_4_CLEM_g0000_t0000_s00012.tif"
      }
    },
    "transforms": {
      "type": "list",
      "specList": [
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "1.0000000000 0.0000000000 0.0000000000 1.0000000000 -2048.0000000000 -1536.0000000000"
        },
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "0.0000000000 1.0000000000 -1.0000000000 0.0000000000 0.0000000000 0.0000000000"
        },
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "1.0000000000 0.0000000000 0.0000000000 1.0000000000 2048.0000000000 1536.0000000000"
        },
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "1.0000000000 0.0000000000 0.0000000000 1.0000000000 70465.2000000000 61726.5000000000"
        }
      ]
    },
    "meshCellSize": 64
  },
  {
    "tileId": "0019.0000.00012",
    "layout": {
      "sectionId": "12",
      "temca": "3View",
      "camera": "3View",
      "stageX": 67778.1,
      "stageY": 39228.4,
      "rotation": 90,
      "pixelsize": 10
    },
    "z": 12,
    "minX": 68290,
    "minY": 38716,
    "maxX": 71362,
    "maxY": 42812,
    "width": 4096,
    "height": 3072,
    "minIntensity": 0,
    "maxIntensity": 255,
    "mipmapLevels": {
      "0": {
        "imageUrl": "file:///g/schwab/Beckwith_MSB/Microscopy_data/SBEM_data/230317_MSB30_4_CLEM/tiles/g0019/t0000/230317_MSB30_4_CLEM_g0019_t0000_s00012.tif"
      }
    },
    "transforms": {
      "type": "list",
      "specList": [
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "1.0000000000 0.0000000000 0.0000000000 1.0000000000 -2048.0000000000 -1536.0000000000"
        },
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "0.0000000000 1.0000000000 -1.0000000000 0.0000000000 0.0000000000 0.0000000000"
        },
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "1.0000000000 0.0000000000 0.0000000000 1.0000000000 2048.0000000000 1536.0000000000"
        },
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "1.0000000000 0.0000000000 0.0000000000 1.0000000000 67778.1000000000 39228.4000000000"
        }
      ]
    },
    "meshCellSize": 64
  },
  {
    "tileId": "0022.0000.00012",
    "layout": {
      "sectionId": "12",
      "temca": "3View",
      "camera": "3View",
      "stageX": 62974.5,
      "stageY": 47543,
      "rotation": 0,
      "pixelsize": 10
    },
    "z": 12,
    "minX": 62974,
    "minY": 47543,
    "maxX": 67070,
    "maxY": 50615,
    "width": 4096,
    "height": 3072,
    "minIntensity": 0,
    "maxIntensity": 255,
    "mipmapLevels": {
      "0": {
        "imageUrl": "file:///g/schwab/Beckwith_MSB/Microscopy_data/SBEM_data/230317_MSB30_4_CLEM/tiles/g0022/t0000/230317_MSB30_4_CLEM_g0022_t0000_s00012.tif"
      }
    },
    "transforms": {
      "type": "list",
      "specList": [
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "1.0000000000 0.0000000000 0.0000000000 1.0000000000 -2048.0000000000 -1536.0000000000"
        },
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "1.0000000000 0.0000000000 -0.0000000000 1.0000000000 0.0000000000 0.0000000000"
        },
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "1.0000000000 0.0000000000 0.0000000000 1.0000000000 2048.0000000000 1536.0000000000"
        },
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "1.0000000000 0.0000000000 0.0000000000 1.0000000000 62974.5000000000 47543.0000000000"
        }
      ]
    },
    "meshCellSize": 64
  },
  {
    "tileId": "0023.0000.00012",
    "layout": {
      "sectionId": "12",
      "temca": "3View",
      "camera": "3View",
      "stageX": 75294.5,
      "stageY": 52909.8,
      "rotation": 90,
      "pixelsize": 10
    },
    "z": 12,
    "minX": 75806,
    "minY": 52397,
    "maxX": 78878,
    "maxY": 56493,
    "width": 4096,
    "height": 3072,
    "minIntensity": 0,
    "maxIntensity": 255,
    "mipmapLevels": {
      "0": {
        "imageUrl": "file:///g/schwab/Beckwith_MSB/Microscopy_data/SBEM_data/230317_MSB30_4_CLEM/tiles/g0023/t0000/230317_MSB30_4_CLEM_g0023_t0000_s00012.tif"
      }
    },
    "transforms": {
      "type": "list",
      "specList": [
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "1.0000000000 0.0000000000 0.0000000000 1.0000000000 -2048.0000000000 -1536.0000000000"
        },
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "0.0000000000 1.0000000000 -1.0000000000 0.0000000000 0.0000000000 0.0000000000"
        },
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "1.0000000000 0.0000000000 0.0000000000 1.0000000000 2048.0000000000 1536.0000000000"
        },
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "1.0000000000 0.0000000000 0.0000000000 1.0000000000 75294.5000000000 52909.8000000000"
        }
      ]
    },
    "meshCellSize": 64
  }
]

We experience the limited area issue during PointMatch generation even when no rotation is applied, just the translations to cope with the change in boundary box position before and after the rotation that in this case add up to zero shift.

[trautmane]

Hi @martinschorb - apologies for my delayed response (been on holiday) ...

There hardly is any change of a tile's rotation within a stack, so it does no matter whether the matches are derived before or after the rotation.

Given this, I think it would be easiest for you to simply rotate entire stack(s) post-alignment. You could use the java TransformSectionClient to do that - but you could also do you your own thing in python. If you think you'll use a java (or spark) client, let me know and I can improve those clients a bit to make rotation easier. Right now they only support a single transformation, so rotation would have to be done in multiple passes.

What is your strategy to make Render handle a more general scenario with tiles changing rotation from slice to slice (as in ssTEM)? Are the tiles loaded with no rotation information and the SIFT PM client catches the rotations?

Yes, our SIFT matching process handles rotations as long as you use an appropriate --matchModelType (e.g. RIGID). We generally generate same z-layer tile pair matches with --matchModelType TRANSLATION and cross z-layer tile pair matches with --matchModelType RIGID. Assuming the solver also allows rotation, tiles will simply be "fixed" by alignment. However, the solver may not orient the aligned stack just the way you like. In those cases, we either create a new translated/rotated stack post-alignment or we handle the translation/rotation outside of render in the exported volume.

We experience the limited area issue during PointMatch generation even when no rotation is applied

This happens because you have multiple unlabelled transforms for each tile. During matching and solving render needs to know which transforms to include or exclude in the process. Unfortunately I did not introduce explicit transform labels until later, so we support a convention that considers the last transform "positional" (to be excluded from matching) and the first one or two transforms as "lens corrections" (to be included in matching). The convoluted logic to determine what is included/excluded is here.

If you prefer to rotate tiles prior to alignment, here are a couple of ways to do that using the first tile from your example above ...

Original attempt that breaks the render match process:

"specList": [
    { "className": "mpicbg.trakem2.transform.AffineModel2D", "dataString": "1 0  0 1 -2048   -1536"   },
    { "className": "mpicbg.trakem2.transform.AffineModel2D", "dataString": "0 1 -1 0     0       0"   },
    { "className": "mpicbg.trakem2.transform.AffineModel2D", "dataString": "1 0  0 1  2048    1536"   },
    { "className": "mpicbg.trakem2.transform.AffineModel2D", "dataString": "1 0  0 1 70465.2 61726.5" }
  ]

Fix A (concatenate 3 rotation transforms into a single one that becomes the default unlabelled lens correction transform):

"specList": [
    { "className": "mpicbg.trakem2.transform.AffineModel2D", "dataString": "0 1 -1 0  3584.0  -512.0" },
    { "className": "mpicbg.trakem2.transform.AffineModel2D", "dataString": "1 0  0 1 70465.2 61726.5" }
  ]

Fix B (explicitly label rotation transforms as "lens"):

"specList": [
    { "className": "mpicbg.trakem2.transform.AffineModel2D", "dataString": "1 0  0 1 -2048   -1536", 
      "metaData": { "labels": ["lens"] }
    },
    { "className": "mpicbg.trakem2.transform.AffineModel2D", "dataString": "0 1 -1 0     0       0", 
      "metaData": { "labels": ["lens"] }
    },
    { "className": "mpicbg.trakem2.transform.AffineModel2D", "dataString": "1 0  0 1  2048    1536",
      "metaData": { "labels": ["lens"] }
    },
    { "className": "mpicbg.trakem2.transform.AffineModel2D", "dataString": "1 0  0 1 70465.2 61726.5" }
  ]

[martinschorb]

Hi Eric,

thanks a lot for the detailed insights into how transforms are treated. This completely explains the behavior we observe. I had a suspicion that some of the (unlabeled) transforms were considered while some weren't. I think I will simply add the "lens" label to the initial rotations and translations and "positional" to the last.
For better comprehensiveness it might be nice to introduce an additional label (maybe called "included", extending TransformSpecMetaData.LENS_CORRECTION_LABEL) for such more general transformations that are supposed to be included in the matching.

[martinschorb]

Hi @trautmane
what branch would you recommend to test that functionality?
Your fix seems to be part of quadratic-fit and ibeam_msem.

[trautmane]

Hi @martinschorb - I'd use the ibeam_msem branch. It is our active development branch at the moment, but should be fine for you in most cases. Everything in quadratic-fit is now in ibeam_msem and we'll likely remove quadratic-fit sometime soon. If you want to play it a little safer, build from a commit on ibeam_msem that is a week or more old. Let me know if you run into any issues and I can try to fix them or come up with a better solution for you.

[martinschorb]

Cool, thanks. I noticed that you changed to Java 11. Does that come with any significant downstream effects? I guess as long as I point every ClientScript to the JDK pulled into./deploy it should all be fine, correct?

[martinschorb]

oh, and I shortly was confused about the server throwing a 404 when connecting on port 8080. Turns out that the automatic redirect to render-ws/view no longer happens.

[trautmane]

I've been running JDK 11 at Janelia since March, so hopefully you won't hit downstream effects - but of course let me know if you do. Are you manually deploying Jetty instances or are you using Docker?

I'll look at the render-ws/view redirect problem separately - can you create a new GitHub issue for that with your deployment details?

[martinschorb]

Hi @trautmane
so far everything seems to work after upgrading. I am using a deployed Jetty on a VM at the moment, but will also fire up some containers for testing and deployment on other hardware.

I have a follow-up issue with the transformations and PM generation. Probably I still don't fully understand everything.

I have generated my stack using the include label.

This is how the tile specs look.

[
  {
    "tileId": "0000.0003.01234",
    "layout": {
      "sectionId": "1234",
      "temca": "3View",
      "camera": "3View",
      "stageX": 21848.7,
      "stageY": -63996.8,
      "rotation": 0,
      "pixelsize": 10
    },
    "z": 1234,
    "minX": 21495,
    "minY": -64495,
    "maxX": 28346,
    "maxY": -58890,
    "width": 6144,
    "height": 4608,
    "minIntensity": 0,
    "maxIntensity": 255,
    "mipmapLevels": {
      "0": {
        "imageUrl": "file://.../tiles/g0000/t0003/20221014_GA_sample21_40nm_Ines_g0000_t0003_s01234.tif"
      }
    },
    "transforms": {
      "type": "list",
      "specList": [
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "1.0000000000 0.0000000000 0.0000000000 1.0000000000 -3072.0000000000 -2304.0000000000"
        },
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "0.9848077530 -0.1736481777 0.1736481777 0.9848077530 0.0000000000 0.0000000000"
        },
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "1.0000000000 0.0000000000 0.0000000000 1.0000000000 3072.0000000000 2304.0000000000"
        },
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "1.0000000000 0.0000000000 0.0000000000 1.0000000000 21848.7000000000 -63996.8000000000"
        }
      ]
    },
    "meshCellSize": 64
  },
  {
    "tileId": "0000.0004.01234",
    "layout": {
      "sectionId": "1234",
      "temca": "3View",
      "camera": "3View",
      "stageX": 27603.9,
      "stageY": -65011.6,
      "rotation": 0,
      "pixelsize": 10
    },
    "z": 1234,
    "minX": 27250,
    "minY": -65510,
    "maxX": 34101,
    "maxY": -59905,
    "width": 6144,
    "height": 4608,
    "minIntensity": 0,
    "maxIntensity": 255,
    "mipmapLevels": {
      "0": {
        "imageUrl": "file://.../tiles/g0000/t0004/20221014_GA_sample21_40nm_Ines_g0000_t0004_s01234.tif"
      }
    },
    "transforms": {
      "type": "list",
      "specList": [
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "1.0000000000 0.0000000000 0.0000000000 1.0000000000 -3072.0000000000 -2304.0000000000"
        },
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "0.9848077530 -0.1736481777 0.1736481777 0.9848077530 0.0000000000 0.0000000000"
        },
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "1.0000000000 0.0000000000 0.0000000000 1.0000000000 3072.0000000000 2304.0000000000"
        },
        {
          "type": "leaf",
          "className": "mpicbg.trakem2.transform.AffineModel2D",
          "dataString": "1.0000000000 0.0000000000 0.0000000000 1.0000000000 27603.9000000000 -65011.6000000000"
        }
      ]
    },
    "meshCellSize": 64
  },
...]

So I explicitly have not used any include label here. However, no matter whether I use them or not for some of the transforms, the transform I like to skip (in this case the rotation) is applied when choosing pointMatches.

Setting "layout": {..., "rotation": 0, ...} also does not help.

Why is it applying the rotation here?
It is perfectly fine for the entire slice, but for the Matches, it should ignore it...

[martinschorb]

Hi @trautmane

our scientists have managed to come up with yet more acquisition schemes, requiring to revive this issue.

I have a grid of tiles that is rotated by a certain angle (like in the last image above). Within that grid the tiles are nicely ordered and oriented.
For the PointMatch generation, I make it ignore (not label the transformation to be included) the rotation and translation ( = not label the transformation to be included) applied to place the tiles in space.
This works nicely and it picks up the correct orientation of the tiles for generating the matches.

However, I now realized that when this rotation is > 90 deg, the pair correspondences that determine the clip position (overlap region) in the tilepair files are wrong.
(LEFT is actually TOP, etc.. depending on the rotation angle)

Obviously, the tilepair algorithm uses the absolute positions of the tiles in space to determine the neighbors. However, it seems that it does not apply the (inverse?) rotation transformation for the relative positioning of tiles.

Would there be a way to use the transformation labels to also affect the spatial orientation of how the tilepairs are detemined?

A workaround I could think of would be to recalculate the global tile coordinates to match their axis orientation. However, I would like to avoid this because it will not allow to handle grids of tiles with different orientations in one stack.

I hope this was clear...

[trautmane]

Hi @martinschorb - can you send me the tile specs and source images for the two adjacent tiles above that are rotated > 90 deg?

[martinschorb]

Hi @trautmane

render-parameters.json

S021_target_g0000_t00.zip

I don't know whether it is related to the non-standard transformation labels that I put for comprehensiveness.
Thanks for checking.

[martinschorb]

I found tilespec parameters imageCol and imageRow in one of the asap importers.
Are these relevant for the tilepair determination? I think I could easily fill them in in most of my cases.

[trautmane]

Hi @martinschorb - I got the chance to look at this today. Thanks for sharing the example data (it helps me :).

I found tilespec parameters imageCol and imageRow in one of the asap importers.
Are these relevant for the tilepair determination? I think I could easily fill them in in most of my cases.

The imageRow and imageCol fields can be relevant for the TilePairClient if you add --useRowColPositions to the command line when you launch the client. I think doing this is likely the easiest way for you to get what you need. I'd recommend always populating imageRow and imageCol in the tileSpecs even though those fields are technically optional - that's what I do.

Your other alternative is to run the TransformSectionClient with a transform that un-rotates the original stack and then use that un-rotated stack to generate the tile pairs.

[martinschorb]

The imageRow and imageCol fields can be relevant for the TilePairClient if you add --useRowColPositions to the command line when you launch the client. I think doing this is likely the easiest way for you to get what you need. I'd recommend always populating imageRow and imageCol in the tileSpecs even though those fields are technically optional - that's what I do.

Perfect, that's what I'll do. Thanks!