Saving state.

- Good place for gainz in liberty - Cutting off PR - known issue: boats' rotations while driving are incorrect.
Saving state, liberty + RW T&L optimizations.
2025-05-03 14:24:46 -05:00 · 2025-04-29 17:05:15 -05:00 · 2025-04-29 10:03:49 -05:00 · 2025-04-27 15:25:54 -05:00 · 2025-04-27 14:59:40 -05:00 · 2025-04-27 14:40:30 -05:00
44 changed files with 3341 additions and 1198 deletions
--- a/README.md
+++ b/README.md
@@ -1,14 +1,14 @@
 ## Intro

-dca3 is a port of GTA III for the Dreamcast made by The Gang, using [re3](https://github.com/halpz/re3/tree/master/) as a base.
+dca3 is a port of GTA III/VC for the Dreamcast made by The Gang, using [re3](https://github.com/halpz/re3/tree/master/) as a base.

-re3 a fully reversed source code for GTA III.
+re3 a fully reversed source code for GTA III/VC.

 This project was started by [Stefanos Kornilios Mitsis Poiitidis](https://x.com/poiitidis) and uses [KallistiOS](https://kos-docs.dreamcast.wiki/).

 ## Baking the CDI
 ### Prerequisites
-You need GTA 3 installed. This version has been tested and works:  https://store.rockstargames.com/game/buy-grand-theft-auto-the-trilogy.
+You need Grand Theft Auto III or Grand Theft Auto: Vice City installed. This version has been tested and works:  https://store.rockstargames.com/game/buy-grand-theft-auto-the-trilogy.

 Please note that *SOME VERSIONS* of the game may not work. It has been reported that `d4_gta.mp3` is corrupted sometimes.

@@ -18,28 +18,33 @@ You will also need the following tools installed
 - git-scm http://git-scm.com/downloads/win
 - dreamsdk r3 https://github.com/dreamsdk/dreamsdk/releases

-### Preparing the gta3 folder
+### Cloning the dca3-game repo (this is required once)
 - Open dreamsdk shell
- type `mkdir gta3` (and press enter)
+- type `git clone --branch beta https://gitlab.com/skmp/dca3-game.git` (and press enter).
+- It should take a moment and successfully clone the repo
+- close the dreamsdk shell and proceed to the next step.
+
+### Grand Theft Auto III build (liberty)
+#### Preparing the liberty folder
+- Open dreamsdk shell
+- type `mkdir liberty` (and press enter)
 - type `explorer .` (and press enter)
- This will open a folder named gta3. Copy your gta3 files in there.
+- This will open a folder named liberty. Copy your gta3 files in there.
  - If you use the 2cdrom version of the game, make sure to also copy the contents of the play disc to this folder.
 - close the folder and the dreamsdk shell and proceed to the next step.

-### Cloning the dca3-game repo and downloading the prebuilt elf
+#### Downloading and extracting the prebuilt elf
 - Open dreamsdk shell
- type `git clone --branch alpha https://gitlab.com/skmp/dca3-game.git` (and press enter)
- It should take a moment and successfully clone the repo
- type `cd dca3-game/dreamcast` (and press enter)
+- type `cd dca3-game/liberty` (and press enter)
 - type `explorer .` (and press enter).
 - A folder named dreamcast with some files should be open. Keep it on the side.
- Download the Alpha Prebuilt Elf from https://gitlab.com/skmp/dca3-game/-/jobs/8725216645
- Open artifacts.zip and extract dca3.elf to the folder that was kept open before.
+- Download the *liberty* beta prebuilt elf from https://gitlab.com/skmp/dca3-game/-/releases/beta
+- Open artifacts.zip and extract dca-liberty.elf to the folder that was kept open before.
 - Close the folder and dreamsdk shell

-### Repacking and making a prebuilt cdi FOR GD-EMU
+#### Repacking and making a prebuilt cdi FOR GD-EMU
 - Open dreamsdk shell
- type `cd dca3-game/dreamcast` (and press enter)
+- type `cd dca3-game/liberty` (and press enter)
 - type `make cdi-prebuilt` (and press enter)
  - This should take a while (5-15 mins)
  - Due to an issue with dreamsdk, this won't fully complete the first time
@@ -47,28 +52,80 @@ You will also need the following tools installed
  - It will continue where it left off before
  - It should run to completion now and show "*** Repack Completed Successfully ***"
 - type `explorer .` (and press enter)
- The dreamcast folder should open up, and it should contain dca3.cdi for you (~ 900 megs)
+- The dreamcast folder should open up, and it should contain dca-liberty.cdi for you (~ 900 megs)

-### Repacking and making a prebuilt cdi FOR burning CD-ROM
+#### Repacking and making a prebuilt cdi FOR burning CD-ROM
 - Open dreamsdk shell
- type `cd dca3-game/dreamcast` (and press enter)
+- type `cd dca3-game/liberty` (and press enter)
 - type `make FOR_DISC=1 cdi-prebuilt` (and press enter)
  - This should take a while (5-15 mins)
  - Due to an issue with dreamsdk, this won't fully complete the first time
 - type `make FOR_DISC=1 cdi-prebuilt` (and press enter)
  - It will continue where it left off before
-  - It should run to completion now and show "*** Repack Completed Successfully ***"
+  - It should run to completion now and show "*** CDI Baked Successfully ***"
 - type `explorer .` (and press enter)
- The dreamcast folder should open up, and it should contain dca3.cdi for you (~ 700 megs or ~260 megs)
+- The dreamcast folder should open up, and it should contain dca-liberty.cdi for you (~ 700 megs or ~260 megs)
 - If the .cdi is not ~ 700 megs (linux/mkdcdisc) or ~260 megs (windows/cdi4dc), then you did something wrong.
  - You can type `rm -rf repack-data` (and press enter)
  - And then start this step from the beggining

+### Grand Theft Auto Vice City build (miami)
+#### Preparing the miami folder
+- Open dreamsdk shell
+- type `mkdir miami` (and press enter)
+- type `explorer .` (and press enter)
+- This will open a folder named liberty. Copy your gtavc files in there.
+- close the folder and the dreamsdk shell and proceed to the next step.
+
+#### Downloading and extracting the prebuilt elf
+- Open dreamsdk shell
+- type `cd dca3-game/miami` (and press enter)
+- type `explorer .` (and press enter).
+- A folder named dreamcast with some files should be open. Keep it on the side.
+- Download the *miami* beta prebuilt elf from https://gitlab.com/skmp/dca3-game/-/releases/beta
+- Open artifacts.zip and extract dca-miami.elf to the folder that was kept open before.
+- Close the folder and dreamsdk shell
+
+#### Repacking and making a prebuilt cdi FOR GD-EMU
+- Open dreamsdk shell
+- type `cd dca3-game/miami` (and press enter)
+- type `FOR_DISC=2 make cdi-prebuilt` (and press enter)
+  - This should take a while (5-15 mins)
+  - Due to an issue with dreamsdk, this won't fully complete the first time
+- type `FOR_DISC=2 make cdi-prebuilt` (and press enter)
+  - It will continue where it left off before
+  - You will have to close the dreamshell window a few times and restart this procesure a few times
+  - Eventually it should run to completion now and show "*** CDI Baked Successfully ***"
+- type `explorer .` (and press enter)
+- The dreamcast folder should open up, and it should contain dca-miami.cdi for you (~ 1.5 gigs)
+
+#### Repacking and making a prebuilt cdi FOR burning CD-ROM
+- Open dreamsdk shell
+- type `cd dca3-game/miami` (and press enter)
+- type `make FOR_DISC=1 cdi-prebuilt` (and press enter)
+  - This should take a while (5-15 mins)
+  - Due to an issue with dreamsdk, this won't fully complete the first time
+- type `make FOR_DISC=1 cdi-prebuilt` (and press enter)
+  - It will continue where it left off before
+  - It should run to completion now and show "*** Repack Completed Successfully ***"
+- type `explorer .` (and press enter)
+- The dreamcast folder should open up, and it should contain dca-miami.cdi for you (~700 or ~550 megs)
+
+## Running on emulators
+You must have the 'trails' options turned off from the graphics settings, or a white overlay may appear over the 3d render
+
+## Fine tuning settings
+Two experimental modes, 24 bpp (640x480x24) and Anti Aliasing are provided in the Graphics settings.
+
+- When using HDMI or VGA out, it is recommended to turn enable 24 bpp mode. Note that you also have to disable the 'trails' effect.
+- Anti Aliasing will work with or without trails, however trails enabled will have a bigger performance hit in that mode.
+
+Enabling any of those modes may result in some missing geometry under heavy scenes (more likely with AA mode).
+They can also be combined.

 ## How to report issues
 - Take a photo of your tv/monitor and vmu
 - open a ticket via https://gitlab.com/skmp/dca3-game/-/issues/new
- state which elf you have used (eg, https://gitlab.com/skmp/dca3-game/-/jobs/8725216645)
 - write something descriptive of what is/went wrong

 ## License
--- a/liberty/Makefile
+++ b/liberty/Makefile
@@ -104,6 +104,8 @@ OBJS_TEXCONV += \
 OBJS_O3 = \
 	../vendor/librw/src/dc/rwdc.o \
 	../src/liberty/core/World.o \
+	../src/liberty/core/Zones.o \
+	../src/liberty/core/ZoneCull.o \
 	../src/liberty/collision/Collision.o \
 	../src/liberty/math/math.o \
 	../src/liberty/math/Matrix.o \
@@ -111,11 +113,22 @@ OBJS_O3 = \
 	../src/liberty/math/Rect.o \
 	../src/liberty/math/Vector.o \
 	../vendor/librw/src/base.o \
-	../src/liberty/renderer/Shadows.o
-
-OBJS_NO_FAST_MATH = \
+	../src/liberty/renderer/Shadows.o \
+	../src/liberty/renderer/Renderer.o \
+	../src/liberty/animation/FrameUpdate.o \
+	../src/liberty/animation/RpAnimblend.o \
+	../src/liberty/control/PathFind.o \
 	../src/liberty/core/Cam.o \
-	../src/liberty/core/Camera.o
+	../src/liberty/peds/Ped.o \
+	../src/liberty/peds/PedAI.o \
+	../src/liberty/vehicles/Automobile.o
+
+# ICE list with with -O3
+OBJS_O2 = \
+	../src/liberty/animation/AnimBlendNode.o
+
+# ICE list with -ffast-math
+OBJS_NO_FAST_MATH =

 KOS_CPPFLAGS += -fbuiltin -ffast-math -ffp-contract=fast \
 	-mfsrra -mfsca
@@ -160,6 +173,7 @@ clean:
 	-rm -f IP.BIN
 	-rm -f $(PROJECT_NAME).iso
 	-rm -f $(PROJECT_NAME).ds.iso
+	-rm -f $(PROJECT_NAME)-prebuilt.ds.iso
 	-rm -f $(PROJECT_NAME).cdi
 	-rm -f $(DEPS)
 	-rm -rf $(REPACK_DIR)
@@ -168,6 +182,9 @@ clean:
 $(OBJS_O3): %.o: %.cpp
 	kos-c++ $(CXXFLAGS) $(CPPFLAGS) -O3 -c $< -o $@

+$(OBJS_O2): %.o: %.cpp
+	kos-c++ $(CXXFLAGS) $(CPPFLAGS) -O2 -c $< -o $@
+
 $(OBJS_NO_FAST_MATH): %.o: %.cpp
 	kos-c++ $(CXXFLAGS) $(CPPFLAGS) -O3 -c $< -o $@ -fno-fast-math

@@ -230,6 +247,12 @@ $(PROJECT_NAME).ds.iso: IP.BIN 1ST_READ.BIN $(REPACK_DIR)/repacked $(REPACK_GTA_
 	$(KOS_BASE)/utils/scramble/scramble $(TARGET)-prebuilt.bin 1ST_READ_PREBUILT.BIN
 	mkdir -p $(REPACK_GTA_DIR)

+$(PROJECT_NAME)-prebuilt.ds.iso: IP.BIN 1ST_READ_PREBUILT.BIN $(REPACK_DIR)/repacked $(REPACK_GTA_DIR)/GTA3SF8.b $(REPACK_GTA_DIR)/0GDTEX.PVR $(REPACK_GTA_DIR)/settings.ico
+	rm -f $(PROJECT_NAME)-prebuilt.ds.iso
+	rm -f $(REPACK_GTA_DIR)/1ST_READ.BIN
+	cp 1ST_READ_PREBUILT.BIN $(REPACK_GTA_DIR)/1ST_READ.BIN
+	mkisofs -V $(PROJECT_NAME) -G IP.BIN -r -J -l -o $(PROJECT_NAME)-prebuilt.ds.iso $(REPACK_GTA_DIR)
+
 $(PROJECT_NAME)-prebuilt.iso: IP.BIN 1ST_READ_PREBUILT.BIN $(REPACK_DIR)/repacked $(REPACK_GTA_DIR)/GTA3SF8.b $(REPACK_GTA_DIR)/0GDTEX.PVR $(REPACK_GTA_DIR)/settings.ico
 	rm -f $(REPACK_GTA_DIR)/1ST_READ.BIN
 	cp 1ST_READ_PREBUILT.BIN $(REPACK_GTA_DIR)/1ST_READ.BIN
@@ -268,6 +291,8 @@ cdi-no-repack: $(PROJECT_NAME)-no-repack.cdi

 dsiso: $(PROJECT_NAME).ds.iso

+dsiso-prebuilt: $(PROJECT_NAME)-prebuilt.ds.iso
+
 cdi-prebuilt: $(PROJECT_NAME)-prebuilt.cdi

 sim: $(REPACK_DIR)/repacked
--- a/liberty/common.mk
+++ b/liberty/common.mk
@@ -381,6 +381,7 @@ INCLUDE = \
 -I../src/liberty/skel/win \
 \
 -I../vendor/librw \
+-I../vendor/librw/src/dc \
 \
 -I../vendor/miniLZO \
 \
--- a/miami/Makefile
+++ b/miami/Makefile
@@ -164,6 +164,7 @@ clean:
 	-rm -f IP.BIN
 	-rm -f $(PROJECT_NAME).iso
 	-rm -f $(PROJECT_NAME).ds.iso
+	-rm -f $(PROJECT_NAME)-prebuilt.ds.iso
 	-rm -f $(PROJECT_NAME).cdi
 	-rm -f $(DEPS)
 	-rm -rf $(REPACK_DIR)
@@ -234,6 +235,12 @@ $(PROJECT_NAME).ds.iso: IP.BIN 1ST_READ.BIN $(REPACK_DIR)/repacked $(REPACK_GTA_
 	$(KOS_BASE)/utils/scramble/scramble $(TARGET)-prebuilt.bin 1ST_READ_PREBUILT.BIN
 	mkdir -p $(REPACK_GTA_DIR)

+$(PROJECT_NAME)-prebuilt.ds.iso: IP.BIN 1ST_READ_PREBUILT.BIN $(REPACK_DIR)/repacked $(REPACK_GTA_DIR)/GTAVCSF8.b $(REPACK_GTA_DIR)/0GDTEX.PVR $(REPACK_GTA_DIR)/settings.ico
+	rm -f $(PROJECT_NAME)-prebuilt.ds.iso
+	rm -f $(REPACK_GTA_DIR)/1ST_READ.BIN
+	cp 1ST_READ_PREBUILT.BIN $(REPACK_GTA_DIR)/1ST_READ.BIN
+	mkisofs -V $(PROJECT_NAME) -G IP.BIN -r -J -l -o $(PROJECT_NAME)-prebuilt.ds.iso $(REPACK_GTA_DIR)
+
 $(PROJECT_NAME)-prebuilt.iso: IP.BIN 1ST_READ_PREBUILT.BIN $(REPACK_DIR)/repacked $(REPACK_GTA_DIR)/GTAVCSF8.b $(REPACK_GTA_DIR)/0GDTEX.PVR $(REPACK_GTA_DIR)/settings.ico
 	rm -f $(REPACK_GTA_DIR)/1ST_READ.BIN
 	cp 1ST_READ_PREBUILT.BIN $(REPACK_GTA_DIR)/1ST_READ.BIN
@@ -272,6 +279,8 @@ cdi-no-repack: $(PROJECT_NAME)-no-repack.cdi

 dsiso: $(PROJECT_NAME).ds.iso

+dsiso-prebuilt:  $(PROJECT_NAME)-prebuilt.ds.iso
+
 cdi-prebuilt: $(PROJECT_NAME)-prebuilt.cdi

 sim: $(REPACK_DIR)/repacked
--- a/miami/common.mk
+++ b/miami/common.mk
@@ -396,6 +396,7 @@ INCLUDE = \
 -I../src/miami/skel/win \
 \
 -I../vendor/librw \
+-I../vendor/librw/src/dc \
 \
 -I../vendor/miniLZO \
 \
--- a/src/common/vmu/vmu.cpp
+++ b/src/common/vmu/vmu.cpp
@@ -94,23 +94,26 @@ void VmuProfiler::run() {
            pvr_stats_t pvrStats;   pvr_get_stats(&pvrStats);
            uint32_t    sramStats = snd_mem_available();
            size_t      pvrAvail  = pvr_mem_available();
-            float       fps       = std::accumulate(std::begin(fps_), std::end(fps_), 0.0f) 
-                                    / static_cast<float>(fpsSamples);

            float sh4Mem = heapUtilization();
            float pvrMem = (8_MB - pvrAvail ) / 8_MB * 100.0f;
            float armMem = (2_MB - sramStats) / 2_MB * 100.0f;
-            float vtxBuf = vertBuffUse_;
-            {
-                std::shared_lock lk(mtx_);

-                vmu_printf("FPS :%6.2f\n"
-                           "SH4 :%6.2f%%\n"
-                           "PVR :%6.2f%%\n"
-                           "ARM :%6.2f%%\n"
-                           "VTX :%6.2f%%",
-                           fps, sh4Mem, pvrMem, armMem, vtxBuf);
+            float vtxBuf;
+            float fps;
+            { /* Critical section with main thread. */
+                std::shared_lock lk(mtx_);
+                vtxBuf = vertBuffUse_;
+                fps    = std::accumulate(std::begin(fps_), std::end(fps_), 0.0f)
+                       / static_cast<float>(fpsSamples);
            }
+
+            vmu_printf(" FPS:%6.2f\n"
+                       " RAM:%6.2f%%\n"
+                       "VRAM:%6.2f%%\n"
+                       "SRAM:%6.2f%%\n"
+                       " VTX:%6.2f%%",
+                       fps, sh4Mem, pvrMem, armMem, vtxBuf);
        }
 #endif

@@ -119,15 +122,21 @@ void VmuProfiler::run() {
 }

 void VmuProfiler::updateVertexBufferUsage() {
+#ifndef DC_SH4
    std::unique_lock lk(mtx_);
    updated_ = true;

-#ifdef DC_SH4
-    vertBuffUse_ = vertexBufferUtilization();
+#else
+    pvr_stats_t pvrStats;
+    pvr_get_stats(&pvrStats);
+    float vtxUtil  = vertexBufferUtilization();

-    pvr_stats_t pvrStats; 
-    pvr_get_stats(&pvrStats);       
-    fps_[fpsFrame_++] = pvrStats.frame_rate;
+    { /* Critical section with VMU thread. */
+        std::unique_lock lk(mtx_);
+        vertBuffUse_ = vtxUtil;
+        updated_ = true;
+        fps_[fpsFrame_++] = pvrStats.frame_rate;
+    }

    if(fpsFrame_ >= fpsSamples)
        fpsFrame_ = 0;
--- a/src/liberty/animation/AnimBlendNode.cpp
+++ b/src/liberty/animation/AnimBlendNode.cpp
@@ -39,7 +39,7 @@ CAnimBlendNode::Update(CVector &trans, CQuaternion &rot, float weight)
 	float blend = association->GetBlendAmount(weight);
 	if(blend > 0.0f){
 		float kfAdt = player->GetNextTimeDelta();
-		float t = kfAdt == 0.0f ? 0.0f : (kfAdt - remainingTime)/kfAdt;
+		float t = kfAdt == 0.0f ? 0.0f : dc::Div<true, false>(kfAdt - remainingTime, kfAdt);
 		if(player->type & CAnimBlendSequence::KF_TRANS){
 			auto kfdAt = player->GetNextTranslationDelta();
 			auto kfBt = player->GetPrevTranslation();
@@ -153,7 +153,7 @@ CAnimBlendNode::CalcDeltas(void)
 	if(cos > 1.0f)
 		cos = 1.0f;
 	theta = Acos(cos);
-	invSin = theta == 0.0f ?  0.0f : 1.0f/Sin(theta);
+	invSin = theta == 0.0f ?  0.0f : dc::Invert<true, false>(Sin(theta));
 }

 void
@@ -164,7 +164,7 @@ CAnimBlendNode::GetCurrentTranslation(CVector &trans, float weight)
 	float blend = association->GetBlendAmount(weight);
 	if(blend > 0.0f){
 		auto kfAdt = player->GetNextTimeDelta();
-		float t = kfAdt == 0.0f ? 0.0f : (kfAdt - remainingTime)/kfAdt;
+		float t = kfAdt == 0.0f ? 0.0f : dc::Div<true, false>(kfAdt - remainingTime, kfAdt);
 		if(player->type & CAnimBlendSequence::KF_TRANS){
 			auto kfdAt = player->GetNextTranslationDelta();
 			auto kfBt = player->GetPrevTranslation();
--- a/src/liberty/collision/Collision.cpp
+++ b/src/liberty/collision/Collision.cpp
@@ -1619,7 +1619,7 @@ CCollision::ProcessLineOfSight(const CColLine &line,
 		point.point = matrix * point.point;
 		point.normal = Multiply3x3(matrix, point.normal);
 #else
-		mat_load(reinterpret_cast<matrix_t *>(const_cast<CMatrix *>(&matrix)));
+		dc::mat_load2(matrix);
 		mat_trans_single3_nodiv(point.point.x,
 		                        point.point.y,
 		                        point.point.z);
@@ -1798,7 +1798,7 @@ CCollision::ProcessVerticalLine(const CColLine &line,
 		point.point = matrix * point.point;
 		point.normal = Multiply3x3(matrix, point.normal);
 #else
-		mat_load(reinterpret_cast<matrix_t *>(const_cast<CMatrix *>(&matrix)));
+		dc::mat_load2(matrix);
 		mat_trans_single3_nodiv(point.point.x,
 		                        point.point.y,
 		                        point.point.z);
@@ -2173,8 +2173,7 @@ CCollision::ProcessColModels(const CMatrix &matrixA, CColModel &modelA,
 #ifndef DC_SH4
 	matAB *= matrixA;
 #else
-	mat_load(reinterpret_cast<const matrix_t*>(&matAB));
-	mat_apply(reinterpret_cast<const matrix_t*>(&matrixA));
+	dc::mat_load_apply(matAB, matrixA);
 #endif

 	CColSphere bsphereAB;	// bounding sphere of A in B space
@@ -2246,8 +2245,7 @@ CCollision::ProcessColModels(const CMatrix &matrixA, CColModel &modelA,
 #ifndef DC_SH4
 	matBA *= matrixB;
 #else
-	mat_load(reinterpret_cast<const matrix_t*>(&matBA));
-	mat_apply(reinterpret_cast<const matrix_t*>(&matrixB));
+	dc::mat_load_apply(matBA, matrixB);
 #endif
 	for(i = 0; i < modelB.numSpheres; i++){
 		s.radius = modelB.spheres[i].radius;
@@ -2309,7 +2307,7 @@ CCollision::ProcessColModels(const CMatrix &matrixA, CColModel &modelA,
 	}

 #ifdef DC_SH4
-	mat_load(reinterpret_cast<matrix_t *>(const_cast<CMatrix *>(&matrixB)));
+	dc::mat_load2(matrixB);
 #endif
 	for(i = 0; i < numCollisions; i++){
 #ifndef DC_SH4
--- a/src/liberty/core/Cam.cpp
+++ b/src/liberty/core/Cam.cpp
@@ -3267,24 +3267,24 @@ CCam::Process_BehindBoat(const CVector &CameraTarget, float TargetOrientation, f
 	if(m_bCollisionChecksOn || ResetStatics){
 		CVector TestPoint;
 		// Weird calculations here, also casting bool to float...
-		c = Cos(TargetOrientation);
-		s = Sin(TargetOrientation);
+		c = Cos<false>(TargetOrientation);
+		s = Sin<false>(TargetOrientation);
 		TestPoint = TheCamera.CarZoomValueSmooth * CVector(-c, -s, 0.0f) +
 			(TheCamera.CarZoomValueSmooth+7.0f) * CVector(-c, -s, 0.0f) +
 			TargetCoors;
 		TestPoint.z = WaterLevel + TheCamera.CarZoomValueSmooth;
 		float Test1 = CWorld::GetIsLineOfSightClear(TestPoint, TargetCoors, true, false, false, true, false, true, true);

-		c = Cos(TargetOrientation + 0.8f);
-		s = Sin(TargetOrientation + DEGTORAD(40.0f));
+		c = Cos<false>(TargetOrientation + 0.8f);
+		s = Sin<false>(TargetOrientation + DEGTORAD(40.0f));
 		TestPoint = TheCamera.CarZoomValueSmooth * CVector(-c, -s, 0.0f) +
 			(TheCamera.CarZoomValueSmooth+7.0f) * CVector(-c, -s, 0.0f) +
 			TargetCoors;
 		TestPoint.z = WaterLevel + TheCamera.CarZoomValueSmooth;
 		float Test2 = CWorld::GetIsLineOfSightClear(TestPoint, TargetCoors, true, false, false, true, false, true, true);

-		c = Cos(TargetOrientation - 0.8);
-		s = Sin(TargetOrientation - DEGTORAD(40.0f));
+		c = Cos<false>(TargetOrientation - 0.8);
+		s = Sin<false>(TargetOrientation - DEGTORAD(40.0f));
 		TestPoint = TheCamera.CarZoomValueSmooth * CVector(-c, -s, 0.0f) +
 			(TheCamera.CarZoomValueSmooth+7.0f) * CVector(-c, -s, 0.0f) +
 			TargetCoors;
@@ -3307,8 +3307,7 @@ CCam::Process_BehindBoat(const CVector &CameraTarget, float TargetOrientation, f
 			DeltaBeta = TargetOrientation - Beta;
 		}

-		c = Cos(Beta);
-		s = Sin(Beta);
+		auto [s, c] = SinCos<false>(Beta);
 		TestPoint.x = TheCamera.CarZoomValueSmooth * -c +
 			(TheCamera.CarZoomValueSmooth + 7.0f) * -c +
 			TargetCoors.x;
@@ -3333,8 +3332,7 @@ CCam::Process_BehindBoat(const CVector &CameraTarget, float TargetOrientation, f
 		// inlined
 		WellBufferMe(TargetWhenChecksWereOn, &Beta, &BetaSpeed, 0.07f, 0.015f, true);

-		s = Sin(Beta);
-		c = Cos(Beta);
+		auto [s, c] = SinCos<false>(Beta);
 		Source = TheCamera.CarZoomValueSmooth * CVector(-c, -s, 0.0f) +
 			(TheCamera.CarZoomValueSmooth+7.0f) * CVector(-c, -s, 0.0f) +
 			TargetCoors;
--- a/src/liberty/core/Camera.cpp
+++ b/src/liberty/core/Camera.cpp
@@ -3694,7 +3694,7 @@ CCamera::IsBoxVisible(CVUVECTOR *box, const CMatrix *mat)
 #ifdef GTA_PS2
 	TransformPoints(box, 8, *mat, box);
 #else
-	#ifdef FIX_BUGS
+	#if defined(FIX_BUGS) && !defined(DC_SH4)
 		for (i = 0; i < 8; i++)
 			box[i] = *mat * box[i];
 	#else
--- a/src/liberty/core/ControllerConfig.cpp
+++ b/src/liberty/core/ControllerConfig.cpp
@@ -2121,14 +2121,14 @@ wchar *CControllerConfigManager::GetControllerSettingTextWithOrderNumber(e_Contr
 				case VEHICLE_ACCELERATE:
 					for (int i = 0; (ActionText[i] = Dreamcast_RightTrigger[i]) != '\0' && i < iLimitCopy; i++);
 					break;
-				case VEHICLE_CHANGE_RADIO_STATION:
-					for (int i = 0; (ActionText[i] = Dreamcast_DPad_Left[i]) != '\0' && i < iLimitCopy; i++);
+				case VEHICLE_CHANGE_RADIO_STATION: // D-Pad Right to switch RADIO on Dreamcast
+					for (int i = 0; (ActionText[i] = Dreamcast_DPad_Right[i]) != '\0' && i < iLimitCopy; i++);
 					break;
 				case VEHICLE_BRAKE:
 					for (int i = 0; (ActionText[i] = Dreamcast_LeftTrigger[i]) != '\0' && i < iLimitCopy; i++);
 					break;
-				case TOGGLE_SUBMISSIONS:
-					for (int i = 0; (ActionText[i] = Dreamcast_DPad_Right[i]) != '\0' && i < iLimitCopy; i++);
+				case TOGGLE_SUBMISSIONS: // D-Pad Left to enable SUBMISSIONS on Dreamcast
+					for (int i = 0; (ActionText[i] = Dreamcast_DPad_Left[i]) != '\0' && i < iLimitCopy; i++);
 					break;
 				case GO_LEFT:
 					for (int i = 0; (ActionText[i] = Dreamcast_Stick_Left[i]) != '\0' && i < iLimitCopy; i++);
@@ -2199,11 +2199,11 @@ wchar *CControllerConfigManager::GetControllerSettingTextWithOrderNumber(e_Contr
 				case VEHICLE_TURRETRIGHT:
 					for (int i = 0; (ActionText[i] = Dreamcast_A[i]) != '\0' && i < iLimitCopy; i++);
 					break;
-				case VEHICLE_TURRETUP:
-					for (int i = 0; (ActionText[i] = Dreamcast_A[i]) != '\0' && i < iLimitCopy; i++);
+				case VEHICLE_TURRETUP: // Remark: VEHICLE_TURRETUP and VEHICLE_TURRETDOWN are used to shift your weight on a bike in Vice City
+					for (int i = 0; (ActionText[i] = Dreamcast_Stick_Up[i]) != '\0' && i < iLimitCopy; i++);
 					break;
-				case VEHICLE_TURRETDOWN:
-					for (int i = 0; (ActionText[i] = Dreamcast_A[i]) != '\0' && i < iLimitCopy; i++);
+				case VEHICLE_TURRETDOWN: // Remark: VEHICLE_TURRETUP and VEHICLE_TURRETDOWN are used to shift your weight on a bike in Vice City
+					for (int i = 0; (ActionText[i] = Dreamcast_Stick_Down[i]) != '\0' && i < iLimitCopy; i++);
 					break;
 				case CAMERA_CHANGE_VIEW_ALL_SITUATIONS:
 					for (int i = 0; (ActionText[i] = Dreamcast_DPad_Up[i]) != '\0' && i < iLimitCopy; i++);
@@ -2259,13 +2259,13 @@ wchar *CControllerConfigManager::GetControllerSettingTextWithOrderNumber(e_Contr
 				case VEHICLE_ACCELERATE:
 					for (int i = 0; (ActionText[i] = Dreamcast_RightTrigger[i]) != '\0' && i < iLimitCopy; i++);
 					break;
-				case VEHICLE_CHANGE_RADIO_STATION:
+				case VEHICLE_CHANGE_RADIO_STATION: // D-Pad Right to switch RADIO on Dreamcast
 					for (int i = 0; (ActionText[i] = Dreamcast_DPad_Right[i]) != '\0' && i < iLimitCopy; i++);
 					break;
 				case VEHICLE_BRAKE:
 					for (int i = 0; (ActionText[i] = Dreamcast_LeftTrigger[i]) != '\0' && i < iLimitCopy; i++);
 					break;
-				case TOGGLE_SUBMISSIONS:
+				case TOGGLE_SUBMISSIONS: // D-Pad Left to enable SUBMISSIONS on Dreamcast
 					for (int i = 0; (ActionText[i] = Dreamcast_DPad_Left[i]) != '\0' && i < iLimitCopy; i++);
 					break;
 				case GO_LEFT:
@@ -2397,14 +2397,14 @@ wchar *CControllerConfigManager::GetControllerSettingTextWithOrderNumber(e_Contr
 				case VEHICLE_ACCELERATE:
 					for (int i = 0; (ActionText[i] = Dreamcast_A[i]) != '\0' && i < iLimitCopy; i++);
 					break;
-				case VEHICLE_CHANGE_RADIO_STATION:
-					for (int i = 0; (ActionText[i] = Dreamcast_DPad_Left[i]) != '\0' && i < iLimitCopy; i++);
+				case VEHICLE_CHANGE_RADIO_STATION: // D-Pad Right to switch RADIO on Dreamcast
+					for (int i = 0; (ActionText[i] = Dreamcast_DPad_Right[i]) != '\0' && i < iLimitCopy; i++);
 					break;
 				case VEHICLE_BRAKE:
 					for (int i = 0; (ActionText[i] = Dreamcast_X[i]) != '\0' && i < iLimitCopy; i++);
 					break;
-				case TOGGLE_SUBMISSIONS:
-					for (int i = 0; (ActionText[i] = Dreamcast_DPad_Right[i]) != '\0' && i < iLimitCopy; i++);
+				case TOGGLE_SUBMISSIONS: // D-Pad Left to enable SUBMISSIONS on Dreamcast
+					for (int i = 0; (ActionText[i] = Dreamcast_DPad_Left[i]) != '\0' && i < iLimitCopy; i++);
 					break;
 				case GO_LEFT:
 					for (int i = 0; (ActionText[i] = Dreamcast_Stick_Left[i]) != '\0' && i < iLimitCopy; i++);
@@ -2542,8 +2542,8 @@ wchar *CControllerConfigManager::GetControllerSettingTextWithOrderNumber(e_Contr
 				case VEHICLE_BRAKE:
 					for (int i = 0; (ActionText[i] = PS2_Square[i]) != '\0' && i < iLimitCopy; i++);
 					break;
-				case TOGGLE_SUBMISSIONS:
-					for (int i = 0; (ActionText[i] = Dreamcast_DPad_Right[i]) != '\0' && i < iLimitCopy; i++);
+				case TOGGLE_SUBMISSIONS: // D-Pad Left to enable SUBMISSIONS on Dreamcast
+					for (int i = 0; (ActionText[i] = Dreamcast_DPad_Left[i]) != '\0' && i < iLimitCopy; i++);
 					break;
 				case GO_LEFT:
 					for (int i = 0; (ActionText[i] = LefAnalog_Left[i]) != '\0' && i < iLimitCopy; i++);
--- a/src/liberty/core/World.cpp
+++ b/src/liberty/core/World.cpp
@@ -1215,6 +1215,9 @@ CWorld::FindObjectsIntersectingAngledCollisionBox(const CColBox &boundingBox, co
 	const int32 nStartY = Max(GetSectorIndexY(fStartY), 0);
 	const int32 nEndX = Min(GetSectorIndexX(fEndX), NUMSECTORS_X - 1);
 	const int32 nEndY = Min(GetSectorIndexY(fEndY), NUMSECTORS_Y - 1);
+#ifdef DC_SH4
+    mat_load_transpose(matrix);
+#endif
 	for(int32 y = nStartY; y <= nEndY; y++) {
 		for(int32 x = nStartX; x <= nEndX; x++) {
 			CSector *pSector = GetSector(x, y);
@@ -1268,14 +1271,19 @@ CWorld::FindObjectsIntersectingAngledCollisionBoxSectorList(CPtrList &list, cons
                                                            int16 *nEntitiesFound, int16 maxEntitiesToFind,
                                                            CEntity **aEntities)
 {
-	for(CPtrNode *pNode = list.first; pNode; pNode = pNode->next) {
+    for(CPtrNode *pNode = list.first; pNode; pNode = pNode->next) {
 		CEntity *pEntity = (CEntity *)pNode->item;
 		if(pEntity->m_scanCode != GetCurrentScanCode()) {
 			pEntity->m_scanCode = GetCurrentScanCode();
 			CColSphere sphere;
 			CVector vecDistance = pEntity->GetPosition() - position;
 			sphere.radius = pEntity->GetBoundRadius();
+#ifndef DC_SH4
 			sphere.center = Multiply3x3(vecDistance, matrix);
+#else   // Transposed matrix was already loaded by the callee, so no need to reload!
+            mat_trans_normal3_nomod(vecDistance.x, vecDistance.y, vecDistance.z,
+                                    sphere.center.x, sphere.center.y, sphere.center.z);
+#endif
 			if(CCollision::TestSphereBox(sphere, boundingBox) && *nEntitiesFound < maxEntitiesToFind) {
 				if(aEntities) aEntities[*nEntitiesFound] = pEntity;
 				++*nEntitiesFound;
@@ -1444,13 +1452,27 @@ CWorld::CallOffChaseForAreaSectorListVehicles(CPtrList &list, float x1, float y1
 				pVehicle->AutoPilot.m_nTimeTempAction = CTimer::GetTimeInMilliseconds() + 2000;
 				CColModel *pColModel = pVehicle->GetColModel();
 				bool bInsideSphere = false;
+#ifdef DC_SH4
+                mat_load2(pVehicle->GetMatrix());
+#endif
 				for(int32 i = 0; i < pColModel->numSpheres; i++) {
+#ifndef DC_SH4
 					CVector pos = pVehicle->GetMatrix() * pColModel->spheres[i].center;
-					float fRadius = pColModel->spheres[i].radius;
+#else
+                    CVector pos;
+                    auto &center = pColModel->spheres[i].center;
+                    mat_trans_single3_nodiv_nomod(center.x, center.y, center.z,
+                                                  pos.x, pos.y, pos.z);
+#endif
+                    float fRadius = pColModel->spheres[i].radius;
 					if(pos.x + fRadius > x1 && pos.x - fRadius < x2 && pos.y + fRadius > y1 &&
-					   pos.y - fRadius < y2)
+					    pos.y - fRadius < y2) {
 						bInsideSphere = true;
 					// Maybe break the loop when bInsideSphere is set to true?
+#ifdef DC_SH4       // Don't see why not!
+                        break;
+#endif
+                    }
 				}
 				if(bInsideSphere) {
 					if(pVehicle->GetPosition().x <= (x1 + x2) * 0.5f)
--- a/src/liberty/core/ZoneCull.cpp
+++ b/src/liberty/core/ZoneCull.cpp
@@ -949,11 +949,15 @@ CCullZone::FindTestPoints()
 	if(ElementsY > 32) ElementsY = 32;
 	if(ElementsZ > 32) ElementsZ = 32;
 	Memsize = ElementsX * ElementsY * ElementsZ;
-	StepX = (maxx-minx)/(ElementsX-1);
-	StepY = (maxy-miny)/(ElementsY-1);
-	StepZ = (maxz-minz)/(ElementsZ-1);
+	StepX = Div<true, false>(maxx-minx, ElementsX-1);
+	StepY = Div<true, false>(maxy-miny, ElementsY-1);
+	StepZ = Div<true, false>(maxz-minz, ElementsZ-1);

+#ifndef DC_SH4
 	pMem = new uint8[Memsize];
+#else
+    pMem = reinterpret_cast<uint8 *>(alloca(Memsize));
+#endif
 	memset(pMem, 0, Memsize);

 	// indices of center
@@ -1496,13 +1500,28 @@ CCullZone::TestEntityVisibilityFromCullZone(CEntity *entity, float extraDist, CE
 	else
 		boundMaxZ += extraDist;

+#ifndef DC_SH4
 	CVector vecMin = entity->GetMatrix() * CVector(boundMinX, boundMinY, boundMinZ);
 	CVector vecMaxX = entity->GetMatrix() * CVector(boundMaxX, boundMinY, boundMinZ);
 	CVector vecMaxY = entity->GetMatrix() * CVector(boundMinX, boundMaxY, boundMinZ);
 	CVector vecMaxZ = entity->GetMatrix() * CVector(boundMinX, boundMinY, boundMaxZ);
-	CVector dirx = vecMaxX - vecMin;
-	CVector diry = vecMaxY - vecMin;
-	CVector dirz = vecMaxZ - vecMin;
+#else
+    mat_load2(entity->GetMatrix());
+
+    CVector vecMin, vecMaxX, vecMaxY, vecMaxZ;
+    mat_trans_single3_nodiv_nomod(boundMinX, boundMinY, boundMinZ,
+                                  vecMin.x, vecMin.y, vecMin.z);
+    mat_trans_single3_nodiv_nomod(boundMaxX, boundMinY, boundMinZ,
+                                  vecMaxX.x, vecMaxX.y, vecMaxX.z);
+    mat_trans_single3_nodiv_nomod(boundMinX, boundMaxY, boundMinZ,
+                                  vecMaxY.x, vecMaxY.y, vecMaxY.z);
+    mat_trans_single3_nodiv_nomod(boundMinX, boundMinY, boundMaxZ,
+                                  vecMaxZ.x, vecMaxZ.y, vecMaxZ.z);
+#endif
+
+    CVector dirx = vecMaxX - vecMin;
+    CVector diry = vecMaxY - vecMin;
+    CVector dirz = vecMaxZ - vecMin;

 	// If building intersects zone at all, it's visible
 	int x, y, z;
@@ -1520,22 +1539,30 @@ CCullZone::TestEntityVisibilityFromCullZone(CEntity *entity, float extraDist, CE

 	float distToZone = CalcDistToCullZone(entity->GetPosition().x, entity->GetPosition().y)/15.0f;
 	distToZone = Max(distToZone, 7.0f);
-	int numX = (boundMaxX - boundMinX)/distToZone + 2.0f;
-	int numY = (boundMaxY - boundMinY)/distToZone + 2.0f;
-	int numZ = (boundMaxZ - boundMinZ)/distToZone + 2.0f;
+    float invDistToZone = Invert<true, false>(distToZone);
+    int numX = (boundMaxX - boundMinX)*invDistToZone + 2.0f;
+	int numY = (boundMaxY - boundMinY)*invDistToZone + 2.0f;
+	int numZ = (boundMaxZ - boundMinZ)*invDistToZone + 2.0f;

-	float stepX = 1.0f/(numX-1);
-	float stepY = 1.0f/(numY-1);
-	float stepZ = 1.0f/(numZ-1);
+	float stepX = Invert<true, false>(numX-1);
+	float stepY = Invert<true, false>(numY-1);
+	float stepZ = Invert<true, false>(numZ-1);
 	float midX = (boundMaxX + boundMinX)/2.0f;
 	float midY = (boundMaxY + boundMinY)/2.0f;
 	float midZ = (boundMaxZ + boundMinZ)/2.0f;

 	// check both xy planes
-	for(int i = 0; i < NumTestPoints; i++){
+#ifndef DC_SH4
+    CVector mid = entity->GetMatrix() * CVector(midX, midY, midZ);
+#else
+    CVector mid;
+    mat_trans_single3_nodiv_nomod(midX, midY, midZ,
+                                  mid.x, mid.y, mid.z);
+#endif
+    mid.z += 0.1f;
+    for(int i = 0; i < NumTestPoints; i++){
 		CVector testPoint = aTestPoints[i];
-		CVector mid = entity->GetMatrix() * CVector(midX, midY, midZ);
-		mid.z += 0.1f;
+
 		if(DoThoroughLineTest(testPoint, mid, entity))
 			return true;

--- a/src/liberty/core/common.h
+++ b/src/liberty/core/common.h
@@ -82,8 +82,6 @@

 #define rwVENDORID_ROCKSTAR 0x0253F2

-__always_inline auto Max(auto a, auto b) { return ((a > b)? a : b); }
-__always_inline auto Min(auto a, auto b) { return ((a < b)? a : b); }

 // Use this to add const that wasn't there in the original code
 #define Const const
@@ -299,15 +297,6 @@ extern int strcasecmp(const char *str1, const char *str2);

 extern wchar *AllocUnicode(const char*src);

-template<typename T>
-__always_inline T Clamp(T v, auto low, auto high) {
-	return std::clamp(v, static_cast<T>(low), static_cast<T>(high));
-}
-
-__always_inline auto Clamp2(auto v, auto center, auto radius) {
-	return (v > center) ? Min(v, center + radius) : Max(v, center - radius);
-}
-
 #define SQR(x) ((x) * (x))
 __always_inline auto sq(auto x) { return SQR(x); }

@@ -418,15 +407,7 @@ template<int s, int t> struct check_size {
 #endif
 #define BIT(num)                         (1<<(num))

-#define ABS(a) std::abs(a)
+#define ABS(a) Abs(a)

-__always_inline auto norm(auto value, auto min, auto max) {
-	return (Clamp(value, min, max) - min) / (max - min);
-}
 // we use std::lerp now
 //#define lerp(norm, min, max) ( (norm) * ((max) - (min)) + (min) )
-
-#define STRINGIFY(x)                    #x
-#define STR(x)                          STRINGIFY(x)
-#define CONCAT_(x,y) x##y
-#define CONCAT(x,y) CONCAT_(x,y)
--- a/src/liberty/math/Matrix.cpp
+++ b/src/liberty/math/Matrix.cpp
@@ -4,7 +4,7 @@ CMatrix::CMatrix(CMatrix const &m)
 {
 	m_attachment = nil;
 	m_hasRwMatrix = false;
-	*this = m;
+	mat_copy(*this, m);
 }

 CMatrix::CMatrix(RwMatrix *matrix, bool owner)
@@ -54,36 +54,46 @@ CMatrix::Detach(void)
 void
 CMatrix::Update(void)
 {
+#ifndef DC_SH4
 	GetRight() = m_attachment->right;
 	GetForward() = m_attachment->up;
 	GetUp() = m_attachment->at;
 	GetPosition() = m_attachment->pos;
+#else
+    mat_copy(*this, *m_attachment);
+#endif
 }

 void
 CMatrix::UpdateRW(void)
 {
 	if (m_attachment) {
+#ifndef DC_SH4
 		m_attachment->right = GetRight();
 		m_attachment->up = GetForward();
 		m_attachment->at = GetUp();
 		m_attachment->pos = GetPosition();
 		RwMatrixUpdate(m_attachment);
+#else
+        mat_copy(*m_attachment, *this);
+#endif
 	}
 }

 void
 CMatrix::operator=(CMatrix const &rhs)
 {
-	memcpy(this, &rhs, sizeof(f));
+	mat_copy(*this, rhs);
+#ifndef DC_SH4
 	if (m_attachment)
 		UpdateRW();
+#endif
 }

 void
 CMatrix::CopyOnlyMatrix(const CMatrix &other)
 {
-	memcpy(this, &other, sizeof(f));
+	mat_copy(*this, other);
 }

 CMatrix &
@@ -99,6 +109,7 @@ CMatrix::operator+=(CMatrix const &rhs)
 void
 CMatrix::SetUnity(void)
 {
+#ifndef DC_SH4
 	rx = 1.0f;
 	ry = 0.0f;
 	rz = 0.0f;
@@ -111,6 +122,10 @@ CMatrix::SetUnity(void)
 	px = 0.0f;
 	py = 0.0f;
 	pz = 0.0f;
+#else
+    dc::mat_identity2();
+    dc::mat_store2(*this);
+#endif
 }

 void
@@ -130,6 +145,7 @@ CMatrix::ResetOrientation(void)
 void
 CMatrix::SetScale(float s)
 {
+#ifndef DC_SH4
 	rx = s;
 	ry = 0.0f;
 	rz = 0.0f;
@@ -145,11 +161,16 @@ CMatrix::SetScale(float s)
 	px = 0.0f;
 	py = 0.0f;
 	pz = 0.0f;
+#else
+    mat_set_scale(s);
+    mat_store2(*this);
+#endif
 }

 void
 CMatrix::SetTranslate(float x, float y, float z)
 {
+#ifndef DC_SH4
 	rx = 1.0f;
 	ry = 0.0f;
 	rz = 0.0f;
@@ -165,6 +186,10 @@ CMatrix::SetTranslate(float x, float y, float z)
 	px = x;
 	py = y;
 	pz = z;
+#else
+    mat_set_translation(x, y, z);
+    mat_store2(*this);
+#endif
 }

 void
@@ -224,34 +249,52 @@ CMatrix::SetRotateZOnly(float angle)
 void
 CMatrix::SetRotateX(float angle)
 {
+#ifndef DC_SH4
 	SetRotateXOnly(angle);
 	px = 0.0f;
 	py = 0.0f;
 	pz = 0.0f;
+#else
+    dc::mat_identity2();
+    dc::mat_apply_rotate_x(angle);
+    dc::mat_store2(*this);
+#endif
 }

-
 void
 CMatrix::SetRotateY(float angle)
 {
+#ifndef DC_SH4
 	SetRotateYOnly(angle);
 	px = 0.0f;
 	py = 0.0f;
 	pz = 0.0f;
+#else
+    dc::mat_identity2();
+    dc::mat_apply_rotate_y(angle);
+    dc::mat_store2(*this);
+#endif
 }

 void
 CMatrix::SetRotateZ(float angle)
 {
+#ifndef DC_SH4
 	SetRotateZOnly(angle);
 	px = 0.0f;
 	py = 0.0f;
 	pz = 0.0f;
+#else
+    dc::mat_identity2();
+    dc::mat_apply_rotate_z(angle);
+    dc::mat_store2(*this);
+#endif
 }

 void
 CMatrix::SetRotate(float xAngle, float yAngle, float zAngle)
 {
+#if 1
 	auto [sX, cX] = SinCos(xAngle);
 	auto [sY, cY] = SinCos(yAngle);
 	auto [sZ, cZ] = SinCos(zAngle);
@@ -271,15 +314,19 @@ CMatrix::SetRotate(float xAngle, float yAngle, float zAngle)
 	px = 0.0f;
 	py = 0.0f;
 	pz = 0.0f;
+#else
+    dc::mat_set_rotate(xAngle, yAngle, zAngle);
+    dc::mat_store2(*this);
+#endif
 }

 void
 CMatrix::RotateX(float x)
 {
-#if 0 && defined(DC_SH4) // this is bugged and does not yield correct results
-	mat_load(reinterpret_cast<matrix_t *>(this));
-	mat_rotate_x(x);
-	mat_store(reinterpret_cast<matrix_t *>(this));
+#if 0// this is bugged and does not yield correct results
+    dc::mat_set_rotate_x(x);
+    mat_apply(*this);
+    dc::mat_store2(*this);
 #else
 	auto [s, c] = SinCos(x);

@@ -306,10 +353,10 @@ CMatrix::RotateX(float x)
 void
 CMatrix::RotateY(float y)
 {
-#if 0 && defined(DC_SH4) // this is bugged and does not yield correct results
-	mat_load(reinterpret_cast<matrix_t *>(this));
-	mat_rotate_y(y);
-	mat_store(reinterpret_cast<matrix_t *>(this));
+#if 0 // this is bugged and does not yield correct results
+    dc::mat_set_rotate_y(y);
+    mat_apply(*this);
+    dc::mat_store2(*this);
 #else
 	auto [s, c] = SinCos(y);

@@ -336,10 +383,10 @@ CMatrix::RotateY(float y)
 void
 CMatrix::RotateZ(float z)
 {
-#if 0 && defined(DC_SH4) // this is bugged and does not yield correct results
-	mat_load(reinterpret_cast<matrix_t *>(this));
-	mat_rotate_z(z);
-	mat_store(reinterpret_cast<matrix_t *>(this));
+#if 0// this is bugged and does not yield correct results
+       dc::mat_set_rotate_z(z);
+       mat_apply(*this);
+       dc::mat_store2(*this);
 #else	
 	auto [s, c] = SinCos(z);

@@ -366,10 +413,10 @@ CMatrix::RotateZ(float z)
 void
 CMatrix::Rotate(float x, float y, float z)
 {
-#if 0 && defined(DC_SH4) // this is bugged and does not yield correct results
-	mat_load(reinterpret_cast<matrix_t *>(this));
-	mat_rotate(x, y, z);
-	mat_store(reinterpret_cast<matrix_t *>(this));
+#if 0 // this is bugged and does not yield correct results
+    dc::mat_set_rotate(x, y, z);
+    mat_apply(*this);
+    dc::mat_store2(*this);
 #else
 	auto [sX, cX] = SinCos(x);
 	auto [sY, cY] = SinCos(y);
@@ -449,65 +496,13 @@ CMatrix::Reorthogonalise(void)
 	f = CrossProduct(u, r);
 }

-#ifdef DC_SH4
-static __always_inline void MATH_Load_Matrix_Product(const matrix_t* matrix1, const matrix_t* matrix2)
-{
-    unsigned int prefetch_scratch;
-
-    asm volatile (
-        "mov %[bmtrx], %[pref_scratch]\n\t" // (MT)
-        "add #32, %[pref_scratch]\n\t" // offset by 32 (EX - flow dependency, but 'add' is actually parallelized since 'mov Rm, Rn' is 0-cycle)
-        "fschg\n\t" // switch fmov to paired moves (note: only paired moves can access XDn regs) (FE)
-        "pref @%[pref_scratch]\n\t" // Get a head start prefetching the second half of the 64-byte data (LS)
-        // back matrix
-        "fmov.d @%[bmtrx]+, XD0\n\t" // (LS)
-        "fmov.d @%[bmtrx]+, XD2\n\t"
-        "fmov.d @%[bmtrx]+, XD4\n\t"
-        "fmov.d @%[bmtrx]+, XD6\n\t"
-        "pref @%[fmtrx]\n\t" // prefetch fmtrx now while we wait (LS)
-        "fmov.d @%[bmtrx]+, XD8\n\t" // bmtrx prefetch should work for here
-        "fmov.d @%[bmtrx]+, XD10\n\t"
-        "fmov.d @%[bmtrx]+, XD12\n\t"
-        "mov %[fmtrx], %[pref_scratch]\n\t" // (MT)
-        "add #32, %[pref_scratch]\n\t" // store offset by 32 in r0 (EX - flow dependency, but 'add' is actually parallelized since 'mov Rm, Rn' is 0-cycle)
-        "fmov.d @%[bmtrx], XD14\n\t"
-        "pref @%[pref_scratch]\n\t" // Get a head start prefetching the second half of the 64-byte data (LS)
-        // front matrix
-        // interleave loads and matrix multiply 4x4
-        "fmov.d @%[fmtrx]+, DR0\n\t"
-        "fmov.d @%[fmtrx]+, DR2\n\t"
-        "fmov.d @%[fmtrx]+, DR4\n\t" // (LS) want to issue the next one before 'ftrv' for parallel exec
-        "ftrv XMTRX, FV0\n\t" // (FE)
-
-        "fmov.d @%[fmtrx]+, DR6\n\t"
-        "fmov.d @%[fmtrx]+, DR8\n\t"
-        "ftrv XMTRX, FV4\n\t"
-
-        "fmov.d @%[fmtrx]+, DR10\n\t"
-        "fmov.d @%[fmtrx]+, DR12\n\t"
-        "ftrv XMTRX, FV8\n\t"
-
-        "fmov.d @%[fmtrx], DR14\n\t" // (LS, but this will stall 'ftrv' for 3 cycles)
-        "fschg\n\t" // switch back to single moves (and avoid stalling 'ftrv') (FE)
-        "ftrv XMTRX, FV12\n\t" // (FE)
-        // Save output in XF regs
-        "frchg\n"
-        : [bmtrx] "+&r" ((unsigned int)matrix1), [fmtrx] "+r" ((unsigned int)matrix2), [pref_scratch] "=&r" (prefetch_scratch) // outputs, "+" means r/w, "&" means it's written to before all inputs are consumed
-        : // no inputs
-        : "fr0", "fr1", "fr2", "fr3", "fr4", "fr5", "fr6", "fr7", "fr8", "fr9", "fr10", "fr11", "fr12", "fr13", "fr14", "fr15" // clobbers (GCC doesn't know about back bank, so writing to it isn't clobbered)
-    );
-}
-#endif
-
 CMatrix
 operator*(const CMatrix &m1, const CMatrix &m2)
 {
 	// TODO: VU0 code
 	CMatrix out;
-#if defined(RW_DC)
-	mat_load(reinterpret_cast<const matrix_t *>(&m1));
-	mat_apply(reinterpret_cast<const matrix_t *>(&m2));
-	mat_store(reinterpret_cast<matrix_t *>(&out));
+#ifdef DC_SH4
+	mat_mult(out, m1, m2);
 #else
 	out.rx = m1.rx * m2.rx + m1.fx * m2.ry + m1.ux * m2.rz;
 	out.ry = m1.ry * m2.rx + m1.fy * m2.ry + m1.uy * m2.rz;
--- a/src/liberty/math/Matrix.h
+++ b/src/liberty/math/Matrix.h
@@ -1,5 +1,7 @@
 #pragma once

+#include "rwdc_common.h"
+
 class alignas(8) CMatrix
 {
 public:
@@ -27,6 +29,8 @@ public:
 		SetScale(scale);
 	}
 	~CMatrix(void);
+	operator matrix_t *() { return reinterpret_cast<matrix_t *>(this); }
+	operator const matrix_t *() const { return reinterpret_cast<const matrix_t *>(this); }
 	void Attach(RwMatrix *matrix, bool owner = false);
 	void AttachRW(RwMatrix *matrix, bool owner = false);
 	void Detach(void);
@@ -102,15 +106,24 @@ CMatrix Invert(const CMatrix &matrix);
 CMatrix operator*(const CMatrix &m1, const CMatrix &m2);
 inline CVector MultiplyInverse(const CMatrix &mat, const CVector &vec)
 {
+#ifndef DC_SH4
 	CVector v(vec.x - mat.px, vec.y - mat.py, vec.z - mat.pz);
 	return CVector(
 		mat.rx * v.x + mat.ry * v.y + mat.rz * v.z,
 		mat.fx * v.x + mat.fy * v.y + mat.fz * v.z,
 		mat.ux * v.x + mat.uy * v.y + mat.uz * v.z);
+#else
+    register float x asm(KOS_FPARG(0)) = vec.x - mat.px;
+    register float y asm(KOS_FPARG(1)) = vec.y - mat.py;
+    register float z asm(KOS_FPARG(2)) = vec.z - mat.pz;
+	return CVector(
+		fipr(x, y, z, 0.0f, mat.rx, mat.ry, mat.rz, 0.0f),
+		fipr(x, y, z, 0.0f, mat.fx, mat.fy, mat.fz, 0.0f),
+		fipr(x, y, z, 0.0f, mat.ux, mat.uy, mat.uz, 0.0f)
+	);
+#endif
 }

-
-
 class CCompressedMatrixNotAligned
 {
 	CVector m_vecPos;
--- a/src/liberty/math/Quaternion.cpp
+++ b/src/liberty/math/Quaternion.cpp
@@ -39,10 +39,14 @@ CQuaternion::Slerp(const CQuaternion &q1, const CQuaternion &q2, float theta, fl
 void
 CQuaternion::Multiply(const CQuaternion &q1, const CQuaternion &q2)
 {
+#ifndef DC_SH4
 	x = (q2.z * q1.y) - (q1.z * q2.y) + (q1.x * q2.w) + (q2.x * q1.w);
 	y = (q2.x * q1.z) - (q1.x * q2.z) + (q1.y * q2.w) + (q2.y * q1.w);
 	z = (q2.y * q1.x) - (q1.y * q2.x) + (q1.z * q2.w) + (q2.z * q1.w);
 	w = (q2.w * q1.w) - (q2.x * q1.x) - (q2.y * q1.y) - (q2.z * q1.z);
+#else
+	quat_mult(*this, q1, q2);
+#endif
 }

 void
@@ -51,9 +55,16 @@ CQuaternion::Get(RwV3d *axis, float *angle)
 	*angle = Acos(w);
 	float s = Sin(*angle);

+#ifndef DC_SH4
 	axis->x = x * (1.0f / s);
 	axis->y = y * (1.0f / s);
 	axis->z = z * (1.0f / s);
+#else
+	float invS = dc::Invert<true, false>(s);
+	axis->x = x * invS;
+	axis->y = y * invS;
+	axis->z = z * invS;
+#endif
 }

 void
@@ -104,7 +115,7 @@ CQuaternion::Set(const RwMatrix &matrix)
 	if (f >= 0.0f) {
 		s = Sqrt(f + 1.0f);
 		w = 0.5f * s;
-		m = 0.5f / s;
+		m = Div<true, false>(0.5f, s);
 		x = (matrix.up.z - matrix.at.y) * m;
 		y = (matrix.at.x - matrix.right.z) * m;
 		z = (matrix.right.y - matrix.up.x) * m;
@@ -115,7 +126,7 @@ CQuaternion::Set(const RwMatrix &matrix)
 	if (f >= 0.0f) {
 		s = Sqrt(f + 1.0f);
 		x = 0.5f * s;
-		m = 0.5f / s;
+		m = Div<true, false>(0.5f, s);
 		y = (matrix.up.x + matrix.right.y) * m;
 		z = (matrix.at.x + matrix.right.z) * m;
 		w = (matrix.up.z - matrix.at.y) * m;
@@ -126,7 +137,7 @@ CQuaternion::Set(const RwMatrix &matrix)
 	if (f >= 0.0f) {
 		s = Sqrt(f + 1.0f);
 		y = 0.5f * s;
-		m = 0.5f / s;
+		m = Div<true, false>(0.5f, s);
 		w = (matrix.at.x - matrix.right.z) * m;
 		x = (matrix.up.x - matrix.right.y) * m;
 		z = (matrix.at.y + matrix.up.z) * m;
--- a/src/liberty/math/Quaternion.h
+++ b/src/liberty/math/Quaternion.h
@@ -1,5 +1,8 @@
 #pragma once

+#include "src/common_defines.h"
+#include "rwdc_common.h"
+
 // TODO: actually implement this
 class CQuaternion
 {
@@ -8,13 +11,23 @@ public:
 	CQuaternion(void) {}
 	CQuaternion(float x, float y, float z, float w) : x(x), y(y), z(z), w(w) {}

-	float Magnitude(void) const { return Sqrt(MagnitudeSqr()); }
+	operator quaternion_t *() { return reinterpret_cast<quaternion_t *>(this); }
+	operator const quaternion_t *() const { return reinterpret_cast<const quaternion_t *>(this); }
+	operator quaternion_t &() { return *reinterpret_cast<quaternion_t *>(this); }
+	operator const quaternion_t &() const { return *reinterpret_cast<const quaternion_t *>(this); }
+	float Magnitude(void) const {
+#ifndef DC_SH4
+		return Sqrt(x*x + y*y + z*z + w*w);
+#else
+		return Sqrt(fipr_magnitude_sqr(x, y, z, w));
+#endif
+	}
 	float MagnitudeSqr(void) const {
-	#ifdef DC_SH4
-		return fipr_magnitude_sqr(x, y, z, w);
-	#else
+#ifndef DC_SH4
 		return x*x + y*y + z*z + w*w;
-	#endif
+#else
+		return fipr_magnitude_sqr(x, y, z, w);
+#endif
 	}
 	void Normalise(void);
 	void Multiply(const CQuaternion &q1, const CQuaternion &q2);
@@ -49,10 +62,11 @@ public:
 	}

 	const CQuaternion &operator/=(float right) {
-		x /= right;
-		y /= right;
-		z /= right;
-		w /= right;
+		right = dc::Invert<false>(right);
+		x *= right;
+		y *= right;
+		z *= right;
+		w *= right;
 		return *this;
 	}

@@ -72,7 +86,11 @@ public:
 inline float
 DotProduct(const CQuaternion &q1, const CQuaternion &q2)
 {
+#ifndef DC_SH4
 	return q1.x*q2.x + q1.y*q2.y + q1.z*q2.z + q1.w*q2.w;
+#else
+	return fipr(q1.x, q1.y, q1.z, q1.w, q2.x, q2.y, q2.z, q2.w);
+#endif
 }

 inline CQuaternion operator+(const CQuaternion &left, const CQuaternion &right)
@@ -97,5 +115,6 @@ inline CQuaternion operator*(float left, const CQuaternion &right)

 inline CQuaternion operator/(const CQuaternion &left, float right)
 {
-	return CQuaternion(left.x / right, left.y / right, left.z / right, left.w / right);
+	right = Invert<false>(right);
+	return CQuaternion(left.x * right, left.y * right, left.z * right, left.w * right);
 }
--- a/src/liberty/math/Vector.cpp
+++ b/src/liberty/math/Vector.cpp
@@ -27,32 +27,34 @@ CrossProduct(const CVector &v1, const CVector &v2)
 CVector
 Multiply3x3(const CMatrix &mat, const CVector &vec)
 {
-#ifdef DC_SH4
-    register float __x __asm__("fr12") = vec.x;
-    register float __y __asm__("fr13") = vec.y;
-    register float __z __asm__("fr14") = vec.z;
-    register float __w __asm__("fr15") = 0.0f;
-	
-	mat_load(reinterpret_cast<matrix_t *>(const_cast<CMatrix *>(&mat)));
-
-	asm volatile( "ftrv  xmtrx, fv12\n"
-                : "=f" (__x), "=f" (__y), "=f" (__z), "=f" (__w)
-                : "0" (__x), "1" (__y), "2" (__z), "3" (__w) );
-    return { __x, __y, __z };
-#else
+#ifndef DC_SH4
 	// TODO: VU0 code
 	return CVector(mat.rx * vec.x + mat.fx * vec.y + mat.ux * vec.z,
 	               mat.ry * vec.x + mat.fy * vec.y + mat.uy * vec.z,
 	               mat.rz * vec.x + mat.fz * vec.y + mat.uz * vec.z);
+#else
+	CVector out;
+	dc::mat_load2(mat);
+	mat_trans_normal3_nomod(vec.x, vec.y, vec.z,
+	                        out.x, out.y, out.z);
+    return out;
 #endif
 }

 CVector
 Multiply3x3(const CVector &vec, const CMatrix &mat)
 {
+#ifndef DC_SH4
 	return CVector(mat.rx * vec.x + mat.ry * vec.y + mat.rz * vec.z,
 	               mat.fx * vec.x + mat.fy * vec.y + mat.fz * vec.z,
 	               mat.ux * vec.x + mat.uy * vec.y + mat.uz * vec.z);
+#else
+	CVector out;
+	dc::mat_load_transpose(mat);
+	mat_trans_normal3_nomod(vec.x, vec.y, vec.z,
+	                        out.x, out.y, out.z);
+	return out;
+#endif
 }

 CVector
@@ -60,7 +62,7 @@ operator*(const CMatrix &mat, const CVector &vec)
 {
 #ifdef DC_SH4
 	CVector out;
-	mat_load(reinterpret_cast<matrix_t *>(const_cast<CMatrix *>(&mat)));
+	dc::mat_load2(mat);
 	mat_trans_single3_nodiv_nomod(vec.x, vec.y, vec.z, out.x, out.y, out.z);
 	return out;
 #else
--- a/src/liberty/math/Vector.h
+++ b/src/liberty/math/Vector.h
@@ -1,5 +1,7 @@
 #pragma once 

+#include "maths.h"
+
 class CVector : public RwV3d
 {
 public:
@@ -19,22 +21,22 @@ public:
  }
  // (0,1,0) means no rotation. So get right vector and its atan
  __always_inline float Heading(void) const { return Atan2(-x, y); }
-  __always_inline float Magnitude(void) const { 
+  __always_inline float Magnitude(void) const {
 #ifdef DC_SH4
    float w;
    vec3f_length(x, y, z, w);
    return w;
 #else
-    return Sqrt(x*x + y*y + z*z); 
+    return Sqrt(x*x + y*y + z*z);
 #endif
 }
-  __always_inline float MagnitudeSqr(void) const { 
+  __always_inline float MagnitudeSqr(void) const {
 #ifdef DC_SH4
-    return fipr_magnitude_sqr(x, y,z, 0.0f); 
+    return fipr_magnitude_sqr(x, y, z, 0.0f);
 #else
    return x*x + y*y + z*z;
 #endif
-} 
+}
  __always_inline float Magnitude2D(void) const { return Sqrt(x*x + y*y); }
  float MagnitudeSqr2D(void) const { return x*x + y*y; }
  void Normalise(void);
@@ -68,7 +70,7 @@ public:
  }

  const CVector &operator/=(float right) {
-    right = Invert(right);
+    right = Invert<true, true>(right);
    x *= right;
    y *= right;
    z *= right;
@@ -112,7 +114,8 @@ inline CVector operator*(float left, const CVector &right)

 inline CVector operator/(const CVector &left, float right)
 {
-  return CVector(left.x / right, left.y / right, left.z / right);
+  right = Invert<true, true>(right);
+  return CVector(left.x * right, left.y * right, left.z * right);
 }

 __always_inline float
--- a/src/liberty/math/VuVector.h
+++ b/src/liberty/math/VuVector.h
@@ -44,7 +44,7 @@ __always_inline void TransformPoint(CVuVector &out, const CMatrix &mat, const CV
 		sqc2    vf06,0x0(%0)\n\
 		": : "r" (&out) , "r" (&mat) ,"r" (&in): "memory");
 #elif defined(DC_SH4)
-	mat_load(reinterpret_cast<matrix_t *>(const_cast<CMatrix *>(&mat)));
+	dc::mat_load2(mat);
 	mat_trans_single3_nodiv_nomod(in.x, in.y, in.z, out.x, out.y, out.z);
 #else
 	out = mat * in;
@@ -71,7 +71,7 @@ __always_inline void TransformPoint(CVuVector &out, const CMatrix &mat, const Rw
 		sqc2    vf06,0x0(%0)\n\
 		": : "r" (&out) , "r" (&mat) ,"r" (&in): "memory");
 #elif defined(DC_SH4)
-	mat_load(reinterpret_cast<matrix_t *>(const_cast<CMatrix *>(&mat)));
+	dc::mat_load2(mat);
 	mat_trans_single3_nodiv_nomod(in.x, in.y, in.z, out.x, out.y, out.z);
 #else
 	out = mat * in;
@@ -108,7 +108,7 @@ __always_inline void TransformPoints(CVuVector *out, int n, const CMatrix &mat,
 		bnez	%1,1b\n\
 		": : "r" (out) , "r" (n), "r" (&mat), "r" (in), "r" (stride): "memory");
 #elif defined(DC_SH4)
-    mat_load(reinterpret_cast<matrix_t *>(const_cast<CMatrix *>(&mat)));
+    dc::mat_load2(mat);
 	while(n--) {
 		mat_trans_single3_nodiv_nomod(in->x, in->y, in->z, out->x, out->y, out->z);
 		in = reinterpret_cast<const RwV3d *>(reinterpret_cast<const uint8_t *>(in) + stride);
--- a/src/liberty/math/maths.h
+++ b/src/liberty/math/maths.h
@@ -1,73 +1,6 @@
 #pragma once

 #include "src/common_defines.h"
+#include "rwdc_common.h"

-#include <tuple>
-#include <dc/matrix.h>
-
-#ifdef DC_SH4
-
-#define mat_trans_nodiv_nomod(x, y, z, x2, y2, z2, w2) do { \
-        register float __x __asm__("fr12") = (x); \
-        register float __y __asm__("fr13") = (y); \
-        register float __z __asm__("fr14") = (z); \
-        register float __w __asm__("fr15") = 1.0f; \
-        __asm__ __volatile__( "ftrv  xmtrx, fv12\n" \
-                              : "=f" (__x), "=f" (__y), "=f" (__z), "=f" (__w) \
-                              : "0" (__x), "1" (__y), "2" (__z), "3" (__w) ); \
-        x2 = __x; y2 = __y; z2 = __z; w2 = __w; \
-    } while(false)
-
-#define mat_trans_w_nodiv_nomod(x, y, z, w) do { \
-        register float __x __asm__("fr12") = (x); \
-        register float __y __asm__("fr13") = (y); \
-        register float __z __asm__("fr14") = (z); \
-        register float __w __asm__("fr15") = 1.0f; \
-        __asm__ __volatile__( "ftrv  xmtrx, fv12\n" \
-                              : "=f" (__x), "=f" (__y), "=f" (__z), "=f" (__w) \
-                              : "0" (__x), "1" (__y), "2" (__z), "3" (__w) ); \
-        w = __w; \
-    } while(false)
-
-__always_inline float Fmac(float a, float b, float c) {
-    asm volatile ("fmac fr0, %[floatb], %[floatc]\n"
-        : [floatc] "+f" (c) : "w" (a), [floatb] "f" (b) : );
-    return c;
-}
-
-#else
-
-#define mat_trans_nodiv_nomod(x_, y_, z_, x2, y2, z2, w2) do { \
-		vector_t tmp = { x_, y_, z_, 1.0f }; \
-		mat_transform(&tmp, &tmp, 1, 0); \
-		x2 = tmp.x; y2 = tmp.y; z2 = tmp.z; w2 = tmp.w; \
-	} while(false)
-
-#define mat_trans_w_nodiv_nomod(x_, y_, z_, w_) do { \
-		vector_t tmp = { x_, y_, z_, 1.0f }; \
-		mat_transform(&tmp, &tmp, 1, 0); \
-		w_ = tmp.w; \
-	} while(false)
-
-__always_inline float Fmac(float a, float b, float c) { return a * b + c; }
-
-#endif
-
-__always_inline float Sin(float x) { return __builtin_sinf(x); }
-__always_inline float Cos(float x) { return __builtin_cosf(x); }
-__always_inline auto SinCos(float x) { return std::pair { Sin(x), Cos(x) }; }
-__always_inline float Tan(float x) { return __builtin_tanf(x); }
-__always_inline float Abs(float x) { return __builtin_fabsf(x); }
-__always_inline float Sqrt(float x) { return __builtin_sqrtf(x); }
-__always_inline float RecipSqrt(float x) { return 1.0f / __builtin_sqrtf(x); }
-__always_inline float Asin(float x) { return __builtin_asinf(x); }
-__always_inline float Acos(float x) { return __builtin_acosf(x); }
-__always_inline float Atan(float x) { return __builtin_atanf(x); }
-__always_inline float Atan2(float y, float x) { return __builtin_atan2f(y, x); }
-__always_inline float RecipSqrt(float x, float y) { return x / __builtin_sqrtf(y); /*y = RecipSqrt(y); return x * y * y;*/ }
-__always_inline float Pow(float x, float y) { return __builtin_powf(x, y); }
-__always_inline float Floor(float x) { return __builtin_floorf(x); }
-__always_inline float Ceil(float x) { return __builtin_ceilf(x); }
-__always_inline float Invert(float x) { return (((x) < 0.0f)? -1.0f : 1.0f) * RecipSqrt((x) * (x)); } 
-__always_inline float Div(float x, float y) { return x * Invert(y); }
-__always_inline float Lerp(float a, float b, float t) { return Fmac(t, (b - a), a); }
+using namespace dc;
--- a/src/liberty/renderer/Coronas.cpp
+++ b/src/liberty/renderer/Coronas.cpp
@@ -300,9 +300,9 @@ CCoronas::Render(void)


 			if(aCoronas[i].fadeAlpha && spriteCoors.z < aCoronas[i].drawDist){
-				float recipz = 1.0f/spriteCoors.z;
+				float recipz = dc::Invert<true, false>(spriteCoors.z);
 				float fadeDistance = aCoronas[i].drawDist / 2.0f;
-				float distanceFade = spriteCoors.z < fadeDistance ? 1.0f : 1.0f - (spriteCoors.z - fadeDistance)/fadeDistance;
+				float distanceFade = spriteCoors.z < fadeDistance ? 1.0f : 1.0f - dc::Div<true, false>((spriteCoors.z - fadeDistance), fadeDistance);
 				int totalFade = aCoronas[i].fadeAlpha * distanceFade;

 				if(aCoronas[i].LOScheck)
@@ -313,6 +313,7 @@ CCoronas::Render(void)
 				// render corona itself
 				if(aCoronas[i].texture){
 					float fogscale = CWeather::Foggyness*Min(spriteCoors.z, 40.0f)/40.0f + 1.0f;
+                    float invFogScale = dc::Invert<true, false>(fogscale);
 					if(CCoronas::aCoronas[i].id == SUN_CORE)
 						spriteCoors.z = 0.95f * RwCameraGetFarClipPlane(Scene.camera);
 					RwRenderStateSet(rwRENDERSTATETEXTURERASTER, RwTextureGetRaster(aCoronas[i].texture));
@@ -328,9 +329,9 @@ CCoronas::Render(void)
 						CSprite::RenderOneXLUSprite(spriteCoors.x, spriteCoors.y, spriteCoors.z,
 							spritew * aCoronas[i].size * wscale,
 							spriteh * aCoronas[i].size * fogscale * hscale,
-							CCoronas::aCoronas[i].red / fogscale,
-							CCoronas::aCoronas[i].green / fogscale,
-							CCoronas::aCoronas[i].blue / fogscale,
+							CCoronas::aCoronas[i].red * invFogScale,
+							CCoronas::aCoronas[i].green * invFogScale,
+							CCoronas::aCoronas[i].blue * invFogScale,
 							totalFade,
 							recipz,
 							255);
@@ -339,9 +340,9 @@ CCoronas::Render(void)
 							spriteCoors.x, spriteCoors.y, spriteCoors.z,
 							spritew * aCoronas[i].size * fogscale,
 							spriteh * aCoronas[i].size * fogscale,
-							CCoronas::aCoronas[i].red / fogscale,
-							CCoronas::aCoronas[i].green / fogscale,
-							CCoronas::aCoronas[i].blue / fogscale,
+							CCoronas::aCoronas[i].red * invFogScale,
+							CCoronas::aCoronas[i].green * invFogScale,
+							CCoronas::aCoronas[i].blue * invFogScale,
 							totalFade,
 							recipz,
 							20.0f * recipz,
@@ -365,7 +366,7 @@ CCoronas::Render(void)
 							(spriteCoors.x - (screenw/2)) * flare->position + (screenw/2),
 							(spriteCoors.y - (screenh/2)) * flare->position + (screenh/2),
 							spriteCoors.z,
-							4.0f*flare->size * spritew/spriteh,
+							4.0f*flare->size * dc::Div<true, false>(spritew, spriteh),
 							4.0f*flare->size,
 							(flare->red * aCoronas[i].red)>>8,
 							(flare->green * aCoronas[i].green)>>8,
@@ -480,9 +481,9 @@ CCoronas::RenderReflections(void)
 					drawDist = Min(drawDist, 55.0f);
 					if(spriteCoors.z < drawDist){
 						float fadeDistance = drawDist / 2.0f;
-						float distanceFade = spriteCoors.z < fadeDistance ? 1.0f : 1.0f - (spriteCoors.z - fadeDistance)/fadeDistance;
+						float distanceFade = spriteCoors.z < fadeDistance ? 1.0f : 1.0f - Div<true, false>((spriteCoors.z - fadeDistance), fadeDistance);
 						distanceFade = Clamp(distanceFade, 0.0f, 1.0f);
-						float recipz = 1.0f/RwCameraGetNearClipPlane(Scene.camera);
+						float recipz = dc::Invert<true, false>(RwCameraGetNearClipPlane(Scene.camera));
 						float heightFade = (20.0f - aCoronas[i].heightAboveRoad)/20.0f;
 						int intensity = distanceFade*heightFade * 230.0 * CWeather::WetRoads;

@@ -606,7 +607,9 @@ CEntity::ProcessLightsForEntity(void)
 	flashTimer1 = 0;
 	flashTimer2 = 0;
 	flashTimer3 = 0;
-
+#ifdef DC_SH4
+    dc:mat_load2(GetMatrix());
+#endif
 	n = CModelInfo::GetModelInfo(GetModelIndex())->GetNum2dEffects();
 	for(i = 0; i < n; i++, flashTimer1 += 0x80, flashTimer2 += 0x100, flashTimer3 += 0x200){
 		effect = CModelInfo::GetModelInfo(GetModelIndex())->Get2dEffect(i);
@@ -614,8 +617,12 @@ CEntity::ProcessLightsForEntity(void)
 		if(effect->type != EFFECT_LIGHT)
 			continue;

+#ifndef DC_SH4
 		pos = GetMatrix() * effect->pos;
-
+#else
+        mat_trans_single3_nodiv_nomod(effect->pos.x, effect->pos.y, effect->pos.z,
+                                      pos.x, pos.y, pos.z);
+#endif
 		lightOn = false;
 		lightFlickering = false;
 		switch(effect->light.lightType){
--- a/src/liberty/renderer/PointLights.cpp
+++ b/src/liberty/renderer/PointLights.cpp
@@ -80,7 +80,7 @@ CPointLights::GenerateLightsAffectingObject(Const CVector *objCoors)
 			distance = dist.Magnitude();
 			if(distance < radius){

-				float distNorm = distance/radius;
+				float distNorm = Div<true, false>(distance, radius);
 				if(aLights[i].type == LIGHT_DARKEN){
 					// darken the object the closer it is
 					ret *= distNorm;
--- a/src/liberty/renderer/Renderer.cpp
+++ b/src/liberty/renderer/Renderer.cpp
@@ -1315,8 +1315,9 @@ CalcNewDelta(RwV2d *a, RwV2d *b)
 #define TOINT(x) ((int)(x))
 #endif

+template<typename F>
 void
-CRenderer::ScanSectorPoly(RwV2d *poly, int32 numVertices, void (*scanfunc)(CPtrList *))
+CRenderer::ScanSectorPoly(RwV2d *poly, int32 numVertices, F&& scanfunc)
 {
 	float miny, maxy;
 	int y, yend;
--- a/src/liberty/renderer/Renderer.h
+++ b/src/liberty/renderer/Renderer.h
@@ -90,7 +90,8 @@ public:
 	static void ConstructRenderList(void);
 	static void ScanWorld(void);
 	static void RequestObjectsInFrustum(void);
-	static void ScanSectorPoly(RwV2d *poly, int32 numVertices, void (*scanfunc)(CPtrList *));
+    template<typename F>
+	static void ScanSectorPoly(RwV2d *poly, int32 numVertices, F &&scanfunc);
 	static void ScanBigBuildingList(CPtrList &list);
 	static void ScanSectorList(CPtrList *lists);
 	static void ScanSectorList_Priority(CPtrList *lists);
--- a/src/miami/collision/Collision.cpp
+++ b/src/miami/collision/Collision.cpp
@@ -1473,7 +1473,7 @@ CCollision::ProcessLineOfSight(const CColLine &line,
 		point.point = matrix * point.point;
 		point.normal = Multiply3x3(matrix, point.normal);
 #else
-		mat_load(reinterpret_cast<matrix_t *>(const_cast<CMatrix *>(&matrix)));
+		dc::mat_load2(matrix);
 		mat_trans_single3_nodiv(point.point.x,
 		                        point.point.y,
 		                        point.point.z);
@@ -1653,7 +1653,7 @@ CCollision::ProcessVerticalLine(const CColLine &line,
 		point.point = matrix * point.point;
 		point.normal = Multiply3x3(matrix, point.normal);
 #else
-		mat_load(reinterpret_cast<matrix_t *>(const_cast<CMatrix *>(&matrix)));
+		dc::mat_load2(matrix);
 		mat_trans_single3_nodiv(point.point.x,
 		                        point.point.y,
 		                        point.point.z);
@@ -2027,8 +2027,7 @@ CCollision::ProcessColModels(const CMatrix &matrixA, CColModel &modelA,
 #ifndef DC_SH4
 	matAB *= matrixA;
 #else
-	mat_load(reinterpret_cast<const matrix_t*>(&matAB));
-	mat_apply(reinterpret_cast<const matrix_t*>(&matrixA));
+	dc::mat_load_apply(matAB, matrixA);
 #endif

 	CColSphere bsphereAB;	// bounding sphere of A in B space
@@ -2099,8 +2098,7 @@ CCollision::ProcessColModels(const CMatrix &matrixA, CColModel &modelA,
 #ifndef DC_SH4
 	matBA *= matrixB;
 #else
-	mat_load(reinterpret_cast<const matrix_t*>(&matBA));
-	mat_apply(reinterpret_cast<const matrix_t*>(&matrixB));
+	dc::mat_load_apply(matBA, matrixB);
 #endif
 	for(i = 0; i < modelB.numSpheres; i++){
 		s.radius = modelB.spheres[i].radius;
@@ -2162,7 +2160,7 @@ CCollision::ProcessColModels(const CMatrix &matrixA, CColModel &modelA,
 	}

 #ifdef DC_SH4
-	mat_load(reinterpret_cast<matrix_t *>(const_cast<CMatrix *>(&matrixB)));
+	dc::mat_load2(matrixB);
 #endif
 	for(i = 0; i < numCollisions; i++){
 #ifndef DC_SH4
--- a/src/miami/core/ControllerConfig.cpp
+++ b/src/miami/core/ControllerConfig.cpp
--- a/src/miami/core/common.h
+++ b/src/miami/core/common.h
@@ -85,8 +85,6 @@

 #define rwVENDORID_ROCKSTAR 0x0253F2

-#define Max(a,b) ((a) > (b) ? (a) : (b))
-#define Min(a,b) ((a) < (b) ? (a) : (b))

 // Use this to add const that wasn't there in the original code
 #define Const const
@@ -303,12 +301,8 @@ extern int strncasecmp(const char *str1, const char *str2, size_t len);

 extern wchar *AllocUnicode(const char*src);

-#define Clamp(v, low, high) ((v)<(low) ? (low) : (v)>(high) ? (high) : (v))
-
-#define Clamp2(v, center, radius) ((v) > (center) ? Min(v, center + radius) : Max(v, center - radius))
-
-inline float sq(float x) { return x*x; }
 #define SQR(x) ((x) * (x))
+__always_inline auto sq(auto x) { return SQR(x); }

 #ifdef __MWERKS__
 #define M_E        2.71828182845904523536   // e
@@ -326,7 +320,10 @@ inline float sq(float x) { return x*x; }
 #define M_SQRT1_2  0.707106781186547524401  // 1/sqrt(2)
 #endif

-#define PI (float)M_PI
+#ifndef DC_SH4
+#define F_PI M_PI
+#endif
+#define PI (float)F_PI
 #define TWOPI (PI*2)
 #define HALFPI (PI/2)
 #define DEGTORAD(x) ((x) * PI / 180.0f)
@@ -380,7 +377,10 @@ __inline__ void TRACE(char *f, ...) { } // this is re3 only, and so the function
 #endif
 #endif

-#ifndef MASTER
+#ifdef assert
+#undef assert
+#endif
+#if !defined(MASTER)
 #define assert(_Expression) (void)( (!!(_Expression)) || (re3_assert(#_Expression, __FILE__, __LINE__, __FUNCTION__), 0) )
 #else
 #define assert(_Expression) (_Expression)
@@ -411,12 +411,7 @@ template<int s, int t> struct check_size {
 #endif
 #define BIT(num)                         (1<<(num))

-#define ABS(a)  (((a) < 0) ? (-(a)) : (a))
-#define norm(value, min, max) (((value) < (min)) ? 0 : (((value) > (max)) ? 1 : (((value) - (min)) / ((max) - (min)))))
-#define Lerp(norm, min, max) ( (norm) * ((max) - (min)) + (min) )
-
-#define STRINGIFY(x)                    #x
-#define STR(x)                          STRINGIFY(x)
-#define CONCAT_(x,y) x##y
-#define CONCAT(x,y) CONCAT_(x,y)
+#define ABS(a) Abs(a)

+// we use std::lerp now
+//#define lerp(norm, min, max) ( (norm) * ((max) - (min)) + (min) )
--- a/src/miami/math/Matrix.cpp
+++ b/src/miami/math/Matrix.cpp
@@ -4,7 +4,7 @@ CMatrix::CMatrix(CMatrix const &m)
 {
 	m_attachment = nil;
 	m_hasRwMatrix = false;
-	*this = m;
+	mat_copy(*this, m);
 }

 CMatrix::CMatrix(RwMatrix *matrix, bool owner)
@@ -75,7 +75,7 @@ CMatrix::UpdateRW(void)
 void
 CMatrix::operator=(CMatrix const &rhs)
 {
-	memcpy(this, &rhs, sizeof(f));
+	mat_copy(*this, rhs);
 	if (m_attachment)
 		UpdateRW();
 }
@@ -83,7 +83,7 @@ CMatrix::operator=(CMatrix const &rhs)
 void
 CMatrix::CopyOnlyMatrix(const CMatrix &other)
 {
-	memcpy(this, &other, sizeof(f));
+	mat_copy(*this, other);
 }

 CMatrix &
@@ -358,12 +358,14 @@ CMatrix::RotateZ(float z)
 void
 CMatrix::Rotate(float x, float y, float z)
 {
-	float cX = Cos(x);
-	float sX = Sin(x);
-	float cY = Cos(y);
-	float sY = Sin(y);
-	float cZ = Cos(z);
-	float sZ = Sin(z);
+#if 0 && defined(DC_SH4) // this is bugged and does not yield correct results
+       dc::mat_load2(reinterpret_cast<matrix_t *>(this));
+       mat_rotate(x, y, z);
+       dc::mat_store2(reinterpret_cast<matrix_t *>(this));
+#else
+	auto [sX, cX] = SinCos(x);
+	auto [sY, cY] = SinCos(y);
+	auto [sZ, cZ] = SinCos(z);
 	
 	float rx = this->rx;
 	float ry = this->ry;
@@ -388,6 +390,20 @@ CMatrix::Rotate(float x, float y, float z)
 	float z2 = sZ * sY - (cZ * sX) * cY;
 	float z3 = cX * cY;

+	#if !defined(DC_TEXCONV) && !defined(DC_SIM)
+	this->rx = fipr(x1, y1, z1, 0,  rx, ry, rz, 0);
+	this->ry = fipr(x2, y2, z2, 0,  rx, ry, rz, 0);
+	this->rz = fipr(x3, y3, z3, 0,  rx, ry, rz, 0);
+	this->fx = fipr(x1, y1, z1, 0,  ux, uy, uz, 0);
+	this->fy = fipr(x2, y2, z2, 0,  ux, uy, uz, 0);
+	this->fz = fipr(x3, y3, z3, 0,  ux, uy, uz, 0);
+	this->ux = fipr(x1, y1, z1, 0,  ax, ay, az, 0);
+	this->uy = fipr(x2, y2, z2, 0,  ax, ay, az, 0);
+	this->uz = fipr(x3, y3, z3, 0,  ax, ay, az, 0);
+	this->px = fipr(x1, y1, z1, 0,  px, py, pz, 0);
+	this->py = fipr(x2, y2, z2, 0,  px, py, pz, 0);
+	this->pz = fipr(x3, y3, z3, 0,  px, py, pz, 0);
+	#else
 	this->rx = x1 * rx + y1 * ry + z1 * rz;
 	this->ry = x2 * rx + y2 * ry + z2 * rz;
 	this->rz = x3 * rx + y3 * ry + z3 * rz;
@@ -400,6 +416,8 @@ CMatrix::Rotate(float x, float y, float z)
 	this->px = x1 * px + y1 * py + z1 * pz;
 	this->py = x2 * px + y2 * py + z2 * pz;
 	this->pz = x3 * px + y3 * py + z3 * pz;
+	#endif
+#endif
 }

 CMatrix &
@@ -429,9 +447,7 @@ operator*(const CMatrix &m1, const CMatrix &m2)
 	// TODO: VU0 code
 	CMatrix out;
 #if defined(RW_DC)
-	mat_load(reinterpret_cast<const matrix_t *>(&m1));
-	mat_apply(reinterpret_cast<const matrix_t *>(&m2));
-	mat_store(reinterpret_cast<matrix_t *>(&out));
+	mat_mult(out, m1, m2);
 #else
 	out.rx = m1.rx * m2.rx + m1.fx * m2.ry + m1.ux * m2.rz;
 	out.ry = m1.ry * m2.rx + m1.fy * m2.ry + m1.uy * m2.rz;
--- a/src/miami/math/Matrix.h
+++ b/src/miami/math/Matrix.h
@@ -43,6 +43,8 @@ public:
 		SetScale(scale);
 	}
 	~CMatrix(void);
+	operator matrix_t *() { return reinterpret_cast<matrix_t *>(this); }
+	operator const matrix_t *() const { return reinterpret_cast<const matrix_t *>(this); }
 	void Attach(RwMatrix *matrix, bool owner = false);
 	void AttachRW(RwMatrix *matrix, bool owner = false);
 	void Detach(void);
@@ -82,13 +84,12 @@ public:
 	void Scale(float sx, float sy, float sz)
 	{
 		for (int i = 0; i < 3; i++){
-			f[i][0] *= sx;
-			f[i][1] *= sy;
-			f[i][2] *= sz;
+				f[i][0] *= sx;
+				f[i][1] *= sy;
+				f[i][2] *= sz;
 		}
 	}

-
 	void SetRotateXOnly(float angle);
 	void SetRotateYOnly(float angle);
 	void SetRotateZOnly(float angle);
@@ -125,11 +126,22 @@ CMatrix Invert(const CMatrix &matrix);
 CMatrix operator*(const CMatrix &m1, const CMatrix &m2);
 inline CVector MultiplyInverse(const CMatrix &mat, const CVector &vec)
 {
+#ifndef DC_SH4
 	CVector v(vec.x - mat.px, vec.y - mat.py, vec.z - mat.pz);
 	return CVector(
 		mat.rx * v.x + mat.ry * v.y + mat.rz * v.z,
 		mat.fx * v.x + mat.fy * v.y + mat.fz * v.z,
 		mat.ux * v.x + mat.uy * v.y + mat.uz * v.z);
+#else
+    register float x asm(KOS_FPARG(0)) = vec.x - mat.px;
+    register float y asm(KOS_FPARG(1)) = vec.y - mat.py;
+    register float z asm(KOS_FPARG(2)) = vec.z - mat.pz;
+	return CVector(
+		fipr(x, y, z, 0.0f, mat.rx, mat.ry, mat.rz, 0.0f),
+		fipr(x, y, z, 0.0f, mat.fx, mat.fy, mat.fz, 0.0f),
+		fipr(x, y, z, 0.0f, mat.ux, mat.uy, mat.uz, 0.0f)
+	);
+#endif
 }


--- a/src/miami/math/Quaternion.cpp
+++ b/src/miami/math/Quaternion.cpp
@@ -39,10 +39,14 @@ CQuaternion::Slerp(const CQuaternion &q1, const CQuaternion &q2, float theta, fl
 void
 CQuaternion::Multiply(const CQuaternion &q1, const CQuaternion &q2)
 {
+#ifndef DC_SH4
 	x = (q2.z * q1.y) - (q1.z * q2.y) + (q1.x * q2.w) + (q2.x * q1.w);
 	y = (q2.x * q1.z) - (q1.x * q2.z) + (q1.y * q2.w) + (q2.y * q1.w);
 	z = (q2.y * q1.x) - (q1.y * q2.x) + (q1.z * q2.w) + (q2.z * q1.w);
 	w = (q2.w * q1.w) - (q2.x * q1.x) - (q2.y * q1.y) - (q2.z * q1.z);
+#else
+	quat_mult(*this, q1, q2);
+#endif
 }

 void
@@ -51,9 +55,16 @@ CQuaternion::Get(RwV3d *axis, float *angle)
 	*angle = Acos(w);
 	float s = Sin(*angle);

+#ifndef DC_SH4
 	axis->x = x * (1.0f / s);
 	axis->y = y * (1.0f / s);
 	axis->z = z * (1.0f / s);
+#else
+	float invS = dc::Invert<true, false>(s);
+	axis->x = x * invS;
+	axis->y = y * invS;
+	axis->z = z * invS;
+#endif
 }

 void
@@ -104,7 +115,7 @@ CQuaternion::Set(const RwMatrix &matrix)
 	if (f >= 0.0f) {
 		s = Sqrt(f + 1.0f);
 		w = 0.5f * s;
-		m = 0.5f / s;
+		m = Div<true, false>(0.5f, s);
 		x = (matrix.up.z - matrix.at.y) * m;
 		y = (matrix.at.x - matrix.right.z) * m;
 		z = (matrix.right.y - matrix.up.x) * m;
@@ -115,7 +126,7 @@ CQuaternion::Set(const RwMatrix &matrix)
 	if (f >= 0.0f) {
 		s = Sqrt(f + 1.0f);
 		x = 0.5f * s;
-		m = 0.5f / s;
+		m = Div<true, false>(0.5f, s);
 		y = (matrix.up.x + matrix.right.y) * m;
 		z = (matrix.at.x + matrix.right.z) * m;
 		w = (matrix.up.z - matrix.at.y) * m;
@@ -126,7 +137,7 @@ CQuaternion::Set(const RwMatrix &matrix)
 	if (f >= 0.0f) {
 		s = Sqrt(f + 1.0f);
 		y = 0.5f * s;
-		m = 0.5f / s;
+		m = Div<true, false>(0.5f, s);
 		w = (matrix.at.x - matrix.right.z) * m;
 		x = (matrix.up.x - matrix.right.y) * m;
 		z = (matrix.at.y + matrix.up.z) * m;
--- a/src/miami/math/Quaternion.h
+++ b/src/miami/math/Quaternion.h
@@ -1,5 +1,8 @@
 #pragma once

+#include "src/common_defines.h"
+#include "rwdc_common.h"
+
 // TODO: actually implement this
 class CQuaternion
 {
@@ -8,8 +11,24 @@ public:
 	CQuaternion(void) {}
 	CQuaternion(float x, float y, float z, float w) : x(x), y(y), z(z), w(w) {}

-	float Magnitude(void) const { return Sqrt(x*x + y*y + z*z + w*w); }
-	float MagnitudeSqr(void) const { return x*x + y*y + z*z + w*w; }
+	operator quaternion_t *() { return reinterpret_cast<quaternion_t *>(this); }
+	operator const quaternion_t *() const { return reinterpret_cast<const quaternion_t *>(this); }
+	operator quaternion_t &() { return *reinterpret_cast<quaternion_t *>(this); }
+	operator const quaternion_t &() const { return *reinterpret_cast<const quaternion_t *>(this); }
+	float Magnitude(void) const {
+#ifndef DC_SH4
+		return Sqrt(x*x + y*y + z*z + w*w);
+#else
+		return Sqrt(fipr_magnitude_sqr(x, y, z, w));
+#endif
+	}
+	float MagnitudeSqr(void) const {
+#ifndef DC_SH4
+		return x*x + y*y + z*z + w*w;
+#else
+		return fipr_magnitude_sqr(x, y, z, w);
+#endif
+	}
 	void Normalise(void);
 	void Multiply(const CQuaternion &q1, const CQuaternion &q2);
 	void Invert(void){	// Conjugate would have been a better name
@@ -43,10 +62,11 @@ public:
 	}

 	const CQuaternion &operator/=(float right) {
-		x /= right;
-		y /= right;
-		z /= right;
-		w /= right;
+		right = dc::Invert<false>(right);
+		x *= right;
+		y *= right;
+		z *= right;
+		w *= right;
 		return *this;
 	}

@@ -66,7 +86,11 @@ public:
 inline float
 DotProduct(const CQuaternion &q1, const CQuaternion &q2)
 {
+#ifndef DC_SH4
 	return q1.x*q2.x + q1.y*q2.y + q1.z*q2.z + q1.w*q2.w;
+#else
+	return fipr(q1.x, q1.y, q1.z, q1.w, q2.x, q2.y, q2.z, q2.w);
+#endif
 }

 inline CQuaternion operator+(const CQuaternion &left, const CQuaternion &right)
@@ -91,5 +115,6 @@ inline CQuaternion operator*(float left, const CQuaternion &right)

 inline CQuaternion operator/(const CQuaternion &left, float right)
 {
-	return CQuaternion(left.x / right, left.y / right, left.z / right, left.w / right);
+	right = Invert<false>(right);
+	return CQuaternion(left.x * right, left.y * right, left.z * right, left.w * right);
 }
--- a/src/miami/math/Vector.cpp
+++ b/src/miami/math/Vector.cpp
@@ -3,6 +3,10 @@
 void
 CVector::Normalise(void)
 {
+#ifdef DC_SH4_BROKEN
+	// TODO: This needs to handle zero vectors here
+	vec3f_normalize(x, y, z);
+#else
 	float sq = MagnitudeSqr();
 	if (sq > 0.0f) {
 		float invsqrt = RecipSqrt(sq);
@@ -11,6 +15,7 @@ CVector::Normalise(void)
 		z *= invsqrt;
 	} else
 		x = 1.0f;
+#endif
 }

 CVector
@@ -22,25 +27,49 @@ CrossProduct(const CVector &v1, const CVector &v2)
 CVector
 Multiply3x3(const CMatrix &mat, const CVector &vec)
 {
+#ifndef DC_SH4
 	// TODO: VU0 code
 	return CVector(mat.rx * vec.x + mat.fx * vec.y + mat.ux * vec.z,
 	               mat.ry * vec.x + mat.fy * vec.y + mat.uy * vec.z,
 	               mat.rz * vec.x + mat.fz * vec.y + mat.uz * vec.z);
+#else
+	CVector out;
+	dc::mat_load2(mat);
+	mat_trans_normal3_nomod(vec.x, vec.y, vec.z,
+	                        out.x, out.y, out.z);
+    return out;
+#endif
 }

 CVector
 Multiply3x3(const CVector &vec, const CMatrix &mat)
 {
+#ifndef DC_SH4
 	return CVector(mat.rx * vec.x + mat.ry * vec.y + mat.rz * vec.z,
 	               mat.fx * vec.x + mat.fy * vec.y + mat.fz * vec.z,
 	               mat.ux * vec.x + mat.uy * vec.y + mat.uz * vec.z);
+#else
+	CVector out;
+	dc::mat_load2(mat);
+	mat_transpose();
+	mat_trans_normal3_nomod(vec.x, vec.y, vec.z,
+	                        out.x, out.y, out.z);
+	return out;
+#endif
 }

 CVector
 operator*(const CMatrix &mat, const CVector &vec)
 {
+#ifdef DC_SH4
+	CVector out;
+	dc::mat_load2(mat);
+	mat_trans_single3_nodiv_nomod(vec.x, vec.y, vec.z, out.x, out.y, out.z);
+	return out;
+#else
 	// TODO: VU0 code
 	return CVector(mat.rx * vec.x + mat.fx * vec.y + mat.ux * vec.z + mat.px,
 	               mat.ry * vec.x + mat.fy * vec.y + mat.uy * vec.z + mat.py,
 	               mat.rz * vec.x + mat.fz * vec.y + mat.uz * vec.z + mat.pz);
+#endif
 }
--- a/src/miami/math/Vector.h
+++ b/src/miami/math/Vector.h
@@ -1,5 +1,7 @@
 #pragma once

+#include "maths.h"
+
 class CVector : public RwV3d
 {
 public:
@@ -18,11 +20,25 @@ public:
 		z = v.z;
 	}
 	// (0,1,0) means no rotation. So get right vector and its atan
-	float Heading(void) const { return Atan2(-x, y); }
-	float Magnitude(void) const { return Sqrt(x*x + y*y + z*z); }
-	float MagnitudeSqr(void) const { return x*x + y*y + z*z; }
-	float Magnitude2D(void) const { return Sqrt(x*x + y*y); }
-	float MagnitudeSqr2D(void) const { return x*x + y*y; }
+__always_inline float Heading(void) const { return Atan2(-x, y); }
+  __always_inline float Magnitude(void) const {
+#ifdef DC_SH4
+    float w;
+    vec3f_length(x, y, z, w);
+    return w;
+#else
+    return Sqrt(x*x + y*y + z*z);
+#endif
+}
+  __always_inline float MagnitudeSqr(void) const {
+#ifdef DC_SH4
+    return fipr_magnitude_sqr(x, y,z, 0.0f);
+#else
+    return x*x + y*y + z*z;
+#endif
+}
+  __always_inline float Magnitude2D(void) const { return Sqrt(x*x + y*y); }
+  float MagnitudeSqr2D(void) const { return x*x + y*y; }
 	void Normalise(void);
 	
 	void Normalise2D(void) {
@@ -54,9 +70,16 @@ public:
 	}

 	const CVector &operator/=(float right) {
+#ifndef DC_SH4
 		x /= right;
 		y /= right;
 		z /= right;
+#else
+        right = dc::Invert<true, true>(right);
+        x *= right;
+        y *= right;
+        z *= right;
+#endif
 		return *this;
 	}

@@ -97,13 +120,22 @@ inline CVector operator*(float left, const CVector &right)

 inline CVector operator/(const CVector &left, float right)
 {
+#ifndef DC_SH4
 	return CVector(left.x / right, left.y / right, left.z / right);
+#else
+    right = dc::Invert<true, true>(right);
+    return CVector(left.x * right, left.y * right, left.z * right);
+#endif
 }

 inline float
 DotProduct(const CVector &v1, const CVector &v2)
 {
-	return v1.x*v2.x + v1.y*v2.y + v1.z*v2.z;
+#ifdef DC_SH4
+  return fipr(v1.x, v1.y, v1.z, 0.0f, v2.x, v2.y, v2.z, 0.0f);
+#else
+  return v1.x*v2.x + v1.y*v2.y + v1.z*v2.z;
+#endif
 }

 CVector CrossProduct(const CVector &v1, const CVector &v2);
@@ -111,7 +143,13 @@ CVector CrossProduct(const CVector &v1, const CVector &v2);
 inline float
 Distance(const CVector &v1, const CVector &v2)
 {
-	return (v2 - v1).Magnitude();
+  float w;
+#ifdef DC_SH4
+  vec3f_distance(v1.x, v1.y, v1.z, v2.x, v2.y, v2.z, w);
+  return w;
+#else
+  return (v2 - v1).Magnitude();
+#endif
 }

 inline float
--- a/src/miami/math/VuVector.h
+++ b/src/miami/math/VuVector.h
@@ -50,7 +50,7 @@ __always_inline void TransformPoint(CVuVector &out, const CMatrix &mat, const CV
 		sqc2    vf06,0x0(%0)\n\
 		": : "r" (&out) , "r" (&mat) ,"r" (&in): "memory");
 #elif defined(DC_SH4)
-	mat_load(reinterpret_cast<matrix_t *>(const_cast<CMatrix *>(&mat)));
+	dc::mat_load2(mat);
 	mat_trans_single3_nodiv_nomod(in.x, in.y, in.z, out.x, out.y, out.z);
 #else
 	out = mat * in;
@@ -77,7 +77,7 @@ __always_inline void TransformPoint(CVuVector &out, const CMatrix &mat, const Rw
 		sqc2    vf06,0x0(%0)\n\
 		": : "r" (&out) , "r" (&mat) ,"r" (&in): "memory");
 #elif defined(DC_SH4)
-	mat_load(reinterpret_cast<matrix_t *>(const_cast<CMatrix *>(&mat)));
+	dc::mat_load2(mat);
 	mat_trans_single3_nodiv_nomod(in.x, in.y, in.z, out.x, out.y, out.z);
 #else
 	out = mat * in;
@@ -86,7 +86,7 @@ __always_inline void TransformPoint(CVuVector &out, const CMatrix &mat, const Rw

 __always_inline void TransformPoints(CVuVector *out, int n, const CMatrix &mat, const RwV3d *in, int stride)
 {
-#ifdef GTA_PS3
+#ifdef GTA_PS2
 	__asm__ __volatile__("\n\
 		paddub	$3,%4,$0\n\
 		lqc2    vf02,0x0(%2)\n\
@@ -114,7 +114,7 @@ __always_inline void TransformPoints(CVuVector *out, int n, const CMatrix &mat,
 		bnez	%1,1b\n\
 		": : "r" (out) , "r" (n), "r" (&mat), "r" (in), "r" (stride): "memory");
 #elif defined(DC_SH4)
-    mat_load(reinterpret_cast<matrix_t *>(const_cast<CMatrix *>(&mat)));
+    dc::mat_load2(mat);
 	while(n--) {
 		mat_trans_single3_nodiv_nomod(in->x, in->y, in->z, out->x, out->y, out->z);
 		in = reinterpret_cast<const RwV3d *>(reinterpret_cast<const uint8_t *>(in) + stride);
@@ -137,7 +137,6 @@ __always_inline void TransformPoints(CVuVector *out, int n, const CMatrix &mat,
 		lqc2    vf03,0x10(%2)\n\
 		lqc2    vf04,0x20(%2)\n\
 		lqc2    vf05,0x30(%2)\n\
-		lqc2    vf01,0x0(%3)\n\
 		nop\n\
 	1:	vmulax.xyz	ACC,   vf02,vf01\n\
 		vmadday.xyz	ACC,   vf03,vf01\n\
--- a/src/miami/math/maths.h
+++ b/src/miami/math/maths.h
@@ -1,46 +1,6 @@
 #pragma once

 #include "src/common_defines.h"
+#include "rwdc_common.h"

-#include <dc/matrix.h>
-
-#ifdef DC_SH4
-
-#define mat_trans_nodiv_nomod(x, y, z, x2, y2, z2, w2) do { \
-        register float __x __asm__("fr12") = (x); \
-        register float __y __asm__("fr13") = (y); \
-        register float __z __asm__("fr14") = (z); \
-        register float __w __asm__("fr15") = 1.0f; \
-        __asm__ __volatile__( "ftrv  xmtrx, fv12\n" \
-                              : "=f" (__x), "=f" (__y), "=f" (__z), "=f" (__w) \
-                              : "0" (__x), "1" (__y), "2" (__z), "3" (__w) ); \
-        x2 = __x; y2 = __y; z2 = __z; w2 = __w; \
-    } while(false)
-
-#else
-
-#define mat_trans_nodiv_nomod(x_, y_, z_, x2, y2, z2, w2) do { \
-		vector_t tmp = { x_, y_, z_, 1.0f }; \
-		mat_transform(&tmp, &tmp, 1, 0); \
-		x2 = tmp.x; y2 = tmp.y; z2 = tmp.z; w2 = tmp.w; \
-	} while(false)
-#endif
-
-
-// wrapper around float versions of functions
-// in gta they are in CMaths but that makes the code rather noisy
-
-inline float Sin(float x) { return sinf(x); }
-inline float Asin(float x) { return asinf(x); }
-inline float Cos(float x) { return cosf(x); }
-inline float Acos(float x) { return acosf(x); }
-inline float Tan(float x) { return tanf(x); }
-inline float Atan(float x) { return atanf(x); }
-inline float Atan2(float y, float x) { return atan2f(y, x); }
-inline float Abs(float x) { return fabsf(x); }
-inline float Sqrt(float x) { return sqrtf(x); }
-inline float RecipSqrt(float x, float y) { return x/Sqrt(y); }
-inline float RecipSqrt(float x) { return RecipSqrt(1.0f, x); }
-inline float Pow(float x, float y) { return powf(x, y); }
-inline float Floor(float x) { return floorf(x); }
-inline float Ceil(float x) { return ceilf(x); }
+using namespace dc;
--- a/vendor/dca3-kos
+++ b/vendor/dca3-kos
--- a/vendor/librw/src/base.cpp
+++ b/vendor/librw/src/base.cpp
@@ -45,12 +45,17 @@ int32 build = 0xFFFF;
 bool32 streamAppendFrames = 0;
 char *debugFile = nil;

-static Matrix identMat = {
-	{ 1.0f, 0.0f, 0.0f }, Matrix::IDENTITY|Matrix::TYPEORTHONORMAL,
-	{ 0.0f, 1.0f, 0.0f }, 0,
-	{ 0.0f, 0.0f, 1.0f }, 0,
-	{ 0.0f, 0.0f, 0.0f }, 0
-};
+static Matrix identMat = {{
+    .right 	= { 1.0f, 0.0f, 0.0f }, 
+    .flags 	= Matrix::IDENTITY|Matrix::TYPEORTHONORMAL,
+	.pad0   = 0,
+    .up 	= { 0.0f, 1.0f, 0.0f }, 
+    .upw 	= 0.0f,
+    .at 	= { 0.0f, 0.0f, 1.0f }, 
+    .atw 	= 0.0f,
+    .pos 	= { 0.0f, 0.0f, 0.0f }, 
+    .posw 	= 1.0f
+}};

 // lazy implementation
 int
@@ -86,16 +91,6 @@ strncmp_ci(const char *s1, const char *s2, int n)
 	return 0;
 }

-Quat
-mult(const Quat &q, const Quat &p)
-{
-	return makeQuat(q.w*p.w - q.x*p.x - q.y*p.y - q.z*p.z,
-	                q.w*p.x + q.x*p.w + q.y*p.z - q.z*p.y,
-	                q.w*p.y + q.y*p.w + q.z*p.x - q.x*p.z,
-	                q.w*p.z + q.z*p.w + q.x*p.y - q.y*p.x);
-}
-
-
 Quat*
 Quat::rotate(const V3d *axis, float32 angle, CombineOp op)
 {
@@ -143,8 +138,9 @@ slerp(const Quat &q, const Quat &p, float32 a)
 	float32 phi = acosf(c);
 	if(phi > 0.00001f){
 		float32 s = sinf(phi);
-		return add(scale(q1, sinf((1.0f-a)*phi)/s),
-		           scale(p,  sinf(a*phi)/s));
+        float invS = dc::Invert<true, false>(s);
+		return add(scale(q1, sinf((1.0f-a)*phi) * invS),
+		           scale(p,  sinf(a*phi) * invS));
 	}
 	return q1;
 }
@@ -152,39 +148,39 @@ slerp(const Quat &q, const Quat &p, float32 a)
 //
 // V3d
 //
-
-V3d
-cross(const V3d &a, const V3d &b)
-{
-	return makeV3d(a.y*b.z - a.z*b.y,
-	               a.z*b.x - a.x*b.z,
-	               a.x*b.y - a.y*b.x);
+void V3d::transformPoints(V3d *out, const V3d *in, int32 n, const Matrix *m) {
+    int32 i;
+    #ifndef DC_SH4
+        V3d tmp;
+        for(i = 0; i < n; i++){
+            tmp.x = in[i].x*m->right.x + in[i].y*m->up.x + in[i].z*m->at.x + m->pos.x;
+            tmp.y = in[i].x*m->right.y + in[i].y*m->up.y + in[i].z*m->at.y + m->pos.y;
+            tmp.z = in[i].x*m->right.z + in[i].y*m->up.z + in[i].z*m->at.z + m->pos.z;
+            out[i] = tmp;
+        }
+    #else
+        dc::mat_load2(*m);
+        for(i = 0; i < n; i++)
+            mat_trans_single3_nodiv_nomod(in[i].x, in[i].y, in[i].z,
+                                          out[i].x, out[i].y, out[i].z);
+    #endif
 }
-
-void
-V3d::transformPoints(V3d *out, const V3d *in, int32 n, const Matrix *m)
-{
-	int32 i;
-	V3d tmp;
-	for(i = 0; i < n; i++){
-		tmp.x = in[i].x*m->right.x + in[i].y*m->up.x + in[i].z*m->at.x + m->pos.x;
-		tmp.y = in[i].x*m->right.y + in[i].y*m->up.y + in[i].z*m->at.y + m->pos.y;
-		tmp.z = in[i].x*m->right.z + in[i].y*m->up.z + in[i].z*m->at.z + m->pos.z;
-		out[i] = tmp;
-	}
-}
-
-void
-V3d::transformVectors(V3d *out, const V3d *in, int32 n, const Matrix *m)
-{
-	int32 i;
-	V3d tmp;
-	for(i = 0; i < n; i++){
-		tmp.x = in[i].x*m->right.x + in[i].y*m->up.x + in[i].z*m->at.x;
-		tmp.y = in[i].x*m->right.y + in[i].y*m->up.y + in[i].z*m->at.y;
-		tmp.z = in[i].x*m->right.z + in[i].y*m->up.z + in[i].z*m->at.z;
-		out[i] = tmp;
-	}
+void V3d::transformVectors(V3d *out, const V3d *in, int32 n, const Matrix *m) {
+    int32 i;
+    #ifndef DC_SH4
+        V3d tmp;
+        for(i = 0; i < n; i++){
+            tmp.x = in[i].x*m->right.x + in[i].y*m->up.x + in[i].z*m->at.x;
+            tmp.y = in[i].x*m->right.y + in[i].y*m->up.y + in[i].z*m->at.y;
+            tmp.z = in[i].x*m->right.z + in[i].y*m->up.z + in[i].z*m->at.z;
+            out[i] = tmp;
+        }
+    #else
+        dc::mat_load2(*m);
+        for(i = 0; i < n; i++)
+            mat_trans_normal3_nomod(in[i].x, in[i].y, in[i].z,
+                                    out[i].x, out[i].y, out[i].z);
+    #endif
 }

 //
@@ -194,6 +190,7 @@ V3d::transformVectors(V3d *out, const V3d *in, int32 n, const Matrix *m)
 void
 RawMatrix::mult(RawMatrix *dst, RawMatrix *src1, RawMatrix *src2)
 {
+#ifndef DC_SH4
 	dst->right.x = src1->right.x*src2->right.x + src1->right.y*src2->up.x + src1->right.z*src2->at.x + src1->rightw*src2->pos.x;
 	dst->right.y = src1->right.x*src2->right.y + src1->right.y*src2->up.y + src1->right.z*src2->at.y + src1->rightw*src2->pos.y;
 	dst->right.z = src1->right.x*src2->right.z + src1->right.y*src2->up.z + src1->right.z*src2->at.z + src1->rightw*src2->pos.z;
@@ -210,11 +207,15 @@ RawMatrix::mult(RawMatrix *dst, RawMatrix *src1, RawMatrix *src2)
 	dst->pos.y   = src1->pos.x*src2->right.y   + src1->pos.y*src2->up.y   + src1->pos.z*src2->at.y + src1->posw*src2->pos.y;
 	dst->pos.z   = src1->pos.x*src2->right.z   + src1->pos.y*src2->up.z   + src1->pos.z*src2->at.z + src1->posw*src2->pos.z;
 	dst->posw    = src1->pos.x*src2->rightw    + src1->pos.y*src2->upw    + src1->pos.z*src2->atw  + src1->posw*src2->posw;
+#else
+    dc::mat_mult(*dst, *src2, *src1);
+#endif
 }

 void
 RawMatrix::transpose(RawMatrix *dst, RawMatrix *src)
 {
+#ifndef DC_SH4
 	dst->right.x = src->right.x;
 	dst->up.x = src->right.y;
 	dst->at.x = src->right.z;
@@ -231,18 +232,27 @@ RawMatrix::transpose(RawMatrix *dst, RawMatrix *src)
 	dst->upw = src->pos.y;
 	dst->atw = src->pos.z;
 	dst->posw = src->posw;
+#else
+    dc::mat_load_transpose(*src);
+    dc::mat_store2(*dst);
+#endif
 }

 void
 RawMatrix::setIdentity(RawMatrix *dst)
 {
-	static RawMatrix identity = {
+#ifndef DC_SH4
+	static RawMatrix identity = {{
 		{ 1.0f, 0.0f, 0.0f }, 0.0f,
 		{ 0.0f, 1.0f, 0.0f }, 0.0f,
 		{ 0.0f, 0.0f, 1.0f }, 0.0f,
 		{ 0.0f, 0.0f, 0.0f }, 1.0f
-	};
+	}};
 	*dst = identity;
+#else
+    dc::mat_identity2();
+    dc::mat_store2(*dst);
+#endif
 }

 //
@@ -301,9 +311,10 @@ Matrix::mult(Matrix *dst, const Matrix *src1, const Matrix *src2)
 		*dst = *src2;
 	else if(src2->flags & IDENTITY)
 		*dst = *src1;
-	else{
+	else {
+        uint8_t flags = src1->flags & src2->flags;
 		mult_(dst, src1, src2);
-		dst->flags = src1->flags & src2->flags;
+		dst->flags = flags;
 	}
 	return dst;
 }
@@ -324,7 +335,8 @@ Matrix::invert(Matrix *dst, const Matrix *src)
 Matrix*
 Matrix::transpose(Matrix *dst, const Matrix *src)
 {
-	if(src->flags & IDENTITY)
+#ifndef DC_SH4
+	if(src->flags & IDENTITY) 
 		*dst = *src;
 	dst->right.x = src->right.x;
 	dst->up.x = src->right.y;
@@ -338,25 +350,31 @@ Matrix::transpose(Matrix *dst, const Matrix *src)
 	dst->pos.x = 0.0;
 	dst->pos.y = 0.0;
 	dst->pos.z = 0.0;
+#else
+    if(src->flags & IDENTITY)
+        *dst = *src;
+    else {
+        dc::mat_load_transpose(*src);
+        dc::mat_store2(*dst);
+    }
+#endif
 	return dst;
 }

 Matrix*
 Matrix::rotate(const V3d *axis, float32 angle, CombineOp op)
 {
-	Matrix tmp, rot;
-	makeRotation(&rot, axis, angle);
+	Matrix rot;
+    makeRotation(&rot, axis, angle);
 	switch(op){
 	case COMBINEREPLACE:
 		*this = rot;
 		break;
 	case COMBINEPRECONCAT:
-		mult(&tmp, &rot, this);
-		*this = tmp;
+		mult(this, &rot, this);
 		break;
 	case COMBINEPOSTCONCAT:
-		mult(&tmp, this, &rot);
-		*this = tmp;
+        mult(this, this, &rot);
 		break;
 	}
 	return this;
@@ -365,41 +383,44 @@ Matrix::rotate(const V3d *axis, float32 angle, CombineOp op)
 Matrix*
 Matrix::rotate(const Quat &q, CombineOp op)
 {
-	Matrix tmp, rot;
-	makeRotation(&rot, q);
+	Matrix rot;
+    makeRotation(&rot, q);
 	switch(op){
 	case COMBINEREPLACE:
 		*this = rot;
 		break;
 	case COMBINEPRECONCAT:
-		mult(&tmp, &rot, this);
-		*this = tmp;
+        mult(this, &rot, this);
 		break;
 	case COMBINEPOSTCONCAT:
-		mult(&tmp, this, &rot);
-		*this = tmp;
+        mult(this, this, &rot);
 		break;
 	}
 	return this;
 }
+
 Matrix*
 Matrix::translate(const V3d *translation, CombineOp op)
 {
-	Matrix tmp;
+#if 1
 	Matrix trans = identMat;
 	trans.pos = *translation;
 	trans.flags &= ~IDENTITY;
+#else
+    Matrix trans;
+    dc::mat_set_translation(translation->x, translation->y, translation->z);
+    dc::mat_store2(trans);
+    trans.flags = TYPEORTHONORMAL;
+#endif
 	switch(op){
 	case COMBINEREPLACE:
 		*this = trans;
 		break;
 	case COMBINEPRECONCAT:
-		mult(&tmp, &trans, this);
-		*this = tmp;
+		mult(this, &trans, this);
 		break;
 	case COMBINEPOSTCONCAT:
-		mult(&tmp, this, &trans);
-		*this = tmp;
+		mult(this, this, &trans);
 		break;
 	}
 	return this;
@@ -408,23 +429,27 @@ Matrix::translate(const V3d *translation, CombineOp op)
 Matrix*
 Matrix::scale(const V3d *scale, CombineOp op)
 {
-	Matrix tmp;
+#ifndef DC_SH4
 	Matrix scl = identMat;
 	scl.right.x = scale->x;
 	scl.up.y = scale->y;
 	scl.at.z = scale->z;
 	scl.flags &= ~IDENTITY;
+#else
+    Matrix scl;
+    dc::mat_set_scale(scale->x, scale->y, scale->z);
+    dc::mat_store2(scl);
+    scl.flags = TYPEORTHONORMAL;
+#endif
 	switch(op){
 	case COMBINEREPLACE:
 		*this = scl;
 		break;
 	case COMBINEPRECONCAT:
-		mult(&tmp, &scl, this);
-		*this = tmp;
+		mult(this, &scl, this);
 		break;
 	case COMBINEPOSTCONCAT:
-		mult(&tmp, this, &scl);
-		*this = tmp;
+		mult(this, this, &scl);
 		break;
 	}
 	return this;
@@ -433,18 +458,15 @@ Matrix::scale(const V3d *scale, CombineOp op)
 Matrix*
 Matrix::transform(const Matrix *mat, CombineOp op)
 {
-	Matrix tmp;
 	switch(op){
 	case COMBINEREPLACE:
 		*this = *mat;
 		break;
 	case COMBINEPRECONCAT:
-		mult(&tmp, mat, this);
-		*this = tmp;
+		mult(this, mat, this);
 		break;
 	case COMBINEPOSTCONCAT:
-		mult(&tmp, this, mat);
-		*this = tmp;
+		mult(this, this, mat);
 		break;
 	}
 	return this;
@@ -459,27 +481,31 @@ Matrix::getRotation(void)
 	if(tr > 0.0f){
 		s = sqrtf(1.0f + tr) * 2.0f;
 		q.w = s / 4.0f;
-		q.x = (up.z - at.y) / s;
-		q.y = (at.x - right.z) / s;
-		q.z = (right.y - up.x) / s;
+        float invS = dc::Invert<true, false>(s);
+		q.x = (up.z - at.y) * invS;
+		q.y = (at.x - right.z) * invS;
+		q.z = (right.y - up.x) * invS;
 	}else if(right.x > up.y && right.x > at.z){
 		s = sqrtf(1.0f + right.x - up.y - at.z) * 2.0f;
-		q.w = (up.z - at.y) / s;
-		q.x = s / 4.0f;
-		q.y = (up.x + right.y) / s;
-		q.z = (at.x + right.z) / s;
+        q.x = s / 4.0f;
+        float invS = dc::Invert<true, false>(s);
+        q.w = (up.z - at.y) * invS;
+		q.y = (up.x + right.y) * invS;
+		q.z = (at.x + right.z) * invS;
 	}else if(up.y > at.z){
 		s = sqrtf(1.0f + up.y - right.x - at.z) * 2.0f;
-		q.w = (at.x - right.z) / s;
-		q.x = (up.x + right.y) / s;
-		q.y = s / 4.0f;
-		q.z = (at.y + up.z) / s;
+        q.y = s / 4.0f;
+        float invS = dc::Invert<true, false>(s);
+        q.w = (at.x - right.z) * invS;
+		q.x = (up.x + right.y) * invS;
+		q.z = (at.y + up.z) * invS;
 	}else{
 		s = sqrtf(1.0f + at.z - right.x - up.y) * 2.0f;
-		q.w = (right.y - up.x) / s;
-		q.x = (at.x + right.z) / s;
-		q.y = (at.y + up.z) / s;
-		q.z = s / 4.0f;
+        q.z = s / 4.0f;
+        float invS = dc::Invert<true, false>(s);
+        q.w = (right.y - up.x) * invS;
+		q.x = (at.x + right.z) * invS;
+		q.y = (at.y + up.z) * invS;
 	}
 	return q;
 }
@@ -501,20 +527,7 @@ Matrix::lookAt(const V3d &dir, const V3d &up)
 void
 Matrix::mult_(Matrix *__restrict__ dst, const Matrix *__restrict__ src1, const Matrix *__restrict__ src2)
 {
-	#if !defined(DC_TEXCONV) && !defined(DC_SIM)
-	dst->right.x = fipr(src1->right.x, src1->right.y,  src1->right.z, 0, 		  src2->right.x, src2->up.x, src2->at.x, 0);
-	dst->right.y = fipr(src1->right.x, src1->right.y,  src1->right.z, 0, 		  src2->right.y, src2->up.y, src2->at.y, 0);
-	dst->right.z = fipr(src1->right.x, src1->right.y,  src1->right.z, 0, 		  src2->right.z, src2->up.z, src2->at.z, 0);
-	dst->up.x    = fipr(src1->up.x,    src1->up.y,  src1->up.z, 0, 				  src2->right.x, src2->up.x, src2->at.x, 0);
-	dst->up.y    = fipr(src1->up.x,    src1->up.y,  src1->up.z, 0, 				  src2->right.y, src2->up.y, src2->at.y, 0);
-	dst->up.z    = fipr(src1->up.x,    src1->up.y,  src1->up.z, 0, 				  src2->right.z, src2->up.z, src2->at.z, 0);
-	dst->at.x    = fipr(src1->at.x,    src1->at.y,  src1->at.z, 0, 				  src2->right.x, src2->up.x, src2->at.x, 0);
-	dst->at.y    = fipr(src1->at.x,    src1->at.y,  src1->at.z, 0, 				  src2->right.y, src2->up.y, src2->at.y, 0);
-	dst->at.z    = fipr(src1->at.x,    src1->at.y,  src1->at.z, 0, 				  src2->right.z, src2->up.z, src2->at.z, 0);
-	dst->pos.x   = fipr(src1->pos.x,   src1->pos.y,  src1->pos.z, 1, 	  		  src2->right.x, src2->up.x, src2->at.x, src2->pos.x);
-	dst->pos.y   = fipr(src1->pos.x,   src1->pos.y,  src1->pos.z, 1, 	  	 	  src2->right.y, src2->up.y, src2->at.y, src2->pos.y);
-	dst->pos.z   = fipr(src1->pos.x,   src1->pos.y,  src1->pos.z, 1, 	  		  src2->right.z, src2->up.z, src2->at.z, src2->pos.z);
-	#else
+#ifndef DC_SH4
 	dst->right.x = src1->right.x*src2->right.x + src1->right.y*src2->up.x + src1->right.z*src2->at.x;
 	dst->right.y = src1->right.x*src2->right.y + src1->right.y*src2->up.y + src1->right.z*src2->at.y;
 	dst->right.z = src1->right.x*src2->right.z + src1->right.y*src2->up.z + src1->right.z*src2->at.z;
@@ -527,12 +540,15 @@ Matrix::mult_(Matrix *__restrict__ dst, const Matrix *__restrict__ src1, const M
 	dst->pos.x   = src1->pos.x*src2->right.x   + src1->pos.y*src2->up.x   + src1->pos.z*src2->at.x + src2->pos.x;
 	dst->pos.y   = src1->pos.x*src2->right.y   + src1->pos.y*src2->up.y   + src1->pos.z*src2->at.y + src2->pos.y;
 	dst->pos.z   = src1->pos.x*src2->right.z   + src1->pos.y*src2->up.z   + src1->pos.z*src2->at.z + src2->pos.z;
-	#endif
+#else
+    dc::mat_mult(*dst, *src2, *src1);
+#endif
 }

 void
 Matrix::invertOrthonormal(Matrix *dst, const Matrix *src)
 {
+#if 1
 	dst->right.x = src->right.x;
 	dst->right.y = src->up.x;
 	dst->right.z = src->at.x;
@@ -551,7 +567,12 @@ Matrix::invertOrthonormal(Matrix *dst, const Matrix *src)
 	dst->pos.z = -(src->pos.x*src->at.x +
 	               src->pos.y*src->at.y +
 	               src->pos.z*src->at.z);
-	dst->flags = TYPEORTHONORMAL;
+#else
+    dc::mat_load_transpose(*src);
+    dc::mat_invert_tranpose();
+    dc::mat_store2(*dst);
+#endif
+    dst->flags = TYPEORTHONORMAL;
 }

 Matrix*
@@ -587,8 +608,13 @@ void
 Matrix::makeRotation(Matrix *dst, const V3d *axis, float32 angle)
 {
 //	V3d v = normalize(*axis);
+#ifndef DC_SH4
 	float32 len = dot(*axis, *axis);
 	if(len != 0.0f) len = 1.0f/sqrtf(len);
+#else
+    float len = fipr_magnitude_sqr(axis->x, axis->y, axis->z, 0.0f);
+    if(len != 0.0f) len = dc::RecipSqrt(len);
+#endif
 	V3d v = rw::scale(*axis, len);
 	angle = angle*(float)M_PI/180.0f;
 	float32 s = sinf(angle);
@@ -646,7 +672,11 @@ Matrix::normalError(void)
 	x = dot(right, right) - 1.0f;
 	y = dot(up, up) - 1.0f;
 	z = dot(at, at) - 1.0f;
+#ifndef DC_SH4
 	return x*x + y*y + z*z;
+#else
+    return fipr_magnitude_sqr(x, y, z, 0.0f);
+#endif
 }

 float32
@@ -656,16 +686,27 @@ Matrix::orthogonalError(void)
 	x = dot(at, up);
 	y = dot(at, right);
 	z = dot(up, right);
+#ifndef DC_SH4
 	return x*x + y*y + z*z;
+#else
+    return fipr_magnitude_sqr(x, y, z, 0.0f);
+#endif
 }

 float32
 Matrix::identityError(void)
 {
-	V3d r = { right.x-1.0f, right.y, right.z };
+    V3d r = { right.x-1.0f, right.y, right.z };
 	V3d u = { up.x, up.y-1.0f, up.z };
 	V3d a = { at.x, at.y, at.z-1.0f };
+#ifndef DC_SH4
 	return dot(r,r) + dot(u,u) + dot(a,a) + dot(pos,pos);
+#else
+    return fipr_magnitude_sqr(r.x, r.y, r.z, 0.0f)    +
+           fipr_magnitude_sqr(u.x, u.y, u.z, 0.0f)    +
+           fipr_magnitude_sqr(at.x, at.y, at.z, 0.0f) +
+           fipr_magnitude_sqr(pos.x, pos.y, pos.z, 0.0f);
+#endif
 }

 void
--- a/vendor/librw/src/camera.cpp
+++ b/vendor/librw/src/camera.cpp
@@ -429,6 +429,7 @@ Camera::frustumTestSphere(const Sphere *s) const
 {
 	int32 res = SPHEREINSIDE;
 	const FrustumPlane *p = this->frustumPlanes;
+#ifndef DC_SH4
 	for(int32 i = 0; i < 6; i++){
 		float32 dist = dot(p->plane.normal, s->center) - p->plane.distance;
 		if(s->radius < dist)
@@ -437,6 +438,88 @@ Camera::frustumTestSphere(const Sphere *s) const
 			res = SPHEREBOUNDARY;
 		p++;
 	}
+#else
+	__builtin_prefetch(p);
+
+	register float sx asm("fr0") = s->center.x;
+	register float sy asm("fr1") = s->center.y;
+	register float sz asm("fr2") = s->center.z;
+	register float sw asm("fr3") = -1.0f;
+
+	// far
+	register float px asm("fr4") = p->plane.normal.x;
+	register float py asm("fr5") = p->plane.normal.y;
+	register float pz asm("fr6") = p->plane.normal.z;
+	register float pw asm("fr7") = p->plane.distance;
+
+	asm volatile("fipr fv0, fv4"
+				: "+f" (pw)
+				: "f" (sx), "f" (sy), "f" (sz), "f" (sw),
+				"f" (px), "f" (py), "f" (pz));
+	if(s->radius < pw)
+		return SPHEREOUTSIDE;
+	else if(s->radius > -pw)
+		res = SPHEREBOUNDARY;
+	p++;
+
+	// near
+	px = p->plane.normal.x;
+	py = p->plane.normal.y;
+	pz = p->plane.normal.z;
+	pw = p->plane.distance;
+	asm volatile("fipr fv0, fv4"
+				: "+f" (pw)
+				: "f" (sx), "f" (sy), "f" (sz), "f" (sw),
+				"f" (px), "f" (py), "f" (pz));
+	if(s->radius < pw)
+		return SPHEREOUTSIDE;
+	if(s->radius > -pw)
+		res = SPHEREBOUNDARY_NEAR;
+	p++;
+
+	const float* base_ptr0 = &p[0].plane.normal.x;
+	const float* base_ptr1 = &p[1].plane.normal.x;
+	const float* base_ptr2 = &p[2].plane.normal.x;
+	const float* base_ptr3 = &p[3].plane.normal.x;
+
+	__builtin_prefetch(base_ptr0);
+
+	static_assert(offsetof (decltype (p[0].plane.normal), y)
+				 -offsetof (decltype (p[0].plane.normal), x) == sizeof (float));
+
+	static_assert(offsetof (decltype (p[0].plane.normal), z)
+				  -offsetof (decltype (p[0].plane.normal), y) == sizeof (float));
+
+	static_assert(offsetof (decltype (p[0].plane), distance)
+				 -offsetof (decltype (p[0].plane.normal), z) == sizeof (float));
+
+    mat_load_rows(base_ptr0, base_ptr1, base_ptr2, base_ptr3);
+
+	float dists[4];
+	mat_trans_vec4_nodiv_nomod(sx, sy, sz, sw,
+	                           dists[0], dists[1], dists[2], dists[3]);
+
+	if(s->radius < dists[0])
+		return SPHEREOUTSIDE;
+	else if(s->radius > -dists[0])
+		res = SPHEREBOUNDARY;
+
+	if(s->radius < dists[1])
+		return SPHEREOUTSIDE;
+	else if(s->radius > -dists[1])
+		res = SPHEREBOUNDARY;
+
+	if(s->radius < dists[2])
+		return SPHEREOUTSIDE;
+	else if(s->radius > -dists[2])
+		res = SPHEREBOUNDARY;
+
+	if(s->radius < dists[3])
+		return SPHEREOUTSIDE;
+	else if(s->radius > -dists[3])
+		res = SPHEREBOUNDARY;
+
+#endif
 	return res;
 }

@@ -445,7 +528,7 @@ Camera::frustumTestSphereNear(const Sphere *s) const
 {
 	int32 res = SPHEREINSIDE;
 	const FrustumPlane *p = this->frustumPlanes;
-
+#ifndef DC_SH4
 	// far
 	float32 dist = dot(p->plane.normal, s->center) - p->plane.distance;
 	if(s->radius < dist)
@@ -481,6 +564,75 @@ Camera::frustumTestSphereNear(const Sphere *s) const
 		return SPHEREOUTSIDE;
 	p++;

+#else
+	__builtin_prefetch(p);
+
+	register float sx asm("fr0") = s->center.x;
+	register float sy asm("fr1") = s->center.y;
+	register float sz asm("fr2") = s->center.z;
+	register float sw asm("fr3") = -1.0f;
+
+	// far
+	register float px asm("fr4") = p->plane.normal.x;
+	register float py asm("fr5") = p->plane.normal.y;
+	register float pz asm("fr6") = p->plane.normal.z;
+	register float pw asm("fr7") = p->plane.distance;
+
+	asm volatile("fipr fv0, fv4"
+				 : "+f" (pw)
+				 : "f" (sx), "f" (sy), "f" (sz), "f" (sw),
+				   "f" (px), "f" (py), "f" (pz));
+	if(s->radius < pw)
+		return SPHEREOUTSIDE;
+	p++;
+
+	// near
+	px = p->plane.normal.x;
+	py = p->plane.normal.y;
+	pz = p->plane.normal.z;
+	pw = p->plane.distance;
+	asm volatile("fipr fv0, fv4"
+				 : "+f" (pw)
+				 : "f" (sx), "f" (sy), "f" (sz), "f" (sw),
+				   "f" (px), "f" (py), "f" (pz));
+	if(s->radius < pw)
+		return SPHEREOUTSIDE;
+	if(s->radius > -pw)
+		res = SPHEREBOUNDARY_NEAR;
+	p++;
+
+	const float* base_ptr0 = &p[0].plane.normal.x;
+	const float* base_ptr1 = &p[1].plane.normal.x;
+	const float* base_ptr2 = &p[2].plane.normal.x;
+	const float* base_ptr3 = &p[3].plane.normal.x;
+
+	__builtin_prefetch(base_ptr0);
+
+	static_assert(offsetof (decltype (p[0].plane.normal), y)
+				 -offsetof (decltype (p[0].plane.normal), x) == sizeof (float));
+
+	static_assert(offsetof (decltype (p[0].plane.normal), z)
+				  -offsetof (decltype (p[0].plane.normal), y) == sizeof (float));
+
+	static_assert(offsetof (decltype (p[0].plane), distance)
+				 -offsetof (decltype (p[0].plane.normal), z) == sizeof (float));
+
+    mat_load_rows(base_ptr0, base_ptr1, base_ptr2, base_ptr3);
+
+	float dists[4];
+	mat_trans_vec4_nodiv_nomod(sx, sy, sz, sw,
+	                           dists[0], dists[1], dists[2], dists[3]);
+
+	if(s->radius < dists[0])
+		return SPHEREOUTSIDE;
+	else if(s->radius < dists[1])
+		return SPHEREOUTSIDE;
+	else if(s->radius < dists[2])
+		return SPHEREOUTSIDE;
+	else if(s->radius < dists[3])
+		return SPHEREOUTSIDE;
+
+#endif
 	return res;
 }

--- a/vendor/librw/src/dc/rwdc.cpp
+++ b/vendor/librw/src/dc/rwdc.cpp
@@ -37,8 +37,6 @@ extern const char* currentFile;
 #include <functional>
 #include <fstream>

-#define ARRAY_SIZE(array)                (sizeof(array) / sizeof(array[0]))
-
 #define errorf(...) dbglog(DBG_CRITICAL, __VA_ARGS__)
 #define logf(...) // printf(__VA_ARGS__)
 bool re3RemoveLeastUsedModel();
@@ -47,21 +45,13 @@ void* re3StreamingAlloc(size_t size);

 // #include "rwdcimpl.h"

-#include <dc/pvr.h>
-#include <dc/matrix.h>
 #include "alloc.h"

-#undef PVR_TXRFMT_STRIDE
-#define PVR_TXRFMT_STRIDE       (1 << 25)
-
-static_assert(PVR_TXRFMT_STRIDE == (1 << 25), "PVR_TXRFMT_STRIDE is bugged in your KOS version");
+using namespace dc;

 // TODO: probably needs a better place to be
 bool doEnvironmentMaps = true;

-#define fclamp0_1(n) ((n) > 1.0f ? 1.0f : n < 0.0f ? 0.0f : n)
-#define fclamp1(n) ((n) > 1.0f ? 1.0f : n)
-
 struct alignas(32) pvr_vertex16_t {
 	uint32_t flags;			/**< \brief TA command (vertex flags) */
 	float	 x;				/**< \brief X coordinate */
@@ -170,177 +160,10 @@ struct alignas(32) pvr_vertex32_ut {
 static_assert(sizeof(pvr_vertex16_t) == 32, "pvr_vertex16_t size mismatch");
 static_assert(alignof(pvr_vertex16_t) == 32, "pvr_vertex16_t alignof mismatch");

-
-#define MATH_Fast_Invert(x) ({ (((x) < 0.0f)? -1.0f : 1.0f) * frsqrt((x) * (x)); }) 
-
 static pvr_dr_state_t drState;

-#include <kos/dbglog.h>
-
 float VIDEO_MODE_SCALE_X;

-#if !defined(DC_TEXCONV) && !defined(DC_SIM)
-#include <kos.h>
-
-#define mat_trans_nodiv_nomod(x, y, z, x2, y2, z2, w2) do { \
-        register float __x __asm__("fr12") = (x); \
-        register float __y __asm__("fr13") = (y); \
-        register float __z __asm__("fr14") = (z); \
-        register float __w __asm__("fr15") = 1.0f; \
-        __asm__ __volatile__( "ftrv  xmtrx, fv12\n" \
-                              : "=f" (__x), "=f" (__y), "=f" (__z), "=f" (__w) \
-                              : "0" (__x), "1" (__y), "2" (__z), "3" (__w) ); \
-        x2 = __x; y2 = __y; z2 = __z; w2 = __w; \
-    } while(false)
-
-#define mat_trans_nodiv_nomod_zerow(x, y, z, x2, y2, z2, w2) do { \
-        register float __x __asm__("fr12") = (x); \
-        register float __y __asm__("fr13") = (y); \
-        register float __z __asm__("fr14") = (z); \
-        register float __w __asm__("fr15") = 0.0f; \
-        __asm__ __volatile__( "ftrv  xmtrx, fv12\n" \
-                              : "=f" (__x), "=f" (__y), "=f" (__z), "=f" (__w) \
-                              : "0" (__x), "1" (__y), "2" (__z), "3" (__w) ); \
-        x2 = __x; y2 = __y; z2 = __z; w2 = __w; \
-    } while(false)
-
-#define mat_trans_w_nodiv_nomod(x, y, z, w) do { \
-        register float __x __asm__("fr12") = (x); \
-        register float __y __asm__("fr13") = (y); \
-        register float __z __asm__("fr14") = (z); \
-        register float __w __asm__("fr15") = 1.0f; \
-        __asm__ __volatile__( "ftrv  xmtrx, fv12\n" \
-                              : "=f" (__x), "=f" (__y), "=f" (__z), "=f" (__w) \
-                              : "0" (__x), "1" (__y), "2" (__z), "3" (__w) ); \
-        w = __w; \
-    } while(false)
-
-	// no declspec naked, so can't do rts / fschg. instead compiler pads with nop?
-
-	inline void rw_mat_load_3x3(const rw::Matrix* mtx) {
-		__asm__ __volatile__ (
-			R"(
-				fschg
-				frchg
-
-				fmov        @%[mtx]+, dr0
-
-				fldi0 		fr12
-				fldi0 		fr13
-
-				fmov        @%[mtx]+, dr2
-				fmov        @%[mtx]+, dr4
-				fmov        @%[mtx]+, dr6
-				fmov        @%[mtx]+, dr8
-				fmov        @%[mtx]+, dr10
-
-				fldi0	    fr3
-				fldi0	    fr7
-				fldi0	    fr11
-				fmov        dr12, dr14
-
-				fschg
-				frchg
-			)"
-			: [mtx] "+r" (mtx)
-		);
-	}
-
-	// sets pos.w to 1
-	inline void rw_mat_load_4x4(const rw::Matrix* mtx) {
-		__asm__ __volatile__ (
-			R"(
-				fschg
-				frchg
-				fmov        @%[mtx]+, dr0
-
-				fmov        @%[mtx]+, dr2
-				fmov        @%[mtx]+, dr4
-				fmov        @%[mtx]+, dr6
-				fmov        @%[mtx]+, dr8
-				fmov        @%[mtx]+, dr10
-				fmov        @%[mtx]+, dr12
-				fmov        @%[mtx]+, dr14
-				fldi1 	 	fr15
-
-				fschg
-				frchg
-			)"
-			: [mtx] "+r" (mtx)
-		);
-	}
-
-#else
-extern matrix_t XMTRX;
-
-void rw_mat_load_3x3(rw::Matrix* mtx) {
-	memcpy(XMTRX, mtx, sizeof(matrix_t));
-	XMTRX[0][3] = 0.0f;
-	XMTRX[1][3] = 0.0f;
-	XMTRX[2][3] = 0.0f;
-
-	XMTRX[3][0] = 0.0f;
-	XMTRX[3][1] = 0.0f;
-	XMTRX[3][2] = 0.0f;
-	XMTRX[3][3] = 0.0f;
-}
-
-void rw_mat_load_4x4(rw::Matrix* mtx) {
-	memcpy(XMTRX, mtx, sizeof(matrix_t));
-	XMTRX[3][3] = 1.0f;
-}
-
-#include <dc/matrix.h>
-#define frsqrt(a) 				(1.0f/sqrt(a))
-#define dcache_pref_block(a)	__builtin_prefetch(a)
-
-#ifndef __always_inline
-#define __always_inline 		__attribute__((always_inline)) inline
-#endif
-
-#ifdef DC_TEXCONV
-#define mat_transform(a, b, c, d)
-#define mat_apply(a)
-#define mat_load(a)
-#define mat_store(a)
-#define mat_identity(a)
-#define pvr_fog_table_color(a,r,g,b)
-#define pvr_fog_table_linear(s,e)
-#define pvr_fog_table_exp(d)
-#define pvr_fog_table_custom(d)
-#endif
-
-#define mat_trans_single3_nomod(x_, y_, z_, x2, y2, z2) do { \
-		vector_t tmp = { x_, y_, z_, 1.0f }; \
-		mat_transform(&tmp, &tmp, 1, 0); \
-		z2 = 1.0f / tmp.w; \
-		x2 = tmp.x * z2; \
-		y2 = tmp.y * z2; \
-	} while(false)
-
-#define mat_trans_nodiv_nomod(x_, y_, z_, x2, y2, z2, w2) do { \
-		vector_t tmp1233123 = { x_, y_, z_, 1.0f }; \
-		mat_transform(&tmp1233123, &tmp1233123, 1, 0); \
-		x2 = tmp1233123.x; y2 = tmp1233123.y; z2 = tmp1233123.z; w2 = tmp1233123.w; \
-	} while(false)
-
-#define mat_trans_nodiv_nomod_zerow(x_, y_, z_, x2, y2, z2, w2) do { \
-		vector_t tmp1233123 = { x_, y_, z_, 0.0f }; \
-		mat_transform(&tmp1233123, &tmp1233123, 1, 0); \
-		x2 = tmp1233123.x; y2 = tmp1233123.y; z2 = tmp1233123.z; w2 = tmp1233123.w; \
-	} while(false)
-
-#define mat_trans_w_nodiv_nomod(x_, y_, z_, w_) do { \
-		vector_t tmp1233123 = { x_, y_, z_, 1.0f }; \
-		mat_transform(&tmp1233123, &tmp1233123, 1, 0); \
-		w_ = tmp1233123.w; \
-	} while(false)
-
-#define memcpy4 memcpy
-
-// END STUBS
-#endif
-
 static pvr_ptr_t fake_tex;

 alignas(4) static const uint16_t fake_tex_data[] = {
@@ -551,9 +374,9 @@ void DCE_MatrixViewport(float x, float y, float width, float height) {

 void DCE_InitMatrices() {
 	// Setup the screenview matrix.  Only need to do once since this matrix does not need to change for single player viewpoint.
-	mat_identity();
+	mat_identity2();
 	
-	mat_store(&DCE_MAT_SCREENVIEW);
+	mat_store2(&DCE_MAT_SCREENVIEW);
 }

 }
@@ -684,7 +507,7 @@ struct atomic_context_t {
 __always_inline void DCE_RenderSubmitVertex(const pvr_vertex_t *v, uint32_t flags) {
    auto *sq  = reinterpret_cast<uint32_t *>(pvr_dr_target(drState));
    auto *src = reinterpret_cast<const uint32_t *>(v);
-    float sz  = MATH_Fast_Invert(v->z);
+    float sz  = Invert<true, false>(v->z);
    float sx  = v->x * sz;
    float sy  = v->y * sz;
    
@@ -711,7 +534,7 @@ __always_inline void DCE_RenderSubmitVertexIM3D(float x, float y, float w,
 {
    auto *sq   = reinterpret_cast<uint32_t *>(pvr_dr_target(drState));
    auto *uv32 = reinterpret_cast<const uint32_t *>(uv);
-    float sz   = MATH_Fast_Invert(w);
+    float sz   = Invert<true, false>(w);
    float sx   = x * sz;
    float sy   = y * sz;

@@ -728,13 +551,6 @@ __always_inline void DCE_RenderSubmitVertexIM3D(float x, float y, float w,

 /* END TA Submission Functions*/

-
-
-
-#if defined(DC_TEXCONV)
-void malloc_stats() { }
-#endif
-
 #if 0
 #define UNIMPL_LOG() printf("TODO: Implement %s @ %s:%d\n", __func__, __FILE__, __LINE__);
 #define UNIMPL_LOGV(fmt, ...) printf("TODO: Implement %s @ %s:%d " fmt "\n", __func__, __FILE__, __LINE__, __VA_ARGS__);
@@ -747,7 +563,7 @@ Camera* rwdcCam;

 void beginUpdate(Camera* cam)  {
 	rwdcCam = cam;
-	float view[16], proj[16];
+	alignas(8) float view[16], proj[16];

 	// View Matrix
 	Matrix inv;
@@ -770,7 +586,7 @@ void beginUpdate(Camera* cam)  {
 	view[13] =  -inv.pos.y;
 	view[14] =  inv.pos.z;
 	view[15] =  1.0f;
-	memcpy4(&cam->devView, view, sizeof(RawMatrix));
+	mat_copy(cam->devView, reinterpret_cast<const matrix_t *>(view));
 //	d3ddevice->SetTransform(D3DTS_VIEW, (D3DMATRIX*)view);

 	// Projection Matrix
@@ -804,13 +620,12 @@ void beginUpdate(Camera* cam)  {
 		proj[15] = 1.0f;
 	}
 	proj[14] = -cam->nearPlane*proj[10];
-	memcpy4(&cam->devProj, proj, sizeof(RawMatrix));
-	
+    mat_copy(cam->devProj, reinterpret_cast<const matrix_t *>(proj));
+
 	DCE_MatrixViewport(0, 0, cam->frameBuffer->width * VIDEO_MODE_SCALE_X, cam->frameBuffer->height);
 	
-	mat_load((matrix_t*)&DCE_MAT_SCREENVIEW);
-	mat_apply((matrix_t*)&cam->devProj);
-	mat_store((matrix_t*)&cam->devProjScreen);
+	mat_load_apply((matrix_t*)&DCE_MAT_SCREENVIEW, (matrix_t*)&cam->devProj);
+	mat_store2((matrix_t*)&cam->devProjScreen);

 }

@@ -996,7 +811,8 @@ struct chunked_vector {
    }

    // Iterate over each element and invoke the callback.
-    void forEach(void(*cb)(T&)) {
+    template<typename F>
+    void forEach(F&& cb) {
        for (chunk* curr = first; curr; curr = curr->header.next) {
            for (size_t i = 0; i < curr->header.used; ++i) {
                cb(curr->items[i]);
@@ -1130,8 +946,6 @@ void dcMotionBlur_v1(uint8_t a, uint8_t r, uint8_t g, uint8_t b) {
 		auto addr2 = (pvr_ptr_t)&emu_vram[addr64b + 640 * 2];
 	#endif

-
-
 		PVR_SET(PVR_TEXTURE_MODULO, 640/32);

 		auto doquad = [=](float x, float y, float w, float h, float tx, float ty, float tw, float th) {
@@ -1838,7 +1652,7 @@ void im2DRenderPrimitive(PrimitiveType primType, void *vertices, int32_t numVert
 			pvrVert->flags = flags;
 			pvrVert->x 	   = gtaVert.x * VIDEO_MODE_SCALE_X;
 			pvrVert->y	   = gtaVert.y;
-			pvrVert->z 	   = MATH_Fast_Invert(gtaVert.w); // this is perfect for almost every case...
+			pvrVert->z 	   = Invert<true, false>(gtaVert.w); // this is perfect for almost every case...
 			pvrVert->u 	   = gtaVert.u;
 			pvrVert->v 	   = gtaVert.v;
 			pvrVert->argb  = (gtaVert.a << 24) |
@@ -1851,13 +1665,13 @@ void im2DRenderPrimitive(PrimitiveType primType, void *vertices, int32_t numVert
 		switch(primType) {
 			case PRIMTYPETRILIST:
 				pvrHeaderSubmit();
-				dcache_pref_block(vtx);
+				__builtin_prefetch(vtx);
 				for(int i = 0; i < numVertices; i += 3) [[likely]] {
-					dcache_pref_block(&vtx[i + 1]);
+					__builtin_prefetch(&vtx[i + 1]);
 					pvrVertexSubmit(vtx[i + 0], PVR_CMD_VERTEX);
-					dcache_pref_block(&vtx[i + 2]);
+					__builtin_prefetch(&vtx[i + 2]);
 					pvrVertexSubmit(vtx[i + 1], PVR_CMD_VERTEX);
-					dcache_pref_block(&vtx[i + 3]);
+					__builtin_prefetch(&vtx[i + 3]);
 					pvrVertexSubmit(vtx[i + 2], PVR_CMD_VERTEX_EOL);
 				}
 			break;
@@ -1865,14 +1679,11 @@ void im2DRenderPrimitive(PrimitiveType primType, void *vertices, int32_t numVert
 				pvrHeaderSubmit();
 				const auto *vtxA = vtx + 0;
 				const auto *vtxB = vtx + 1;
-				dcache_pref_block(vtxA);
+				__builtin_prefetch(vtxA);
 				for(int i = 2; i < numVertices; ++i) [[likely]] {
 					const auto *vtxC = vtx + i;
-					dcache_pref_block(vtxB);
 					pvrVertexSubmit(*vtxA, PVR_CMD_VERTEX);
-					dcache_pref_block(vtxC);
 					pvrVertexSubmit(*vtxB, PVR_CMD_VERTEX);
-					dcache_pref_block(&vtx[i]);
 					pvrVertexSubmit(*vtxC, PVR_CMD_VERTEX_EOL);
 					vtxB = vtxC;
 				}
@@ -1959,7 +1770,7 @@ void im2DRenderIndexedPrimitive(PrimitiveType primType, void *vertices, int32 nu
 			pvrVert->flags = flags;
 			pvrVert->x 	   = gtaVert.x * VIDEO_MODE_SCALE_X;
 			pvrVert->y	   = gtaVert.y;
-			pvrVert->z 	   = MATH_Fast_Invert(gtaVert.w); // this is perfect for almost every case...
+			pvrVert->z 	   = Invert<true, false>(gtaVert.w); // this is perfect for almost every case...
 			pvrVert->u 	   = gtaVert.u;
 			pvrVert->v 	   = gtaVert.v;
 			pvrVert->argb  = (gtaVert.a << 24) |
@@ -1972,13 +1783,13 @@ void im2DRenderIndexedPrimitive(PrimitiveType primType, void *vertices, int32 nu
 		switch(primType) {
 			case PRIMTYPETRILIST:
 				pvrHeaderSubmit();
-				dcache_pref_block(vtx);
+				__builtin_prefetch(vtx);
 				for(int i = 0; i < numIndices; i += 3) [[likely]] {
-					dcache_pref_block(&vtx[idx[i + 1]]);
+					__builtin_prefetch(&vtx[idx[i + 1]]);
 					pvrVertexSubmit(vtx[idx[i + 0]], PVR_CMD_VERTEX);
-					dcache_pref_block(&vtx[idx[i + 2]]);
+					__builtin_prefetch(&vtx[idx[i + 2]]);
 					pvrVertexSubmit(vtx[idx[i + 1]], PVR_CMD_VERTEX);
-					dcache_pref_block(&vtx[idx[i + 3]]);
+					__builtin_prefetch(&vtx[idx[i + 3]]);
 					pvrVertexSubmit(vtx[idx[i + 2]], PVR_CMD_VERTEX_EOL);
 				}
 			break;
@@ -2015,15 +1826,21 @@ void im3DTransform(void *vertices, int32 numVertices, Matrix *worldMat, uint32 f
 		worldMat = &ident;
 	}
 	
-	rw::RawMatrix mtx, proj, world, worldview;
 	rw::Camera *cam = engine->currentCamera;
-
-	rw::convMatrix(&world, worldMat);
+#ifndef DC_SH4
+    rw::RawMatrix mtx, proj, world, worldview;
+    rw::convMatrix(&world, worldMat);
 	rw::RawMatrix::mult(&worldview, &world, &cam->devView);
 	rw::RawMatrix::mult(&proj, &worldview, &cam->devProj);
 	rw::RawMatrix::mult(&mtx, &proj, (RawMatrix*)&DCE_MAT_SCREENVIEW);
-	// mat_load(&DCE_MAT_SCREENVIEW);     // ~11 cycles.
-	mat_load(( matrix_t*)&mtx.right);  // Number of cycles: ~32.
+	// mat_load2(&DCE_MAT_SCREENVIEW);     // ~11 cycles.
+	mat_load2(( matrix_t*)&mtx.right);  // Number of cycles: ~32.
+#else
+    mat_load_apply(&DCE_MAT_SCREENVIEW, cam->devProj);
+    mat_apply(cam->devView);
+    mat_apply(*worldMat);
+#endif
+
    if (im3dVertices) {
 		free(im3dVertices);
 	}
@@ -2117,7 +1934,7 @@ void im3DRenderIndexedPrimitive(PrimitiveType primType,

 			// assuming near plane is 0.0f
 			// gv1 is visible (posi), and gv2 is behind the plane (negative)
-			float t = (1.0f - gv1.position.z) * MATH_Fast_Invert(gv2.position.z - gv1.position.z);
+			float t = (1.0f - gv1.position.z) * Invert<true, true>(gv2.position.z - gv1.position.z);

 			pvr_vertex_t pvrVert; 

@@ -2321,39 +2138,36 @@ static_assert(sizeof(MeshletInfo) == 40); // or 32 if !skin

 inline __attribute__((always_inline))  void setLights(Atomic *atomic, WorldLights *lightData, UniformObject &uniformObject)
 {
-	int n = 0;
+    int i = 0;

 	uniformObject.ambLight = lightData->ambient;

 	if (lightData->numDirectionals) {
 		Matrix mat;
 		Matrix matsrc = *atomic->getFrame()->getLTM();
-		matsrc.pos = V3d {0,0,0};
+		matsrc.pos = V3d {0.0f, 0.0f, 0.0f};
+        matsrc.posw = 0.0f;
 		Matrix::invert(&mat, &matsrc);
-		
-		n = 0;
-		for(int i = 0; i < lightData->numDirectionals && i < MAX_LIGHTS; i++){
+		mat_load2(mat);
+
+		for(; i < lightData->numDirectionals && i < MAX_LIGHTS; i++){
 			Light *l = lightData->directionals[i];
-			uniformObject.col[n] = scale(l->color, 255);
+			uniformObject.col[i] = scale(l->color, 255);
 			V3d at = l->getFrame()->getLTM()->at;

 			V3d dir;

-			V3d::transformVectors(&dir, &at, 1, &mat);
+			mat_trans_normal3_nomod(at.x, at.y, at.z,
+                                    dir.x, dir.y, dir.z);

-			uniformObject.dir[n>>2][0][n&3] = -dir.x / 127.0f;
-			uniformObject.dir[n>>2][1][n&3] = -dir.y / 127.0f;
-			uniformObject.dir[n>>2][2][n&3] = -dir.z / 127.0f;
-			uniformObject.dir[n>>2][3][n&3] = 0;
-
-			n++;
-			if(n >= MAX_LIGHTS)
-				goto out;
+			uniformObject.dir[i>>2][0][i&3] = -dir.x / 127.0f;
+			uniformObject.dir[i>>2][1][i&3] = -dir.y / 127.0f;
+			uniformObject.dir[i>>2][2][i&3] = -dir.z / 127.0f;
+			uniformObject.dir[i>>2][3][i&3] = 0;
 		}
 	}

-out:
-	uniformObject.lightCount = n;
+	uniformObject.lightCount = i;
 }


@@ -3020,7 +2834,7 @@ void* interpolateAndSubmit(void* dst, const void* src1, const void* src2, uint32
 	float y = v1->o_r + t * (v2->o_r - v1->o_r);
 	float w = v1->o_g + t * (v2->o_g - v1->o_g);

-	w = frsqrt(w * w);
+	w = Invert<true, false>(w);

 	v->x = x * w;
 	v->y = y * w;
@@ -3394,15 +3208,14 @@ void tnlMeshletSkinVertices(uint8_t *OCR, uint8_t *OCR_normal, const uint8_t* ve
 		auto skinningIndexData = (int16_t*)skinIndexes;
 		auto skinningWeightData = (uint8_t*)skinWeights;

-		if (!matrix0Identity) {
-			rw_mat_load_4x4(&skinMatrices[0]);
-			if (small_xyz) {
-				mat_apply(&DCE_MESHLET_MAT_DECODE);
-			}
-		} else {
-			if (small_xyz) {
-				mat_load(&DCE_MESHLET_MAT_DECODE);
-			}
+		if constexpr (!matrix0Identity) {
+			if (!small_xyz)
+                mat_load2(skinMatrices[0]);
+			else 
+                mat_load_apply(skinMatrices[0], &DCE_MESHLET_MAT_DECODE);
+			
+		} else if constexpr (small_xyz) {
+			mat_load2(&DCE_MESHLET_MAT_DECODE);
 		}

 		for(;;) {
@@ -3412,7 +3225,7 @@ void tnlMeshletSkinVertices(uint8_t *OCR, uint8_t *OCR_normal, const uint8_t* ve
 				int count = *skinningIndexData++;
 				uint8_t* dstVertexBytes = dest + *skinningIndexData++;

-				if (matrix0Identity && !small_xyz) {
+				if constexpr (matrix0Identity && !small_xyz) {
 					do {
 						const V3d* srcVtx = (const V3d*)(srcVtxBytes);
 						srcVtxBytes += vertexSize;
@@ -3431,11 +3244,8 @@ void tnlMeshletSkinVertices(uint8_t *OCR, uint8_t *OCR_normal, const uint8_t* ve
 						srcVtxBytes += vertexSize;
 						V3d* dstVertex = (V3d*)(dstVertexBytes);
 						dstVertexBytes += 64;
-						float x, y, z, w;
-						mat_trans_nodiv_nomod(srcVtx->x, srcVtx->y, srcVtx->z, x, y, z, w);
-						dstVertex->x = x;
-						dstVertex->y = y;
-						dstVertex->z = z;
+						mat_trans_single3_nodiv_nomod(srcVtx->x, srcVtx->y, srcVtx->z,
+                                                      dstVertex->x, dstVertex->y, dstVertex->z);
 					} while(--count != 0);
 				}
 			} else if (!(flags & 0x80)) {
@@ -3465,10 +3275,10 @@ void tnlMeshletSkinVertices(uint8_t *OCR, uint8_t *OCR_normal, const uint8_t* ve
 				break;
 			}

-			rw_mat_load_4x4(currentMatrix);
-			if (small_xyz){
-				mat_apply(&DCE_MESHLET_MAT_DECODE);
-			}
+			if constexpr(!small_xyz)
+                mat_load2(*currentMatrix);
+            else
+				mat_load_apply(*currentMatrix, &DCE_MESHLET_MAT_DECODE);

 			do {
 				auto srcOffset = *skinningIndexData++;
@@ -3482,9 +3292,9 @@ void tnlMeshletSkinVertices(uint8_t *OCR, uint8_t *OCR_normal, const uint8_t* ve
 					srcVtx = &tmpSrc;
 				}
 				auto dstVtx = (V3d*)(dest + dstOffset);
-				float x, y, z, w;
-				mat_trans_nodiv_nomod(srcVtx->x, srcVtx->y, srcVtx->z, x, y, z, w);
-				V3d tmp = { x, y, z };
+				V3d tmp;
+				mat_trans_single3_nodiv_nomod(srcVtx->x, srcVtx->y, srcVtx->z, 
+                                              tmp.x, tmp.y, tmp.z);
 				tmp = scale(tmp, *skinningWeightData++ / 255.0f);
 				*dstVtx = add(*dstVtx, tmp);
 			} while (--count != 0);
@@ -3497,8 +3307,8 @@ void tnlMeshletSkinVertices(uint8_t *OCR, uint8_t *OCR_normal, const uint8_t* ve
 		auto skinningIndexData = (int16_t*)skinIndexes;
 		auto skinningWeightData = (uint8_t*)skinWeights;

-		if (!matrix0Identity) {
-			rw_mat_load_3x3(&skinMatrices[0]);
+		if constexpr(!matrix0Identity) {
+			mat_load2(skinMatrices[0]);
 		}

 		for(;;) {
@@ -3508,7 +3318,7 @@ void tnlMeshletSkinVertices(uint8_t *OCR, uint8_t *OCR_normal, const uint8_t* ve
 				int count = *skinningIndexData++;
 				uint8_t* dstNormalBytes = destNormal + *skinningIndexData++;

-				if (matrix0Identity) {
+				if constexpr (matrix0Identity) {
 					do {
 						V3d srcNormal = { static_cast<float32>(srcNormalBytes[0]), static_cast<float32>(srcNormalBytes[1]), static_cast<float32>(srcNormalBytes[2]) };
 						
@@ -3524,9 +3334,9 @@ void tnlMeshletSkinVertices(uint8_t *OCR, uint8_t *OCR_normal, const uint8_t* ve
 						srcNormalBytes += vertexSize;
 						V3d* dstNormal = (V3d*)(dstNormalBytes);
 						dstNormalBytes += 64;
-						float x, y, z, w;
-						mat_trans_nodiv_nomod_zerow(srcNormal.x, srcNormal.y, srcNormal.z, x, y, z, w);
-						*dstNormal = { x, y, z };
+						float x, y, z;
+						mat_trans_normal3_nomod(srcNormal.x, srcNormal.y, srcNormal.z,
+                                                dstNormal->x, dstNormal->y, dstNormal->z);
 					} while(--count != 0);
 				}
 			} else if (!(flags & 0x80)) {
@@ -3556,7 +3366,7 @@ void tnlMeshletSkinVertices(uint8_t *OCR, uint8_t *OCR_normal, const uint8_t* ve
 				break;
 			}

-			rw_mat_load_3x3(currentMatrix);
+			mat_load2(*currentMatrix);

 			do {
 				auto srcOffset = *skinningIndexData++;
@@ -3568,8 +3378,8 @@ void tnlMeshletSkinVertices(uint8_t *OCR, uint8_t *OCR_normal, const uint8_t* ve
 				auto dstNormal = (V3d*)(destNormal + dstOffset);

 				V3d tmp;
-				float w;
-				mat_trans_nodiv_nomod_zerow(srcNormal.x, srcNormal.y, srcNormal.z, tmp.x, tmp.y, tmp.z, w);
+				mat_trans_normal3_nomod(srcNormal.x, srcNormal.y, srcNormal.z,
+                                        tmp.x, tmp.y, tmp.z);
 				tmp = scale(tmp, *skinningWeightData++ / 255.0f);
 				*dstNormal = add(*dstNormal, tmp);
 			} while (--count != 0);
@@ -3581,7 +3391,7 @@ void tnlMeshletSkinVertices(uint8_t *OCR, uint8_t *OCR_normal, const uint8_t* ve
 __attribute__((noinline))
 void tnlMeshletEnvMap(uint8_t* OCR, uint8_t* normal, int vertexCount, int vertexSize, matrix_t* matfxMatrix, float matfxCoefficient) {

-	mat_load(matfxMatrix);
+	mat_load2(matfxMatrix);

 	do {
 		pvr_vertex64_t* v = (pvr_vertex64_t*)OCR;
@@ -3683,22 +3493,20 @@ uploadSkinMatrices(Atomic *a, Matrix* skinMatrices)

 	if(hier){
 		Matrix *invMats = (Matrix*)skin->inverseMatrices;
-		Matrix tmp;

 		assert(skin->numBones == hier->numNodes);
 		if(hier->flags & HAnimHierarchy::LOCALSPACEMATRICES){
 			for(i = 0; i < hier->numNodes; i++){
-				invMats[i].flags = 0;
-				Matrix::mult(m, &invMats[i], &hier->matrices[i]);
+				mat_mult(*m, invMats[i], hier->matrices[i]);
 				m++;
 			}
 		}else{
 			Matrix invAtmMat;
 			Matrix::invert(&invAtmMat, a->getFrame()->getLTM());
 			for(i = 0; i < hier->numNodes; i++){
-				invMats[i].flags = 0;
-				Matrix::mult(&tmp, &hier->matrices[i], &invAtmMat);
-				Matrix::mult(m, &invMats[i], &tmp);
+				mat_load_apply(invAtmMat, hier->matrices[i]);
+				mat_apply(invMats[i]);
+                mat_store2(*m);
 				m++;
 			}
 		}
@@ -3715,12 +3523,12 @@ uploadSkinMatrices(Atomic *a, Matrix* skinMatrices)
 	return skinMatrices[0].identityError() < 0.01f;
 }

-static RawMatrix normal2texcoord = {
+static RawMatrix normal2texcoord = {{
 	{ 0.5f / 127,  0.0f, 0.0f }, 0.0f,
 	{ 0.0f, -0.5f / 127, 0.0f }, 0.0f,
 	{ 0.0f,  0.0f, 1.0f }, 0.0f,
 	{ 0.5f,  0.5f, 0.0f }, 1.0f
-};
+}};

 void
 uploadEnvMatrix(Frame *frame, RawMatrix *world, matrix_t* envMatrix)
@@ -3732,13 +3540,13 @@ uploadEnvMatrix(Frame *frame, RawMatrix *world, matrix_t* envMatrix)
 	RawMatrix *envMtx = (RawMatrix*)envMatrix;
 	{

-		RawMatrix invMtx;
+		//RawMatrix invMtx;
 		Matrix::invert(&invMat, frame->getLTM());
-		convMatrix(&invMtx, &invMat);
-		invMtx.pos.set(0.0f, 0.0f, 0.0f);
+		//convMatrix(&invMtx, &invMat);
+		//invMtx.pos.set(0.0f, 0.0f, 0.0f);
 		float uscale = fabs(normal2texcoord.right.x);
 		normal2texcoord.right.x = MatFX::envMapFlipU ? -uscale : uscale;
-		
+#if 0
 		RawMatrix tmpMtx;

 		RawMatrix::mult(&tmpMtx, &invMtx, &normal2texcoord);
@@ -3748,6 +3556,11 @@ uploadEnvMatrix(Frame *frame, RawMatrix *world, matrix_t* envMatrix)
 		world->upw = 0;
 		world->atw = 0;
 		RawMatrix::mult(envMtx, world, &tmpMtx);
+#else
+        mat_load_apply(normal2texcoord, invMat);
+        mat_apply(*world);
+        mat_store2(envMatrix);
+#endif
 	}
 }

@@ -3997,16 +3810,10 @@ void defaultRenderCB(ObjPipeline *pipe, Atomic *atomic) {
 	ac->skinMatrix0Identity = skinMatrix0Identity;

 	lightingCB(atomic, ac->uniform);
-
-
-	rw::RawMatrix world;
-	rw::convMatrix(&world, atomic->getFrame()->getLTM());
 	
-
-	mat_load((matrix_t*)&cam->devProjScreen);
-	mat_apply((matrix_t*)&cam->devView);
-	mat_apply((matrix_t*)&world);
-	mat_store((matrix_t*)&atomicContexts.back().mtx);
+	mat_load_apply((matrix_t*)&cam->devProjScreen, (matrix_t*)&cam->devView);
+	mat_apply(*atomic->getFrame()->getLTM());
+	mat_store2((matrix_t*)&atomicContexts.back().mtx);

 	auto meshes = geo->meshHeader->getMeshes();

@@ -4029,7 +3836,7 @@ void defaultRenderCB(ObjPipeline *pipe, Atomic *atomic) {
 			matfxContextPointer = &matfxContexts.back();
 			// N.B. world here gets converted to a 3x3 matrix
 			// 		this is fine, as we only use it for env mapping from now on
-			uploadEnvMatrix(matfx->fx[0].env.frame, &world, &matfxContexts.back().mtx);
+			uploadEnvMatrix(matfx->fx[0].env.frame, reinterpret_cast<rw::RawMatrix*>(atomic->getFrame()->getLTM()), &matfxContexts.back().mtx);
 			matfxContextPointer->coefficient = matfxCoefficient;
 			
 			pvr_poly_cxt_t cxt;
@@ -4245,15 +4052,14 @@ void defaultRenderCB(ObjPipeline *pipe, Atomic *atomic) {
 						unsigned skinSelector = small_xyz + acp->skinMatrix0Identity*2;
 						tnlMeshletSkinVerticesSelector[skinSelector](OCR_SPACE, normalDst, &dcModel->data[meshlet->vertexOffset],  normalSrc, &dcModel->data[meshlet->skinWeightOffset], &dcModel->data[meshlet->skinIndexOffset], meshlet->vertexCount, meshlet->vertexSize, &acp->skinContextPointer->mtx);
 						
-						mat_load(&mtx);
+						mat_load2(&mtx);
 						tnlMeshletTransformSelector[clippingRequired * 2](OCR_SPACE, OCR_SPACE + 4, meshlet->vertexCount, 64);
 					} else {
 						
 						if (selector & 8) {
-							mat_load(&mtx);
-							mat_apply(&DCE_MESHLET_MAT_DECODE);
+							mat_load_apply(&mtx, &DCE_MESHLET_MAT_DECODE);
 						} else {
-							mat_load(&mtx);
+							mat_load2(&mtx);
 						}
 						tnlMeshletTransformSelector[smallSelector](OCR_SPACE, &dcModel->data[meshlet->vertexOffset], meshlet->vertexCount, meshlet->vertexSize);
 					}
@@ -4280,7 +4086,7 @@ void defaultRenderCB(ObjPipeline *pipe, Atomic *atomic) {

 						unsigned dstColOffset = textured ? offsetof(pvr_vertex64_t, a) : offsetof(pvr_vertex32_ut, a);
 						dce_set_mat_vertex_color(&residual, &material);
-						mat_load(&DCE_MESHLET_MAT_VERTEX_COLOR);
+						mat_load2(&DCE_MESHLET_MAT_VERTEX_COLOR);
 						tnlMeshletVertexColorSelector[0](OCR_SPACE + dstColOffset, (int8_t*)&dcModel->data[meshlet->vertexOffset] + colOffset, meshlet->vertexCount, meshlet->vertexSize);
 					} else {
 						unsigned dstColOffset = textured ? offsetof(pvr_vertex64_t, a) : offsetof(pvr_vertex32_ut, a);
@@ -4302,7 +4108,7 @@ void defaultRenderCB(ObjPipeline *pipe, Atomic *atomic) {

 						
 						unsigned normalSelector = (pass1 - 1) + (skin != 0) * 4;
-						mat_load((matrix_t*)&uniformObject.dir[0][0][0]);
+						mat_load2((matrix_t*)&uniformObject.dir[0][0][0]);
 						auto normalPointer = &dcModel->data[meshlet->vertexOffset] + normalOffset;
 						auto vtxSize = meshlet->vertexSize;
 						if (skin) {
@@ -4317,7 +4123,7 @@ void defaultRenderCB(ObjPipeline *pipe, Atomic *atomic) {
 					
 						if (pass2) {
 							unsigned normalSelector = (pass2 - 1) + (skin != 0) * 4;
-							mat_load((matrix_t*)&uniformObject.dir[1][0][0]);
+							mat_load2((matrix_t*)&uniformObject.dir[1][0][0]);
 							tnlMeshletDiffuseColorSelector[normalSelector](OCR_SPACE + dstColOffset, normalPointer, meshlet->vertexCount, vtxSize, &lightDiffuseColors[4]);
 						}
 					}
@@ -4380,7 +4186,7 @@ void defaultRenderCB(ObjPipeline *pipe, Atomic *atomic) {
 					indices.back() |= 0x80;

 					pvr_vertex64_t *vd = (pvr_vertex64_t *)OCR_SPACE;
-					mat_load(&mtx);  // Number of cycles: ~11
+					mat_load2(&mtx);  // Number of cycles: ~11

 					for (int idx = 0; idx < geo->numVertices; idx++) {
 						auto& vert = vertices[idx];
@@ -6436,6 +6242,18 @@ writeNativeSkin(Stream *stream, int32 len, void *object, int32 offset)

 	stream->write8(&skin->numBones, 4);

+	for(int32 i = 0; i < skin->numBones; i++){
+		Matrix &m = *reinterpret_cast<Matrix *>(
+			&skin->inverseMatrices[i * 16]);
+
+		if(m.flags & MatrixBase::IDENTITY_OLD)
+			m.flags |= MatrixBase::IDENTITY;
+		m.pad0 = 0;
+		m.upw = 0.0f;
+		m.atw = 0.0f;
+		m.posw = 1.0f;
+	}
+
 	stream->write32(skin->inverseMatrices, skin->numBones*64);
 	return stream;
 }
--- a/vendor/librw/src/dc/rwdc_common.h
+++ b/vendor/librw/src/dc/rwdc_common.h
--- a/vendor/librw/src/rwbase.h
+++ b/vendor/librw/src/rwbase.h
@@ -1,6 +1,9 @@
 #pragma once

 #include "common_defines.h"
+#ifdef RW_DC
+#include "rwdc_common.h"
+#endif

 #ifndef RW_PS2
 #include <stdint.h>
@@ -235,8 +238,8 @@ inline V2d neg(const V2d &a) { return makeV2d(-a.x, -a.y); }
 inline V2d add(const V2d &a, const V2d &b) { return makeV2d(a.x+b.x, a.y+b.y); }
 inline V2d sub(const V2d &a, const V2d &b) { return makeV2d(a.x-b.x, a.y-b.y); }
 inline V2d scale(const V2d &a, float32 r) { return makeV2d(a.x*r, a.y*r); }
-inline float32 length(const V2d &v) { return sqrtf(v.x*v.x + v.y*v.y); }
-inline V2d normalize(const V2d &v) { return scale(v, 1.0f/length(v)); }
+inline float32 length(const V2d &v) { return dc::Sqrt(v.x*v.x + v.y*v.y); }
+inline V2d normalize(const V2d &v) { return scale(v, dc::RecipSqrt(v.x*v.x + v.y*v.y)); }

 struct V3d
 {
@@ -253,11 +256,37 @@ inline V3d neg(const V3d &a) { return makeV3d(-a.x, -a.y, -a.z); }
 inline V3d add(const V3d &a, const V3d &b) { return makeV3d(a.x+b.x, a.y+b.y, a.z+b.z); }
 inline V3d sub(const V3d &a, const V3d &b) { return makeV3d(a.x-b.x, a.y-b.y, a.z-b.z); }
 inline V3d scale(const V3d &a, float32 r) { return makeV3d(a.x*r, a.y*r, a.z*r); }
-inline float32 length(const V3d &v) { return sqrtf(v.x*v.x + v.y*v.y + v.z*v.z); }
-inline V3d normalize(const V3d &v) { return scale(v, 1.0f/length(v)); }
-inline V3d setlength(const V3d &v, float32 l) { return scale(v, l/length(v)); }
-V3d cross(const V3d &a, const V3d &b);
-inline __attribute__((always_inline)) float32 dot(const V3d &a, const V3d &b) { return a.x*b.x + a.y*b.y + a.z*b.z; }
+inline float32 length(const V3d &v) {
+#ifndef DC_SH4
+	return sqrtf(v.x*v.x + v.y*v.y + v.z*v.z);
+#else
+	float len;
+	vec3f_length(v.x, v.y, v.z, len);
+	return len;
+#endif
+}
+inline V3d normalize(const V3d &v) {
+    float invLen; 
+#ifndef DC_SH4
+    invLen = 1.0f / length(v);
+#else
+    invLen = dc::RecipSqrt(fipr_magnitude_sqr(v.x, v.y, v.z, 0.0f));
+#endif
+    return scale(v, invLen); 
+}
+inline V3d setlength(const V3d &v, float32 l) { return scale(v, dc::Div<true, false>(l, length(v))); }
+inline V3d cross(const V3d &a, const V3d &b) {
+    return makeV3d(a.y*b.z - a.z*b.y,
+        a.z*b.x - a.x*b.z,
+        a.x*b.y - a.y*b.x);
+}
+inline float32 dot(const V3d &a, const V3d &b) {
+#ifdef DC_SH4
+	return fipr(a.x, a.y, a.z, 0.0f, b.x, b.y, b.z, 0.0f);
+#else
+	return a.x*b.x + a.y*b.y + a.z*b.z;
+#endif
+}
 inline V3d lerp(const V3d &a, const V3d &b, float32 r){
 	return makeV3d(a.x + r*(b.x - a.x),
 	               a.y + r*(b.y - a.y),
@@ -300,19 +329,52 @@ inline Quat makeQuat(float32 w, const V3d &vec) { Quat q = { vec.x, vec.y, vec.z
 inline Quat add(const Quat &q, const Quat &p) { return makeQuat(q.w+p.w, q.x+p.x, q.y+p.y, q.z+p.z); }
 inline Quat sub(const Quat &q, const Quat &p) { return makeQuat(q.w-p.w, q.x-p.x, q.y-p.y, q.z-p.z); }
 inline Quat negate(const Quat &q) { return makeQuat(-q.w, -q.x, -q.y, -q.z); }
-inline float32 dot(const Quat &q, const Quat &p) { return q.w*p.w + q.x*p.x + q.y*p.y + q.z*p.z; }
+inline float32 dot(const Quat &q, const Quat &p) {
+#ifdef DC_SH4
+	return fipr(q.x, q.y, q.z, q.w, p.x, p.y, p.z, p.w);
+#else
+	return q.w*p.w + q.x*p.x + q.y*p.y + q.z*p.z;
+#endif
+}
 inline Quat scale(const Quat &q, float32 r) { return makeQuat(q.w*r, q.x*r, q.y*r, q.z*r); }
-inline float32 length(const Quat &q) { return sqrtf(q.w*q.w + q.x*q.x + q.y*q.y + q.z*q.z); }
-inline Quat normalize(const Quat &q) { return scale(q, 1.0f/length(q)); }
+inline float32 length(const Quat &q) {
+#ifndef DC_SH4
+	return sqrtf(q.w*q.w + q.x*q.x + q.y*q.y + q.z*q.z);
+#else
+	return dc::Sqrt(fipr_magnitude_sqr(q.x, q.y, q.z, 0.0f));
+#endif
+}
+inline Quat normalize(const Quat &q) {
+    float invLen; 
+#ifndef DC_SH4
+    invLen = 1.0f / length(q);
+#else
+    invLen = dc::RecipSqrt(fipr_magnitude_sqr(q.x, q.y, q.z, 0.0f));
+#endif
+    return scale(q, invLen);
+}
 inline Quat conj(const Quat &q) { return makeQuat(q.w, -q.x, -q.y, -q.z); }
-Quat mult(const Quat &q, const Quat &p);
+inline Quat mult(const Quat &q, const Quat &p) {
+#ifndef DC_SH4
+	return makeQuat(q.w*p.w - q.x*p.x - q.y*p.y - q.z*p.z,
+	                q.w*p.x + q.x*p.w + q.y*p.z - q.z*p.y,
+	                q.w*p.y + q.y*p.w + q.z*p.x - q.x*p.z,
+	                q.w*p.z + q.z*p.w + q.x*p.y - q.y*p.x);
+#else
+    Quat o;
+    dc::quat_mult(reinterpret_cast<dc::quaternion_t *>(&o),
+	              reinterpret_cast<const dc::quaternion_t &>(q),
+                  reinterpret_cast<const dc::quaternion_t &>(p));
+	return o;
+#endif
+}
 inline V3d rotate(const V3d &v, const Quat &q) { return mult(mult(q, makeQuat(0.0f, v)), conj(q)).vec(); }
 Quat lerp(const Quat &q, const Quat &p, float32 r);
 Quat slerp(const Quat &q, const Quat &p, float32 a);

-struct __attribute__((aligned(8))) RawMatrix 
+struct alignas(8) RawMatrixBase
 {
-	V3d right;
+    V3d right;
 	float32 rightw;
 	V3d up;
 	float32 upw;
@@ -320,6 +382,32 @@ struct __attribute__((aligned(8))) RawMatrix
 	float32 atw;
 	V3d pos;
 	float32 posw;
+};
+
+struct RawMatrix: public RawMatrixBase
+{
+    RawMatrix() {}
+
+    RawMatrix(RawMatrixBase &&aggregate):
+	    RawMatrixBase{aggregate}
+    {}
+
+    RawMatrix(const RawMatrix &rhs) {
+        *this = rhs;
+    }
+
+    operator matrix_t *() {
+        return reinterpret_cast<matrix_t *>(this);
+    }
+
+    operator const matrix_t *() const {
+        return reinterpret_cast<const matrix_t *>(this);
+    }
+
+    RawMatrix &operator=(const RawMatrix &rhs) {
+        dc::mat_copy(*this, rhs);
+        return *this;
+    }

 	// NB: this is dst = src2*src1, i.e. src1 is applied first, then src2
 	static void mult(RawMatrix *dst, RawMatrix *src1, RawMatrix *src2);
@@ -327,7 +415,7 @@ struct __attribute__((aligned(8))) RawMatrix
 	static void setIdentity(RawMatrix *dst);
 };

-struct Matrix
+struct alignas(8) MatrixBase
 {
 	enum Type {
 		TYPENORMAL	= 1,
@@ -336,22 +424,65 @@ struct Matrix
 		TYPEMASK = 3
 	};
 	enum Flags {
-		IDENTITY = 0x20000
+		IDENTITY = 0x4,
+        IDENTITY_OLD = 0x20000
 	};
+
+	V3d right;
+	union {
+		struct {
+			uint32_t flags: 3 = TYPEORTHONORMAL|IDENTITY;
+			uint32_t pad0: 29 = 0;
+		};
+		float rightw;
+	};
+	V3d up;
+	union {
+		uint32 pad1;
+		float  upw = 0.0f;
+	};
+	V3d at;
+	union {
+		uint32 pad2;
+		float  atw = 0.0f;
+	};
+	V3d pos;
+	union {
+		uint32 pad3;
+		float  posw = 1.0f;
+	};
+
+    operator matrix_t *() { return reinterpret_cast<matrix_t *>(this); }
+    operator const matrix_t *() const { return reinterpret_cast<const matrix_t *>(this); }
+};
+
+struct Matrix: public MatrixBase
+{
 	struct Tolerance {
 		float32 normal;
 		float32 orthogonal;
 		float32 identity;
 	};

-	V3d right;
-	uint32 flags;
-	V3d up;
-	uint32 pad1;
-	V3d at;
-	uint32 pad2;
-	V3d pos;
-	uint32 pad3;
+    Matrix() {}
+
+    Matrix(MatrixBase &&aggregate){
+        *this = aggregate;
+    }
+
+    Matrix(const Matrix &rhs) {
+        *this = rhs;
+    }
+
+    Matrix &operator=(const RawMatrix &rhs) {
+        dc::mat_copy(*this, rhs);
+        return *this;
+    }
+
+    Matrix &operator=(const MatrixBase &rhs) {
+        dc::mat_copy(*this, rhs);
+        return *this;
+    }

 	static Matrix *create(void);
 	void destroy(void);
@@ -388,10 +519,12 @@ inline void convMatrix(Matrix *dst, RawMatrix *src){

 inline void convMatrix(RawMatrix *dst, Matrix *src){
 	*dst = *(RawMatrix*)src;
+#ifndef DC_SH4
 	dst->rightw = 0.0;
 	dst->upw = 0.0;
 	dst->atw = 0.0;
 	dst->posw = 1.0;
+#endif
 }

 struct Line
--- a/vendor/pvrtex/Makefile
+++ b/vendor/pvrtex/Makefile
@@ -14,8 +14,8 @@ else
 endif

 MYFLAGS=-Wall -Wextra -Wno-unused-parameter -Wno-sign-compare -Ilibavutil -I. -DCONFIG_MEMORY_POISONING=0 -DHAVE_FAST_UNALIGNED=0
-MYCPPFLAGS=$(MYFLAGS)
-MYCFLAGS=$(MYFLAGS) -Wno-pointer-sign
+MYCPPFLAGS=$(MYFLAGS) -std=gnu++17
+MYCFLAGS=$(MYFLAGS) -Wno-pointer-sign -std=gnu17


 .PHONY: all clean
Author	SHA1	Message	Date
Falco Girgis	f0b38ff86f	Saving state. - Good place for gainz in liberty - Cutting off PR - known issue: boats' rotations while driving are incorrect.	2025-05-03 14:24:46 -05:00
Falco Girgis	92f8e04aa7	Saving state, liberty + RW T&L optimizations.	2025-04-29 17:05:15 -05:00
Falco Girgis	f59f84c133	Accelerated lots of RW math + Coronas (liberty) - lot of the RW matrix stuff has become accelerated - went through and accelerated liberty's coronas/reflections ! apparently introduced a bug somewhere along the lines that cause boats to freak out and do summersaults when trying to drive. Will resolve later.	2025-04-29 10:03:49 -05:00
Falco Girgis	c3454ac7ec	Basic rw::Matrix acceleration.	2025-04-27 15:25:54 -05:00
Falco Girgis	056fe39567	Accelerated rw::RawMatrix layer.	2025-04-27 14:59:40 -05:00
Falco Girgis	2361fcd882	Collision + Vec optimizations (liberty only).	2025-04-27 14:40:30 -05:00
Stefanos Kornilios Mitsis Poiitidis	1f2f270da9	Merge branch 'falco/gainz_phase_3' into 'main' Liberty/Miami Perf Gainz Phase 3 See merge request skmp/dca3-game!97	2025-04-26 17:00:40 +00:00
Falco Girgis	1f4dace511	Fixing simulator builds. Ugh. Stupid typo. Accidentally defined mat_store2() to mat_store2(). :(	2025-04-26 11:36:49 -05:00
Falco Girgis	fc4a7e3791	Fixing simulator builds. - Had forgotten to define mat_load2(), mat_store2(), mat_identity2() for simulator builds in rwdc_common.h. - Just defining them to use the regular versions.	2025-04-26 11:28:08 -05:00
Falco Girgis	04b11dfb9d	Synchronized dca3-kos with main branch. Which is where KOS's Gainz Phase 3 work now resides.	2025-04-26 11:15:07 -05:00
Falco Girgis	2716147db4	Everything tested/working: liberty/miami + gainz!	2025-04-26 10:52:36 -05:00
Falco Girgis	dc96ffc551	FIXED LIBERTY BOAT ROTATION ISSUES.	2025-04-15 13:53:19 -05:00
Falco Girgis	4d13e821b5	Liberty/Miami Perf Gainz Phase 3 1) synced dca3-kos repo which has some gainzy commits 2) rwdc_common.h - all low-level and matrix/vector routines for SH4 are now shared in this common file, included in both RW and Liberty/Miami engines 3) CMatrix a. assignment operator: now uses asm-optimized mat_copy() b. multiplication operator: now use mat_mult() SH4 routine c. Scale(): applies a scale matrix via mat_scale d. MultiplyInverse: fipr-optimizations 4) CQuaternion a. multiplication: SH4 ASM FIPR optimized b. Get(V3d& axis, float &angle): fast inversion/division c. Set(RWMatrix&): fast division 5) CVector a. Multiply3x3() now accelerated with mat_transpose 5) RwQuat a. mult(): FIPR accelerated b. length(): FIPR/FSRRA accelerated	2025-04-15 12:06:12 -05:00
Stefanos Kornilios Mitsis Poiitidis	cca7b3c6fa	Merge branch 'falco/gainz_phase_2' into 'main' Phase 2 of Liberty/Mmiami engine+RW SH4 math See merge request skmp/dca3-game!96	2025-04-04 17:26:28 +00:00
Stefanos Kornilios Mitsis Poiitidis	dbcc46b774	Merge branch 'frogbull/fix-not-assigned-hints' into 'main' Fix for "Not assigned" in two Hints (Answer cell phone and replace melee weapon) See merge request skmp/dca3-game!95	2025-04-04 17:25:15 +00:00
Falco Girgis	92be9cfbdd	Phase 2 of liberty/miami engine+RW SH4 math. This PR is the second batch of SH4-specific optimizations for key linear algebra routines within both the RW layer and liberty + miami engines. It also includes a bunch of routines that were accelerated within liberty but were never added to miami, which has helped contribute to the lower perf in Miami. 1) CQuaternion - magnitude, squared magnitude, and dot product use SH4 instructions now 2) CVector (Miami port from Liberty) - Multiply3x3, magnitude, squared magnitude, dot product, distance use SH4 instructions 3) rwbase.h - vector3 length, vector3 dot product, and quaternion dot product use SH4 instructions	2025-04-02 22:39:53 -05:00
Frogbull	7bb5e1640e	Fix for "Not assigned" in two Hints (Answer cell phone and replace melee weapon)	2025-04-01 23:18:19 +02:00
Stefanos Kornilios Mitsis Poiitidis	67cd1106ab	Edit README.md for beta instructions	2025-04-01 05:58:55 +00:00
Stefanos Kornilios Mitsis Poiitidis	a31d73d2e4	Merge branch 'MastaG/c23_host_cc_workaround' into 'main' Build pvrtex using GNU(++) 17 standard See merge request skmp/dca3-game!90	2025-04-01 05:13:53 +00:00
Stefanos Kornilios Mitsis Poiitidis	d630a329a8	Merge branch 'MastaG/dreamshell_prebuilt' into 'main' Add a "dsiso-prebuilt" target See merge request skmp/dca3-game!91	2025-04-01 05:13:24 +00:00
Stefanos Kornilios Mitsis Poiitidis	7cc91b9a76	Merge branch 'frogbull/liberty-hints' into 'main' Small Fix for the Hints of Liberty See merge request skmp/dca3-game!94	2025-03-31 20:30:08 +00:00
Stefanos Kornilios Mitsis Poiitidis	29c346c9ae	Merge branch 'frogbull/miami-hints' into 'main' Miami Hints (Dreamcast Hints for Miami like in Liberty) See merge request skmp/dca3-game!93	2025-03-31 20:29:56 +00:00
Frogbull	186339b854	Fix for the wrong Hint about Submissions and Radio in Liberty (Right/Left was inversed for them)	2025-03-31 16:41:26 +02:00
Frogbull	823a283caa	Dreamcast Hints like in Liberty	2025-03-31 16:32:53 +02:00
MastaG	6dc8ed77bc	Add a "dsiso-prebuilt" target This way DreamShell users can easily build and test using prebuilt elf binaries.	2025-03-31 15:05:06 +02:00
MastaG	a3ce1e4ed5	Build pvrtex using GNU(++) 17 standard Recent distros feature GCC 15 which defaults to C23. ATOMIC_VAR_INIT has been removed from the C23 standard. This causes the following build failure: mem.c:72:39: fout: implicit declaration of function ‘ATOMIC_VAR_INIT’; did you mean ‘ATOMIC_FLAG_INIT’? [-Wimplicit-function-declaration] 72 \| static atomic_size_t max_alloc_size = ATOMIC_VAR_INIT(INT_MAX); \| ^~~~~~~~~~~~~~~ \| ATOMIC_FLAG_INIT Workaround this by telling the compiler to build using GNU 17 instead.	2025-03-30 20:49:27 +02:00